# Unitarity in the non-relativistic regime and implications for dark matter

Marcos M Flores, Kalliopi Petraki

► To cite this version:

Marcos M Flores, Kalliopi Petraki. Unitarity in the non-relativistic regime and implications for dark matter. Physics Letters B, 2024, 858, pp.139022. 10.1016/j.physletb.2024.139022. hal-04576674

HAL Id: hal-04576674

https://hal.science/hal-04576674v1

Submitted on 2 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Letter

# Unitarity in the non-relativistic regime and implications for dark matter

Marcos M. Flores*, Kalliopi Petraki

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>

Laboratoire de Physique de l'École Normale Supérieure, ENS, Université PSL, CNRS, Sorbonne Université, Université Paris Cité, Paris, F-75005, France

<figure><img src="image_5.png" /></figure>

Check for
updates

## ARTICLE INFO

Editor: A. Ringwald

**Keywords:**

Unitarity

Resummation

Imaginary potential

Long-range interactions

Dark matter

## ABSTRACT

Unitarity sets upper limits on partial-wave elastic and inelastic cross-sections, which are often violated by perturbative computations. We discuss the dynamics underlying these limits in the non-relativistic regime, namely long-range interactions, and show how the resummation of the 2-particle-irreducible diagrams arising from squaring inelastic processes unitarizes both elastic and inelastic cross-sections. We provide a simple prescription to obtain the unitarized cross-sections from those that do not include resummation of the squared inelastic processes. Our results are model-independent, apply to all partial waves, and affect elastic and inelastic cross-sections, with extensive implications for new physics scenarios, such as dark-matter freeze-out, indirect detection and self-interactions.

## 1. Introduction

Unitarity has long stood as a fundamental principle in particle physics, crucially invoked in the past to argue for the existence of then-unknown physics, including the electroweak interactions and the Higgs mechanism [1-3]. More recently, it has been extensively employed to circumscribe the range of validity of effective theories aiming to parametrize the still-unknown physics beyond the Standard Model. In this respect, the implications of unitarity have been widely considered in the high-energy limit, relevant for new physics searches at colliders.

In this letter, we consider unitarity in the antipodal limit, the non-relativistic regime [4]. In the forefront of our motivation lies dark matter (DM), which makes up 85% of the non-relativistic matter in our Universe. Besides the DM phenomenology in today's Universe, many important cosmological events, such as the DM production in the early universe, may also fall within the non-relativistic realm. Our results pertain to the cosmology and astrophysics of DM in a variety of theories at the frontiers of research, but are also more generally applicable to any non-relativistic system.

Unitarity sets upper bounds on partial-wave elastic and inelastic cross-sections that indicate the saturation of the corresponding scattering probability to one. In a theory that is by construction unitary, the ramifications of these bounds are two-fold: (a) They may limit unknown physical parameters if phenomenological requirements are imposed that amount to certain processes being sufficiently efficient. (b) They indicate the point where approximations employed to compute scattering

amplitudes fail. While this failure may be simply due to the break-down of a perturbative expansion in a parameter that has become large, in certain cases it hints towards the onset of new effects.

In the present work, we reframe these points, posing the questions: (a') What theories give rise to physical processes that can approach or attain the unitarity limits on scattering cross-sections, thus allowing to approach or attain the bounds on physical parameters set by phenomenological requirements? (b') How can we improve perturbative calculations in weakly-coupled unitary theories to ensure consistency with unitarity? We begin with reviewing the unitarity bounds, before addressing these points and discussing examples in the DM context.

## 2. Off-shell partial-wave optical theorem

To minimize technicalities that are tangential to the main goal of this work, we neglect spin in the following. We consider 2-to-2 transition amplitudes that may in general be off-shell, an essential ingredient in our analysis. We expand them into partial waves as follows,

$$ \mathcal{M}^{ab}(s, \mathbf{k}^a, \mathbf{k}^b) = 16\pi \sum_{\ell} (2\ell + 1) P_{\ell}(\hat{\mathbf{k}}^a \cdot \hat{\mathbf{k}}^b) \mathcal{M}_{\ell}^{ab}(s, |\mathbf{k}^a|, |\mathbf{k}^b|), \quad (1) $$

where a and b denote the initial and final state, s is the first Mandelstam variable, and $\mathbf{k}_a$ stands for the 3-momentum of each particle in the center-of-momentum (CM) frame in the state a. Note that the two interacting particles may in general be different. Next, we define the rescaled partial-wave amplitude

<footnote>* Corresponding author.</footnote>

<footnote>E-mail addresses: marcos.flores@phys.ens.fr (M.M. Flores), kalliopi.petraki@phys.ens.fr (K. Petraki).</footnote>
$$M_{\ell}^{ab}(s, |\mathbf{k}^a|, |\mathbf{k}^b|) = \sqrt{\frac{4|\mathbf{k}^a||\mathbf{k}^b|}{2^{\delta_a} 2^{\delta_b} s}} M_{\ell}^{ab}(s, |\mathbf{k}^a|, |\mathbf{k}^b|), \quad (2)$$

where we have introduced the symmetry parameter $\delta_a = 0$ or 1 if the two particles in the state a are distinguishable or identical respectively.

Following standard methods generalized to off-shell amplitudes [5], S-matrix unitarity implies

$$\mathrm{Im}[M_{\ell}^{aa}(s, |\mathbf{k}|, |\mathbf{k}'|)] = \sum_{b: \text{on-shell}} X_{\ell}^{ab}(s, |\mathbf{k}|, |\mathbf{k}'|), \quad (3)$$

where b runs over all on-shell states into which a can scatter, and

$$X_{\ell}^{ab}(s, |\mathbf{k}|, |\mathbf{k}'|) = \sqrt{\frac{4|\mathbf{k}||\mathbf{k}'|}{2^{2\delta_a} s}} \times \frac{1}{2(8\pi)^2} \int d\Pi_b (2\pi)^4 \delta^4(\{\mathbf{k}\} - \{\mathbf{k}^b\}) \\ \times \left[ \int d\Omega_k Y_{\ell m}^*(\Omega_k) M^{ab}(s, \mathbf{k}, \{\mathbf{k}^b\}) \right] \\ \times \left[ \int d\Omega_{k'} Y_{\ell m}(\Omega_{k'}) M^{ab*}(s, \mathbf{k}', \{\mathbf{k}^b\}) \right]. \quad (4)$$

Equation (4) holds for any $m \in [-\ell, \ell]$, and can also be averaged over $m$. Curly brackets denote collectively the momenta of multi-particle states, and $d\Pi_b$ is the phase-space integration measure for the n-particle state b,

$$d\Pi_b = f_b \prod_{j=1}^{n} \frac{d^3 k_j^b}{(2\pi)^3 2E_j^b}, \quad (5)$$

where $f_b$ is the symmetry factor of the state b that ensures the phase space is not multiply counted. For 2-particle states b, $f_b = 1/2^{\delta_b}$, and Eq. (4) reduces to

$$X_{\ell}^{ab}(s, |\mathbf{k}|, |\mathbf{k}'|) = M_{\ell}^{ab*}(s, |\mathbf{k}|, |\mathbf{k}^b|) M_{\ell}^{ab}(s, |\mathbf{k}'|, |\mathbf{k}^b|), \quad (6)$$

with $|\mathbf{k}^b|$ determined by s due to b being on-shell. Equations (3) to (6) generalize the optical theorem, and we shall employ them in different ways below.

Notably, $X_{\ell}^{ab}(s, |\mathbf{k}|, |\mathbf{k}'|) \ge 0$, which leads to the important observation that the sign of the imaginary part of the elastic amplitude is not arbitrary, $\mathrm{Im}[M_{\ell}^{aa}(s, |\mathbf{k}|, |\mathbf{k}'|)] \ge 0$. This fact will prove critical in the following.

### 3. Unitarity bounds

Next, we focus on on-shell in and out states. We define the cross-section

$$\sigma_{\ell}^{U} \equiv 2^{\delta} 4\pi(2\ell + 1)/k^2, \quad (7)$$

with the symmetry parameter $\delta$ and the CM momentum k referring to the incoming state. Then, all cross-sections can be expressed in the convenient form

$$\sigma_{\ell}^{a \to b} = \sigma_{\ell}^{U} |M_{\ell}^{ab}|^2, \quad (8)$$

where for 2-particle states, $M_{\ell}^{ab}$ is given by Eq. (2), while for states of different multiplicity we set $|M_{\ell}^{ab}|^2 \to X_{\ell}^{ab}(s, |\mathbf{k}|, |\mathbf{k}'|)$, with $|\mathbf{k}|$ fixed by s due to the incoming state being on-shell.

For on-shell states, Eq. (3) implies

$$\mathrm{Im}(M_{\ell}^{aa}) = \sum_{b} |M_{\ell}^{ab}|^2 \ge |M_{\ell}^{aa}|^2. \quad (9)$$

This can be recast for $M_{\ell}^{\mathrm{elas}} = M_{\ell}^{aa}$ as

$$|M_{\ell}^{\mathrm{elas}} - i/2| \le 1/2, \quad (10)$$

which is the familiar unitarity circle centered at i/2 on the complex plane. We may obtain a more restrictive combined bound on elastic

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 1. Upper: The unitarity circle bounding the on-shell partial-wave elastic matrix elements $M_{\ell}^{\mathrm{elas}}$. The border of the circle corresponds to vanishing in-elasticity. Lower: Combined bound on the partial-wave elastic and inelastic cross-sections. The symbols identify points between the two plots. The filled red circle corresponds to maximum inelastic scattering.</figcaption></figure>

and inelastic cross-sections as follows. The inequality Eq. (10) implies $\mathrm{Im}(M_{\ell}^{\mathrm{elas}}) \le |M_{\ell}^{\mathrm{elas}}|$. Then, recalling Eq. (8),

$$\begin{aligned} \left(\sigma_{\ell}^{\mathrm{elas}} + \sigma_{\ell}^{\mathrm{inel}}\right)/\sigma_{\ell}^{U} &= |M_{\ell}^{\mathrm{elas}}|^2 + \sum_{j} |M_{\ell}^{\mathrm{inel}, j}|^2 \\ &= \mathrm{Im}|M_{\ell}^{\mathrm{elas}}| \le |M_{\ell}^{\mathrm{elas}}| = (\sigma_{\ell}^{\mathrm{elas}}/\sigma_{\ell}^{U})^{1/2}, \end{aligned}$$

which can be recast as

$$\sigma_{\ell}^{\mathrm{inel}} / \sigma_{\ell}^{U} \le \sqrt{\sigma_{\ell}^{\mathrm{elas}} / \sigma_{\ell}^{U}} \left(1 - \sqrt{\sigma_{\ell}^{\mathrm{elas}} / \sigma_{\ell}^{U}}\right). \quad (11)$$

The above encompasses the individual bounds $\sigma_{\ell}^{\mathrm{inel}} \le \sigma_{\ell}^{U}/4$ and $\sigma_{\ell}^{\mathrm{elas}} \le \sigma_{\ell}^{U}$. It agrees with [6] except for the symmetry factor $2^{\delta}$ in $\sigma_{\ell}^{U}$ that has been missed in the past. The unitarity circle (10), the combined bound on elastic and inelastic cross-sections (11), and the mapping between them are depicted in Fig. 1.

### 4. The physics of the unitarity limit

#### 4.1. Attaining the unitarity limit

The dependence of the cross-section (7) on the kinematic variables allows us to address question (a'). Let us first motivate this phenomenologically with an example.
Reference [7] argued that the inelastic unitarity bounds imply an upper bound on the mass of thermal-relic DM. The calculation amounted to bounding the so-called canonical DM annihilation cross-section by the inelastic s-wave unitarity limit, $\sigma_{\text{can}} v_{\text{rel}} = 2.2 \times 10^{-26} \text{cm}^3/\text{s} < 4\pi/(m_{\text{DM}}^2 v_{\text{rel}})$, where we set $|k| = m_{\text{DM}} v_{\text{rel}}/2$ in the non-relativistic regime, with $m_{\text{DM}}$ and $v_{\text{rel}}$ being the mass and relative velocity of the annihilating DM particles, and we used the current measurement of the DM density to determine $\sigma_{\text{can}}$. In this comparison, there is clearly a discrepancy between the velocity scalings of $\sigma_{\text{can}}$ and $\sigma_0^U$. Reference [7] claimed that the momentum scaling of $\sigma_\ell^U$ cannot be obtained in the non-relativistic regime (and simply set $v_{\text{rel}}^{\text{freeze-out}} \sim 0.5$ in order to evaluate $\sigma_0^U$), and that s-wave annihilation dominates. This however, is not so, as unitarity implies [4,8] and we shall now discuss. This, in turn, has implications for the upper limit $m_{\text{DM}} < 116$ (82) TeV, deduced according to the considerations of Ref. [7]$^1$ for (non-)self-conjugate thermal-relic DM, but, more importantly, for the underlying physics of viable thermal-relic DM scenarios at the multi-TeV mass scale [4].

In the relativistic regime, $|k| \to \sqrt{s}/2$, thus $\sigma_\ell^U \approx 2^\delta 16\pi(2\ell + 1)/s$, leading to the well-known conclusion that cross-sections should decrease at least as fast as $\sigma \propto 1/s$ at high energies, in order to be consistent with unitarity up to arbitrary high scales. However, for processes computed within an effective theory where a large scale $\Lambda > m$ has been integrated out, with $m$ being the mass of the interacting particles, the cross-sections scale as $\sigma \propto (1/s)(s/\Lambda^2)^\rho$, where $\rho > 0$ depends on the effective theory. Unitarity then implies that the effective theory breaks down at $s \gtrsim \Lambda^2$, where the high-energy physics must be resolved and incorporated in the computation. If doing so realizes the scaling $\sigma \propto 1/s$, then approaching or attaining the unitarity limit in a continuum of energies $s > \Lambda^2$ amounts essentially to the couplings involved being sufficiently large.

If $\Lambda \lesssim m$, then the breakdown of the effective theory and the transition described above occur in the quasi-relativistic or non-relativistic regime. As previously, the scaling $\sigma \propto 1/s$ can retain consistency with unitarity at increasing energies since $s > |k|^2$, thus $\sigma_\ell^U \propto 1/|k|^2 > 1/s$. However, now, this scaling does not allow to approach or attain the unitarity limit in a continuum of low momenta, $|k| \lesssim m$. To identify what can realize the unitarity limit in this case, we repeat the previous considerations in the non-relativistic regime.

For concreteness, we consider a toy model, and follow in part the discussion of Refs. [4,9]. Similar considerations to what follows can be employed in a variety of theories, including gauge theories. Let $\chi$ be a scalar of mass $m$, coupled to a scalar force mediator $\varphi$ of mass $M$ via the trilinear coupling $\mathcal{L} \supset -(y/2)m \chi^2\varphi$, with $y$ being a dimensionless coupling. We define $\alpha \equiv y^2/(16\pi)$ for convenience. The mediator $\varphi$ may also couple in a similar fashion to other lighter species. $\chi\chi$ pairs can scatter elastically or annihilate via $\varphi$ exchanges. If $\varphi$ is sufficiently heavy, then it can be effectively integrated out, giving rise to the $\chi\chi$ elastic and annihilation cross-sections $\sigma_{\text{elas}} \sim \sigma_{\text{ann}} v_{\text{rel}} \sim \alpha^2 s/M^4$, with $v_{\text{rel}}$ being the $\chi\chi$ relative velocity. Setting $s \simeq 4m^2$ and $|k| \simeq m v_{\text{rel}}/2$, unitarity implies, as previously, that this effective description breaks down at sufficiently high energies – with the energy now more conveniently parametrized by $|k|$ or $v_{\text{rel}}$ – in particular at $|k| > M^2/(m\alpha)$. This can be recast in the more indicative form, $M < m\sqrt{\alpha v_{\text{rel}}} \ll m$.

The comparison between scales that arises from this consideration hints to the reason for the apparent violation of unitarity: $\varphi$ mediates a $\chi\chi$ interaction of range $\lambda_{\text{inter}} \sim M^{-1}$ that is comparable to or longer than the geometric mean of the de Broglie wavelength of the interacting particles, $\lambda_{dB} \sim |k|^{-1} = (m v_{\text{rel}}/2)^{-1}$ and their Bohr radius, $\lambda_B \sim (m\alpha/2)^{-1}$. Consequently, virtual $\varphi$ exchanges cause long-range elas-

tic $\chi\chi$ interactions that distort the wavefunction of the incoming $\chi\chi$ pair, and must thus be resummed in order to properly compute both the elastic and inelastic cross-sections. (More details on resummation will be provided in the next section.) In the context of inelastic scattering, this is called the Sommerfeld effect [10,11]. As we discuss next, the resummation of an (attractive) long-range interaction reproduces the same energy dependence as that of the unitarity cross-section (7).

The resummation of the $1/u$-channel $\varphi$ exchanges generates an attractive Yukawa potential between the $\chi$ particles, $V(r) = -\alpha e^{-Mr}/r$ (see e.g. [12]). In the Coulomb limit (which is valid for $M < m v_{\text{rel}}/2$, see e.g. [13]), this results in the cross-sections

$$ \sigma_{\ell}^{\text{elas}} = 2 \frac{4\pi(2\ell + 1)}{\mathbf{k}^2} \sin^2 \left( \frac{1}{2\ell} \ln \left[ \frac{\Gamma(1 + \ell + \ell\zeta)}{\Gamma(1 + \ell - \ell\zeta)} \right] \right), \quad (12\text{a}) $$

$$ \sigma_{\ell}^{\text{ann}} v_{\text{rel}} \sim \frac{\pi\alpha^2}{m^2} v_{\text{rel}}^{2\ell} S_{\ell}(\zeta), \quad (12\text{b}) $$

where $\zeta \equiv \alpha/v_{\text{rel}}$. $S_\ell(\zeta) = S_0(\zeta) \prod_{j=1}^\ell (1 + \zeta^2/j^2)$ with $S_0(\zeta) = 2\pi\zeta/(1 - e^{-2\pi\zeta})$ are the Sommerfeld factors for a Coulomb potential. Expanding for small values of the coupling, $\zeta \ll 1 + \ell$, the above reduce to the well-known Rutherford scattering, $\sigma_\ell^{\text{elas}} \propto \alpha^2/(m^2 v_{\text{rel}}^4)$, and the tree-level annihilation cross-sections, $\sigma_\ell^{\text{ann}} v_{\text{rel}} \sim (\pi\alpha^2/m^2) v_{\text{rel}}^{2\ell}$. On the other hand, for large values of the coupling, $\zeta \gg 1 + \ell$, $\sigma_\ell^{\text{ann}} \propto \alpha^{2\ell+3} \sigma_\ell^U$ while the elastic cross-section oscillates rapidly in the interval $\sigma_\ell^{\text{elas}} \in [0, \sigma_\ell^U]$ with a $\zeta$-dependent phase. The annihilation cross-sections thus scale with momentum in the same way as the unitarity limit,$^2$ which they may approach or attain in a continuum of momenta provided that the coupling is sufficiently large. Conversely, we conclude that the inelastic unitarity limit in the non-relativistic regime can be realized in a continuum of energies only if the incoming particles interact via a long-range force [4,8]. Indeed, an interaction of range comparable to or smaller than $|k|^{-1}$ introduces a new scale – the range of the interaction – that can only spoil the $1/|k|^2$ scaling of the cross-sections.$^3$

Returning to the phenomenological example of DM freeze-out, we may conclude that thermal-relic DM can be as heavy as unitarity allows only if it possesses attractive long-range interactions [4]. Such interactions imply the existence of bound states, and it has been shown that the formation and decay of metastable bound states in the early Universe can deplete DM very efficiently [8], also via higher partial waves [4]. The significant contribution of higher partial waves is in contrast to the direct annihilation processes, for which higher $\ell$ modes are suppressed by $v_{\text{rel}}^{2\ell}$ or $\alpha^{2\ell}$ at tree level and in the Sommerfeld-enhanced regime respectively, as seen above. In the case of bound-state formation, the bound-state wavefunction introduces negative powers of the coupling, compensating for this suppression [15,16]. How many and which partial waves are important depends on the model, there is thus no model-independent upper bound on the mass of thermal-relic DM. Nevertheless, considering carefully the unitarity limit has allowed us to understand and predict, in a model-independent fashion, the underlying physics of viable heavy thermal-relic DM models.

<footnote>$^2$ In Eq. (12b), we have considered only the leading order contribution in $v_{\text{rel}}^2$ for every partial wave. Higher order contributions receive different Sommerfeld corrections that result in the same asymptotic behavior at large $\alpha/v_{\text{rel}}$ [4,14].</footnote>

<footnote>$^3$ It is reasonable to wonder whether other, perhaps exotic, interactions could result in cross-sections that scale as the unitarity limit in the non-relativistic regime. Any exotic mechanism that could generate a Coulomb potential would be encompassed by what we call a long-range force. On the other hand, a different potential cannot result in Coulomb wavefunctions, which is what gives rise to cross-sections with the same momentum dependence as the unitarity limit. This reinforces the above statement. Nevertheless, the unitarity limit may be approached or attained by other types of interactions within a limited energy range, e.g. interactions of finite range at resonant points.</footnote>

<footnote>$^1$ We note that the computation of Ref. [7] of the upper limit on $m_{\text{DM}}$ is in conflict with unitarity. The DM particles are presumed to have a velocity distribution in the early universe. Setting $\sigma_{\text{can}}$ equal to the s-wave inelastic unitarity limit evaluated at a fixed $v_{\text{rel}}$ implies that particles with higher velocity violate unitarity, while particles with lower velocity do not saturate the unitarity limit.</footnote>
## 4.2. Unitarity violation

Our discussion of the dynamics that can approach or attain the (in-elastic) unitarity limit in a continuum of energies in the non-relativistic regime, reframes the issue of (apparent) unitarity violation. In the toy model of Section 4.1, the elastic cross-section (12a) always respects the global limit $\sigma_{\ell}^{\text{elas}} \le \sigma_{\ell}^{U}$. However, the inelastic cross-sections (12b) violate $\sigma_{\ell}^{\text{inel}} \le \sigma_{\ell}^{U}/4$ at large $\alpha$. Moreover the cross-sections (12) violate the combined unitarity bound (11). This is made apparent, for example, in the limit $\sigma_{\ell}^{\text{ann}} \to \sigma_{\ell}^{U}/4$, where the elastic cross-section oscillates in the interval $[0, \sigma_{\ell}^{U}]$, in contradiction with the constraint (11) which mandates that $\sigma_{\ell}^{\text{elas}}$ must approach $\sigma_{\ell}^{U}/4$ as well (cf. Fig. 1). This example represents one of several qualitatively distinct cases of unitarity violation observed in state-of-the-art computations in theories involving long-range interactions:

(i) Attractive Coulomb(-like) potentials violate inelastic unitarity at large couplings. In a QED-like theory, this occurs at $\alpha \gtrsim 0.85$, considering the dominant s- and p-wave inelastic processes for a fermion-antifermion pair [13].

(ii) Radiative transitions between states of different Coulomb(-like) potentials, e.g. bound-state formation with emission of a charged scalar [15,17–19] or a non-Abelian gauge boson [16,20], exhibit non-monotonic dependence on the incoming momentum, with peaks that can violate unitarity for a finite velocity range, even for fairly small values of the couplings. The problem is exacerbated if the potential in the scattering state is repulsive [15,16].

(iii) Potentials generated by the exchange of light but massive media-tors exhibit parametric resonances at the thresholds for the existence of bound levels, violating inelastic unitarity at low $v_{\text{rel}}$, even for small couplings [13,21].

These issues highlight a limitation in our understanding of this regime; they also hamper phenomenological studies of, for example, heavy thermal-relic DM or self-interacting DM, both of which feature long-range interactions. This brings us to question (b') that we now address.

The inequalities (9) cannot be satisfied if the amplitudes are calculated at any finite order in perturbation theory, as the two sides would then be of different order in the couplings of the theory. This suggests that in order to ensure consistency with unitarity, we must resum all interactions of interest. Indeed, returning to the cross-sections (12), the resummation of the elastic $\chi\chi$ interaction (the one-boson exchange), resulted in the elastic cross-section (12a) being consistent with the elastic unitarity bound for all momenta and values of the coupling. On the other hand, the inelastic cross-sections (12b), despite including the resummed one-boson exchanges, also involve a purely inelastic interaction that has not been resummed. However, the optical theorem implies that squaring inelastic interactions generates a contribution to elastic scattering; this must be self-consistently resummed in any computation that invokes these inelastic interactions, in order to ensure consistency with unitarity.

In the following we shall show that the proper resummation of the squared inelastic interactions indeed renders elastic and inelastic cross-sections consistent with the combined unitarity limit (11). We begin with reviewing the resummation process, establishing the appropriate framework.

## 5. Resummation

The resummation of 2-particle interactions amounts to a Dyson-Schwinger equation for the 4-point function of the two interacting particles, with the kernel being the sum of the 2-particle-irreducible (2PI) diagrams. The latter are defined as those that cannot be cut to give two individual contributions to the same 4-point function. The Dyson-

Schwinger equation yields the Bethe-Salpeter equation for the two-particle wavefunction, which reduces to the Schrödinger equation under the instantaneous and non-relativistic approximations (see e.g. [12] for precise definitions and derivations). In momentum space, it reads

$$ \left( \frac{\mathbf{p}^2}{2\mu} - \varepsilon_k \right) \bar{\psi}_k(\mathbf{p}) = \frac{1}{4m_t \mu} \int \frac{d^3 \mathbf{p}'}{(2\pi)^3} \mathcal{K}^{2\text{PI}}(s, \mathbf{p}, \mathbf{p}') \bar{\psi}_k(\mathbf{p}'), \quad (13) $$

with $\varepsilon_k = \mathbf{k}^2/(2\mu)$ being the kinetic energy of the system in the CM frame. Here, $m_t$ and $\mu$ stand for the total and reduced masses of the interacting particles. Equation (13) can be Fourier transformed to position space using

$$ \psi_k(\mathbf{r}) = \int \frac{d^3 p}{(2\pi)^3} e^{i\mathbf{p}\cdot\mathbf{r}} \bar{\psi}_k(\mathbf{p}), \quad (14) $$

$$ \mathcal{V}(\mathbf{r}, \mathbf{r}') = \frac{-1}{4m_t \mu} \int \frac{d^3 p}{(2\pi)^3} \frac{d^3 p'}{(2\pi)^3} e^{i\mathbf{p}\cdot\mathbf{r}} \mathcal{K}^{2\text{PI}}(s, \mathbf{p}, \mathbf{p}') e^{-i\mathbf{p}'\cdot\mathbf{r}'}, \quad (15) $$

where we left the dependence of $\mathcal{V}(\mathbf{r}, \mathbf{r}')$ on $s$ implicit.

Two important remarks are in order. First, no on-shell condition should be imposed on the 2PI kernel, $\mathcal{K}^{2\text{PI}}(s, \mathbf{p}, \mathbf{p}')$, since resummation involves its off-shell iterations of the 2PI diagrams. Note that when on-shell, $|\mathbf{p}| = |\mathbf{p}'|$ are determined by $s$. In fact, in the instantaneous approximation, the dependence of $\mathcal{K}^{2\text{PI}}$ on the energy transfer is neglected (cf. ref. [12]). In position space this yields, in general, non-local potentials $\mathcal{V}(\mathbf{r}, \mathbf{r}')$. Second, guided by the generalized optical theorem (3), we do not restrict our analysis to the case of $\mathcal{K}^{2\text{PI}}$ depending on $\mathbf{p} - \mathbf{p}'$ only, rather than $\mathbf{p}$ and $\mathbf{p}'$ separately. 2PI kernels that depend only on $\mathbf{p} - \mathbf{p}'$ lead to central potentials, where all partial waves are related. We treat instead all partial waves independently, and do not specify the momentum dependence of $\mathcal{K}^{2\text{PI}}$.

In a weakly coupled system, the potential vanishes at $r \to \infty$ as $\mathcal{V}(\mathbf{r}, \mathbf{r}') \sim r^{-(1+\epsilon)}$ with $\epsilon \ge 0$. Then, scattering is described by wavefunctions that asymptote to an incoming plane wave plus an outgoing spherical wave at spatial infinity; in position space

$$ \psi_k(\mathbf{r}) \stackrel{r \to \infty}{\approx} e^{i\mathbf{k}\cdot\mathbf{r}} + f_k(\Omega_r) \frac{e^{ikr}}{r}, \quad (16) $$

where $k = |\mathbf{k}|$, and $f_k(\Omega_r)$ encodes the scattering dynamics. (For $\mathcal{V}(\mathbf{r}, \mathbf{r}') \propto 1/r$, $f_k$ has a mild $r$ dependence at $r \to \infty$, which does not spoil any of our conclusions.)

As is customary, $\psi_k$ and $f_k$ may be expanded into partial waves and re-expressed in terms of phase-shifts,

$$ \psi_k(\mathbf{r}) = \sum_{\ell} (2\ell + 1) P_{\ell}(\hat{\mathbf{k}} \cdot \hat{\mathbf{r}}) \psi_{|\mathbf{k}|, \ell}(r), \quad (17) $$

$$ f_k(\Omega_r) = \sum_{\ell} (2\ell + 1) P_{\ell}(\hat{\mathbf{k}} \cdot \hat{\mathbf{r}}) f_{|\mathbf{k}|, \ell}, \quad (18) $$

$$ f_{|\mathbf{k}|, \ell} = \frac{e^{i2\Delta_{\ell}(k)} - 1}{i2k}, \quad (19) $$

where the phase shifts $\Delta_{\ell}(k)$ are in general complex. The partial-wave cross-sections are (see e.g. [22])

$$ \frac{\sigma_{\ell}^{\text{tot}}}{\sigma_{\ell}^{U}} = k \operatorname{Im} f_{|\mathbf{k}|, \ell} = \frac{1}{2} [1 - \operatorname{Re}(e^{i2\Delta_{\ell}})], \quad (20a) $$

$$ \frac{\sigma_{\ell}^{\text{elas}}}{\sigma_{\ell}^{U}} = k^2 |f_{|\mathbf{k}|, \ell}|^2 = \frac{1}{4} [1 + e^{-4\text{Im}\Delta_{\ell}} - 2\text{Re}(e^{i2\Delta_{\ell}})], \quad (20b) $$

$$ \frac{\sigma_{\ell}^{\text{inel}}}{\sigma_{\ell}^{U}} = \frac{\sigma_{\ell}^{\text{tot}} - \sigma_{\ell}^{\text{elas}}}{\sigma_{\ell}^{U}} = \frac{1}{4} [1 - e^{-4\text{Im}\Delta_{\ell}}]. \quad (20c) $$

The total cross-section is obtained from the standard form of the optical theorem, and Eq. (20c) includes all inelastic channels. The factors $2^{\delta}$ incorporated in $\sigma_{\ell}^{U}$ arise from the (anti)symmetrization of the wavefunction in the case of identical particles that doubles the strength of the $\ell$ modes that survive; for all other $\ell$ modes, $\sigma_{\ell}$ vanish.
<figure><img src="image_8.png" /><figcaption>Fig. 2. The imaginary part of the 2PI kernel arises from squaring inelastic processes (cf. Eq. (21)): the intermediate particles must be on-shell and different from the incoming / outgoing particles. Here we show only 2-to-2 processes; multiparticle intermediate states require also integration over the allowed intermediate phase space.</figcaption></figure>

To ensure that the combined unitarity bound (11) is respected by the above solution, we must still show that $Im\Delta_{\ell} \ge 0$. This is related to the imaginary part of $\mathcal{K}^{2\text{PI}}$ that we discuss next.

## 6. Imaginary potential from generalized optical theorem

The generalized optical theorem (3) suggests that inelastic processes generate an imaginary $\mathcal{K}^{2\text{PI}}$ component. Resumming this contribution renders inelastic cross-sections consistent with unitarity, as we show next.

An imaginary potential due to inelastic processes has been previously considered in Ref. [23], for the purpose of unitarizing s-wave annihilation. However, the connection with the generalized optical theorem was not identified, resulting in ambiguity on the sign of the potential, which is an important issue. Moreover, a central potential was adopted, which does not allow to treat different partial waves independently [24], and the chosen potential, a $\delta$-function in position space, amounted to a specific choice for the underlying inelastic interaction. (See also [25] for a study of similar scope, and [26,27] for works that do not resum inelastic processes but propose related ansätze to improve perturbative amplitudes.) Our treatment relies on resummation, handles all partial waves independently, and does not make model assumptions for the momentum dependence of the resummed process. The sole constraint is the generalized optical theorem (3), as it should be.

To obtain $Im(\mathcal{K}^{2\text{PI}})$, we must remove from the right-hand side of Eq. (3) all elastic contributions (which yield 2-particle-reducible diagrams), and consider only inelastic contributions with any 2PI factors amputated. We denote the latter as $\mathcal{A}^{\text{inel.}j}$. Note that $\mathcal{A}^{\text{inel.}j}$ is not the full amplitude for the inelastic channel $j$; the latter requires resummation of the 2PI diagrams of the incoming particles (cf. e.g. [12]). For simplicity, we will neglect contributions from $2 \to n \ge 3$ scattering.

We partial-wave analyze $\mathcal{K}^{2\text{PI}}$ and $\mathcal{A}^{\text{inel.}j}$ according to Eq. (1), and define their corresponding rescaled versions, that we shall call $\mathcal{K}_{\ell}^{2\text{PI}}$ and $a_{\ell}^{j}$, according to Eq. (2). Then, the generalized optical theorem, Eqs. (3) and (6), imply

$$ \mathrm{Im}[\mathcal{K}_{\ell}^{2\mathrm{PI}}(s, |\mathbf{k}|, |\mathbf{k}'|)] = \sum_{j} a_{\ell}^{j*}(s, |\mathbf{k}|) a_{\ell}^{j}(s, |\mathbf{k}'|), \quad (21) $$

where we do not specify the final-state momenta in the inelastic factors $a_{\ell}^{j}$ since they are on-shell and determined by $s$. Equation (21) is a key element in our analysis. We depict it pictorially in Fig. 2. Examples of contributions to $\mathcal{K}^{2\text{PI}}$ from purely elastic and squared inelastic diagrams are contrasted in Fig. 3.

To proceed, we transcribe Eq. (21) in position space, according to the Fourier transform (15). We analyze $\mathcal{V}(\mathbf{r}, \mathbf{r}')$ in partial waves,

$$ \mathcal{V}(\mathbf{r}, \mathbf{r}') = \sum_{\ell} \frac{2\ell' + 1}{4\pi} P_{\ell}(\hat{\mathbf{r}} \cdot \hat{\mathbf{r}}') \mathcal{V}_{\ell}(\mathbf{r}, \mathbf{r}'), \quad (22) $$

and obtain

$$ \mathrm{Im}[\mathcal{V}_{\ell}(\mathbf{r}, \mathbf{r}')] = - \sum_{j} \nu_{\ell}^{j*}(\mathbf{r}) \nu_{\ell}^{j}(\mathbf{r}'), \quad (23) $$

<figure><img src="image_9.png" /><figcaption>Fig. 3. Examples of contributions to $\mathcal{K}^{2\text{PI}}$ in the toy model discussed in Section 4. Top: The one-boson exchange diagram is a purely elastic interaction that generates the leading order contribution to $\mathrm{Re}(\mathcal{K}^{2\text{PI}})$. Bottom: 2PI diagrams arising from squaring the 2-to-2 inelastic processes, $\chi\chi \to \phi\phi$ and $\chi\chi \to \mathcal{B}(\chi\chi) + \phi$, where $\mathcal{B}(\chi\chi)$ is bound state, described at leading order by an infinite ladder of one-boson exchange diagrams. These diagrams generate subleading contributions to $\mathrm{Re}(\mathcal{K}^{2\text{PI}})$ from off-shell intermediate states, and leading contributions to $\mathrm{Im}(\mathcal{K}^{2\text{PI}})$ from on-shell intermediate states. The latter are described by Eq. (21) (see text for explanations).</figcaption></figure>

with

$$ \nu_{\ell}^{j}(\mathbf{r}) = \left(\frac{2^{1+\delta}}{\pi^2 \mu}\right)^{1/2} \int_{0}^{\infty} dp \, p^{3/2} j_{\ell}(pr) a_{\ell}^{j}(s, p), \quad (24) $$

where we set $\sqrt{s} \approx m_t$ in the prefactor, consistently with the non-relativistic approximation leading to Eq. (15) [12].

We shall employ Eq. (23) in Section 7, solve Schrödinger's equation and obtain the corrected wavefunction that, as we show, renders elastic and inelastic cross-sections consistent with unitarity. Before doing so, we offer here a simpler proof of this unitarization procedure, valid at leading order.

We consider the current $\mathbf{j}_k(\mathbf{r}) \equiv \mathrm{Im}[\psi_k^*(\mathbf{r})\nabla\psi_k(\mathbf{r})]/\mu$. By virtue of Stokes' theorem and the continuity equation,

$$ \begin{aligned} \sigma_{\text{inel}} &= -\frac{\mu}{k} \int d^3\mathbf{r} \nabla \cdot \mathbf{j}_k(\mathbf{r}) \\ &= 2 \operatorname{Im} \left[ \int d^3\mathbf{r} d^3\mathbf{r}' \psi_k^*(\mathbf{r}) \mathcal{V}(\mathbf{r}, \mathbf{r}') \psi_k(\mathbf{r}') \right]. \end{aligned} \quad (25) $$

Using the first Born approximation for the right-hand side, and expanding correspondingly the left-hand side in $\mathrm{Im}\Delta_{\ell}$ (cf. Eq. (20c)), keeping lowest order terms, we find

$$ \sum_{\ell} \sigma_{\ell}^{U} \mathrm{Im}\Delta_{\ell}(k) \approx \frac{1}{2m_{\ell}k} \mathrm{Im}[\mathcal{K}^{2\mathrm{PI}}(s, k, k)], \quad (26) $$

which implies $\mathrm{Im}\Delta_{\ell} \ge 0$ as per Eq. (21). (This remains so if we include $2 \to n \ge 3$ contributions to $\mathrm{Im}[\mathcal{K}^{2\text{PI}}]$.)

The above analysis, distilled in Eqs. (20), (21) and (26), shows that the resummation of the imaginary part of the 2PI diagrams — which arises from squaring inelastic processes, an essential point — ensures that the elastic and inelastic cross-sections respect the combined unitarity bound (11). While this proof is at leading order in the couplings that give rise to inelastic scattering, and concerns the sum of all inelastic processes, we show in the following how the resummation of Eq. (21) ensures that unitarity is respected at all orders, and obtain unitarized solutions for the exclusive inelastic cross-sections.

## 7. Unitarization

Next, we compute how the consistent resummation of the inelastic contributions to elastic scattering, and in particular of Eq. (21), unitarizes the partial-wave elastic and exclusive inelastic processes.

To begin, we note that the elastic scattering cross-section will be computed according to Eq. (20b), using the phase shifts that take into account the imaginary part of the potential given by Eq. (23). The cross-section for an inelastic channel $j$ is $\sigma_{\ell}^{\text{inel.}j} = \sigma_{\ell}^{U} |M_{\ell}^{\text{inel.}j}|^2$, where $M_{\ell}^{\text{inel.}j}$
stands for the corresponding amplitude rescaled according to Eq. (2), and incorporates both the hard-scattering inelastic amplitude, $a_{\ell}^{j}$, and the effect of the non-relativistic potential on the incoming particles,

$$ 
\begin{aligned}
M_{\ell}^{\text{inel},j}(k) &= \frac{k^{1/2}}{2\pi^2} \int_0^{\infty} dp \, p^{3/2} \bar{\psi}_{|k|,\ell}(p) a_{\ell}^{j}(s,p) \\
&= \sqrt{\frac{2\mu k}{2\delta}} \int_0^{\infty} dr \, r^2 \, \psi_{|k|,\ell}(r) v_{\ell}^{j}(r).
\end{aligned} 
\quad (27) $$

We must now determine the wavefunctions $\psi_{|k|,\ell}$.

### 7.1. Schrödinger's equation with complex potential

We shall consider a complex potential whose real part we take to be central for simplicity, Re[$V(\mathbf{r}, \mathbf{r}')$] = $V(r)\delta^3(\mathbf{r} - \mathbf{r}')$, while its imaginary part is given by Eq. (23). Then, the partial-wave Schrödinger equation for $u_{|k|,\ell}(r) \equiv kr\psi_{|k|,\ell}(r)$ reads

$$ S_{\ell}(r)u_{k,\ell}(r) - i \sum_{j} r v_{\ell}^{j*}(r) \int_{0}^{\infty} dr' r' v_{\ell}^{j}(r') u_{k,\ell}(r') = \mathcal{E}_{k} u_{k,\ell}(r), \quad (28) $$

where $S_{\ell}(r)$ is the differential operator that includes the real part of the potential only,

$$ S_{\ell}(r) = -\frac{1}{2\mu} \frac{d^2}{dr^2} + \frac{\ell(\ell+1)}{2\mu r^2} + V(r), \quad (29) $$

and $\mathcal{E}_k = k^2/(2\mu)$, as before. We define two versions of the amplitudes and cross-sections: the *unregulated* ones, for which we consider the real but neglect the imaginary part of the potential, and the *regulated* ones, for which we take into account the entire complex potential. We denote them by the subscripts 'unreg' and 'reg' respectively, and aim to express the latter in terms of the former.

To do so, we consider the solution of Schrödinger's equation in the former case and the Green's function,

$$ [S_{\ell}(r) - \mathcal{E}_{k}] \mathcal{F}_{k,\ell}(r) = 0, \quad (30) $$

$$ [S_{\ell}(r) - \mathcal{E}_{k}] G_{k,\ell}(r) = 0, \quad (31) $$

$$ [S_{\ell}(r) - \mathcal{E}_{k}] G_{k,\ell}(r, r') = \delta(r - r'), \quad (32) $$

where $\mathcal{F}_{k,\ell}$ and $G_{k,\ell}$ are the regular and irregular families of solutions that feature the following asymptotic behaviors at $r \to \infty$,

$$ \mathcal{F}_{k,\ell}(r \to \infty) \mapsto +\frac{1}{i2} (e^{ikr} e^{i2\theta_{\ell}(k)} - e^{-i(kr-\ell\pi)}), \quad (33a) $$

$$ G_{k,\ell}(r \to \infty) \mapsto -\frac{1}{2} (e^{ikr} e^{i2\theta_{\ell}(k)} + e^{-i(kr-\ell\pi)}), \quad (33b) $$

with $\theta_{\ell}(k) \in \mathbb{R}$. For some well-known real potentials, such as the Coulomb potential, $\mathcal{F}_{k,\ell}$ and $G_{k,\ell}$ are known functions [28,29]. The functions $\mathcal{F}_{k,\ell}(r)$ describe states of an incident plane wave plus an outgoing spherical wave, and are the solutions of interest. $G_{k,\ell}(r)$ will be useful for computational purposes. We shall require that $G_{k,\ell}(r, r')$ acts as an outgoing spherical wave at $r \to \infty$.

We may now write an implicit solution of Eq. (28) for $u_{k,\ell}(r)$ for the full potential,

$$ u_{k,\ell}(r) = \mathcal{F}_{k,\ell}(r) + i \sqrt{\frac{2\delta k}{2\mu}} \sum_j M_{\ell,\text{reg}}^{\text{inel},j}(k) \left( \int_0^\infty dr' r' G_{k,\ell}(r, r') v_{\ell}^{j*}(r') \right), \quad (34) $$

where we used

$$ M_{\ell,\text{reg}}^{\text{inel},j}(k) = \sqrt{\frac{2\mu}{2\delta k}} \int_0^\infty dr r u_{k,\ell}(r) v_{\ell}^{j}(r), \quad (35a) $$

$$ M_{\ell,\text{unreg}}^{\text{inel},j}(k) = \sqrt{\frac{2\mu}{2\delta k}} \int_0^\infty dr r \mathcal{F}_{k,\ell}(r) v_{\ell}^{j}(r), \quad (35b) $$

to express the r-independent factor of the second term in terms of the rescaled regulated inelastic amplitudes, which depend on the wavefunction $u_{k,\ell}$. The unregulated amplitudes will become useful shortly.

Inserting Eq. (34) into (35a), and introducing the matrix

$$ [\mathbb{N}_{\ell}(k)]^{ij} \equiv \delta^{ij} - i \int_0^{\infty} dr r \int_0^{\infty} dr' r' v_{\ell}^{j*}(r) G_{k,\ell}(r, r') v_{\ell}^{j}(r'), \quad (36) $$

we obtain a matrix equation between regulated and unregulated inelastic amplitudes. Upon inversion, it gives

$$ M_{\ell,\text{reg}}^{\text{inel},j}(k) = \sum_j [\mathbb{N}_{\ell}^{-1}(k)]^{ij} M_{\ell,\text{unreg}}^{\text{inel},j}(k). \quad (37) $$

Equation (37) is the key result of this work: it expresses the regulated inelastic amplitudes in terms of the unregulated ones. Similarly, using Eq. (37), the solution (34) can be re-expressed in terms of input parameters only: the imaginary potential and unregulated wavefunctions.

$$ u_{k,\ell}(r) = \mathcal{F}_{k,\ell}(r) + i \sqrt{\frac{2\delta k}{2\mu}} \sum_{i,j} M_{\ell,\text{unreg}}^{\text{inel},j}(k) [\mathbb{N}_{\ell}^{-1}(k)]^{ij} \times \int_0^{\infty} dr' r' G_{k,\ell}(r, r') v_{\ell}^{j*}(r'). \quad (38) $$

### 7.2. Green's function

To further simplify the results (37) and (38), we focus on the Green's function $G_{k,\ell}(r, r')$. We seek an expression for $G_{k,\ell}(r, r')$ that is regular at $r \to 0$ and has the asymptotic behavior of an outgoing spherical wave at $r \to \infty$. Following Refs. [24,29–32], we expand the r dependence of $G_{k,\ell}(r, r')$ in terms of $\mathcal{F}_{q,\ell}(r)$, which, being eigenfunctions of the hermitian operator $S_{\ell}(r)$, constitute a complete orthonormal set of states with the desired asymptotic behavior,<fnref n="4" />

$$ G_{k,\ell}(r, r') = \int_0^\infty dq \, g_{k,\ell}(q, r') \mathcal{F}_{q,\ell}(r). \quad (39) $$

With the above decomposition, Eq. (32) becomes

$$ \begin{aligned} \delta(r-r') &= \int_0^\infty dq \, g_{k,\ell}(q,r') [S_\ell(r) - \mathcal{E}_k] \mathcal{F}_{q,\ell}(r), \\ &= \int_0^\infty dq \, g_{k,\ell}(q,r') \left( \frac{q^2-k^2}{2\mu} \right) \mathcal{F}_{q,\ell}(r). \end{aligned} \quad (40) $$

We must now invert the above to obtain $g_{k,\ell}(q, r')$. For this, we need the orthonormality relation for $\mathcal{F}_{q,\ell}(r)$,

$$ \int_0^\infty dr \, \mathcal{F}_{q,\ell}(r) \mathcal{F}_{q',\ell}^*(r) = \frac{2^\delta \pi}{2} \delta(q-q'). \quad (41) $$

<footnote n="4">⁴ This requires that the operator $S_{\ell}(r)$ does not depend on the energy eigenvalue, which in turn implies that any dependence on s has been eliminated by setting $s \approx m_{\ell}^2$. Despite this non-relativistic approximation, the solutions to Eqs. (30) to (32) can be defined for arbitrarily large energies. While such solutions may not well approximate the actual wavefunctions in the relativistic regime, they will be useful for computational purposes (cf. Footnote 5).</footnote>
Considering this, Eq. (40) yields

$$g_{k,\ell}(q,r) = \frac{4\mu}{2\delta\pi} \frac{\mathcal{F}_{q,\ell}^*(r)}{q^2 - k^2}. \qquad (42)$$

Then, Eq. (39) becomes

$$G_{k,\ell}(r, r') = \frac{4\mu}{2\delta\pi} \int_0^\infty dq \frac{\mathcal{F}_{q,\ell}(r) \mathcal{F}_{q,\ell}^*(r')}{q^2 - k^2 - \frac{1}{2}\epsilon}. \qquad (43)$$

In order to evaluate $G_{k,\ell}(r, r')$, it will be convenient to extend the range of integration over $dq$ to $(-\infty, \infty)$. To analytically continue to $q < 0$, we consider the asymptotic behaviors (33), and set self-consistently

$$\theta_{\ell}(-q) = -\theta_{\ell}(q) + \ell\pi, \qquad (44a)$$

$$\mathcal{F}_{-q,\ell}(r) = -\mathcal{F}_{q,\ell}^*(r), \qquad (44b)$$

$$G_{-q,\ell}(r) = G_{q,\ell}^*(r). \qquad (44c)$$

With this, Eq. (43) becomes

$$G_{k,\ell}(r, r') = \frac{2\mu}{2\delta\pi} \int_{-\infty}^{\infty} dq \frac{\mathcal{F}_{q,\ell}(r) \mathcal{F}_{q,\ell}^*(r')}{q^2 - k^2 - \frac{1}{2}\epsilon}, \qquad (45)$$

with the prescription $\epsilon \to 0^+$ chosen such that $G_{k,\ell}(r, r')$ behaves as an outgoing spherical wave at $r \to \infty$, as will be shown below. Since $\mathcal{F}_{q,\ell}(r)$ and $G_{q,\ell}(r)$ are real up to the same $\ell$- and $q$-dependent but $r$-independent phase, it follows that

$$\mathcal{F}_{q,\ell}(r) \mathcal{F}_{q,\ell}^*(r') = \mathcal{F}_{q,\ell}^*(r) \mathcal{F}_{q,\ell}(r'), \qquad (46a)$$

$$G_{q,\ell}(r) G_{q,\ell}^*(r') = G_{q,\ell}^*(r) G_{q,\ell}(r'), \qquad (46b)$$

$$\mathcal{F}_{q,\ell}(r) \mathcal{F}_{q,\ell}^*(r') = \mathcal{F}_{q,\ell}^*(r) \mathcal{F}_{q,\ell}(r'). \qquad (46c)$$

Here, this implies $G_{k,\ell}(r, r') = G_{k,\ell}(r', r)$, as expected.

To evaluate Eq. (45), we consider the two independent eigenfunction families of $S_\ell$, characterized by the asymptotic behaviors of Eqs. (33). We then define

$$\mathcal{H}_{k,\ell}^{(\pm)}(r) = \mathcal{F}_{k,\ell}(r) \pm i \mathcal{G}_{k,\ell}(r). \qquad (47)$$

The functions $\mathcal{H}_{k,\ell}^{(\pm)}(r)$ have the properties

$$\mathcal{H}_{k,\ell}^{(+)}(r \to \infty) \to -i e^{+i(kr+2\theta_{\ell})}, \qquad (48a)$$

$$\mathcal{H}_{k,\ell}^{(-)}(r \to \infty) \to +i e^{-i(kr-\ell\pi)}, \qquad (48b)$$

$$\mathcal{H}_{-k,\ell}^{(\pm)}(r) = -\mathcal{H}_{k,\ell}^{(\pm)*}(r). \qquad (48c)$$

We invert Eq. (47) to express $\mathcal{F}_{k,\ell}(r)$ in terms of $\mathcal{H}_{k,\ell}^{(\pm)}(r)$, which allows to separate Eq. (45) into two integrals,

$$G_{k,\ell}(r, r') = \frac{\mu}{2\delta\pi} \int_{-\infty}^{\infty} dq \frac{\mathcal{F}_{q,\ell}^*(r) \mathcal{H}_{q,\ell}^{(+)}(r')}{q^2 - k^2 - \frac{1}{2}\epsilon} \\ + \frac{\mu}{2\delta\pi} \int_{-\infty}^{\infty} dq \frac{\mathcal{F}_{q,\ell}^*(r) \mathcal{H}_{q,\ell}^{(-)}(r')}{q^2 - k^2 - \frac{1}{2}\epsilon}. \qquad (49)$$

We evaluate the above via contour integration, using the residue theorem, and choosing the contours based on the asymptotic behavior of $\mathcal{H}_{k,\ell}^{(\pm)}(r)$. For $r < r'$, the first integral must be evaluated in the upper-half plane while the second in the lower half plane. The case $r > r'$ follows a similar logic. Considering also the analytic continuations (44b) and (48c), as well as Eqs. (46), we find,

$$G_{k,\ell}(r, r') = +\frac{2\mu i}{2\delta k} \mathcal{F}_{k,\ell}^*(r_<(r)) \mathcal{H}_{k,\ell}^{(+)}(r_>), \qquad (50)$$

where $r_< \equiv \min\{r, r'\}$ and $r_> \equiv \max\{r, r'\}$. The necessary asymptotic behavior of $G_{k,\ell}(r, r')$ is made apparent by the appearance of $\mathcal{H}_{k,\ell}^{(+)}$ in the above expression.

If the spectrum of $S_\ell$ includes bound states, they should be included in the Green's function spectral decomposition of Eq. (39). However, in this case, the scattering state wavefunctions $\mathcal{F}_{q,\ell}$ exhibit poles at the imaginary values of $q$ that yield the bound-state energies. Upon contour integration in Eq. (49), the residues of those poles cancel exactly the bound-state contributions from the spectral decomposition (cf. e.g. [33]), thereby leaving the final result, Eq. (50), unaffected.

### 7.3. Regulated cross sections

The results of the previous section can be used to derive simple expressions for the regulated elastic and inelastic cross sections in terms of the corresponding unregulated quantities. To begin, using Eq. (45), we can re-express the normalization matrix, Eq. (36), as

$$[\mathbb{N}_{\ell}(k)]^{ij} = \delta^{ij} - \frac{i}{\pi} \int_{-\infty}^{\infty} \frac{h_{\ell}^{i}(q) \bar{h}_{\ell}^{j}(q)}{q^2 - k^2 - \frac{1}{2}\epsilon} dq, \qquad (51)$$

where, for $q \in \mathbb{R}$, we define

$$h_{\ell}^{i}(q) \equiv \sqrt{\frac{2\mu}{2\delta}} \int_{0}^{\infty} dr r \mathcal{F}_{q,\ell}(r) v_{\ell}^{i}(r), \qquad (52a)$$

$$\bar{h}_{\ell}^{i}(q) \equiv \mathrm{Re}[h_{\ell}^{i}(q)] - i \mathrm{Im}[h_{\ell}^{i}(q)]. \qquad (52b)$$

Equations (44b) and (46a) imply that $h_{\ell}^{i}(-q)\bar{h}_{\ell}^{j}(-q) = h_{\ell}^{i}(q)\bar{h}_{\ell}^{j}(q)$.

With these definitions, the integrand in Eq. (51) can be analytically continued in the complex $q$ plane, as it is an analytic function of $q$ except perhaps for isolated poles (meromorphic). If the functions $h_{\ell}^{i}(q)$ are holomorphic in the upper half of the complex plane and increase at complex infinity more slowly than $|q|^{1/2}$, then we can evaluate the contour integral to obtain,<fnref n="5" />

$$[\mathbb{N}_{\ell}(k)]^{ij} = \delta^{ij} + M_{\ell,\text{unreg}}^{\text{inel},i}(k) M_{\ell,\text{unreg}}^{\text{inel},j*}(k). \qquad (53)$$

From here, it is easy to demonstrate that the amplitudes $M_{\ell,\text{unreg}}^{\text{inel},j}(k)$ form an eigenvector of $\mathbb{N}_\ell(k)$,

$$\sum_j [\mathbb{N}_\ell]^{ij} M_{\ell, \text{unreg}}^{\text{inel},j} = M_{\ell, \text{unreg}}^{\text{inel},i} \left( 1 + \sigma_{\ell, \text{unreg}}^{\text{inel}} / \sigma_{\ell}^{U} \right). \qquad (54)$$

The above implies that the regulated amplitude (37) is

$$M_{\ell, \text{reg}}^{\text{inel},i} = \frac{M_{\ell, \text{unreg}}^{\text{inel},i}}{1 + \sigma_{\ell, \text{unreg}}^{\text{inel}} / \sigma_{\ell}^{U}}, \qquad (55)$$

from where which we obtain the regulated exclusive partial-wave inelastic cross-sections, which will be presented shortly.

Next, we want to determine the elastic cross-section. For this purpose, we need to deduce the complex phase shift $\Delta_\ell(k)$ that determines

<footnote n="5">⁵ For physical, non-relativistic momenta, $0 < q \ll \mu$, the functions $h_{\ell}^{i}(q)$ are related to the unregulated inelastic amplitudes given by Eq. (35b), $h_{\ell}^{i}(q) = \sqrt{q} M_{\ell, \text{unreg}}^{\text{inel},i}(q)$. However, $M_{\ell, \text{unreg}}^{\text{inel},i}(q)$ cannot be analytically continued on the entire real $q$ axis, due to the rescaling by the non-analytic factor $\sqrt{q}$ relative to the full amplitudes in Eq. (2). In Eq. (51), we use $h_{\ell}^{i}(q)$ instead of $M_{\ell, \text{unreg}}^{\text{inel},i}(q)$ to make analyticity manifest. To ensure analyticity in the complex plane, we also use $\bar{h}_{\ell}^{i}(q)$, which depends only on $q$, instead of using $[h_{\ell}^{i}(q)]^* = \bar{h}_{\ell}^{i}(q^*)$, which introduces dependence on $q^*$. Finally, $h_{\ell}^{i}(q)$ is defined via Eq. (52a) for arbitrarily large $|q|$, while Eq. (35b) holds only for non-relativistic momenta (cf. Footnote 4). The large-momentum modes are needed for the spectral analysis of the Green's function, employed in Eq. (51).</footnote>
the asymptotic behavior of the partial-wave modes at $r \to \infty$, according to

$$u_{k,\ell}(r \to \infty) \mapsto \frac{1}{\bar{i}2} (e^{i[kr+2\Delta_{\ell}(k)]} - e^{-i(kr-\ell\pi)}). \qquad (56)$$

We expand our solution, Eq. (34), in the $r \to \infty$ limit, using the Green's function (50) and considering the asymptotic forms (33a) and (48). We find the asymptotic form (56) with the phase shift being $\Delta_{\ell} = \theta_{\ell} + \delta_{\ell}$, where

$$\begin{aligned} e^{i2\delta_{\ell}(k)} &= 1 - 2 \sum_{l,j} M_{\ell, \text{unreg}}^{\text{inel}, i*}(k) [\mathbb{N}_{\ell}^{-1}(k)]^{ij} M_{\ell, \text{unreg}}^{\text{inel}, j}(k) \\ &= 1 - 2 \left( \frac{\sigma_{\ell, \text{unreg}}^{\text{inel}} / \sigma_{\ell}^{U}}{1 + \sigma_{\ell, \text{unreg}}^{\text{inel}} / \sigma_{\ell}^{U}} \right), \end{aligned} \qquad (57)$$

where the second line has been obtained under the analyticity and convergence assumptions leading to Eq. (53). Note that $\text{Re}(\delta_{\ell}) = 0$ or $\pi/2$ for $\sigma_{\ell, \text{unreg}}^{\text{inel}} < \sigma_{\ell}^{U}$ or $> \sigma_{\ell}^{U}$ respectively. From Eq. (57), and given that $\sigma_{\ell}^{\text{tot}} / \sigma_{\ell}^{U} = k \text{Im}(f_{|\mathbf{k}|, \ell})$, or

$$\sigma_{\ell, \text{unreg}}^{\text{elas}} / \sigma_{\ell}^{U} = \frac{1}{2} (1 - \text{Re}[e^{i2\theta_{\ell}(k)}]), \qquad (58a)$$

$$(\sigma_{\ell, \text{reg}}^{\text{elas}} + \sigma_{\ell, \text{reg}}^{\text{inel}}) / \sigma_{\ell}^{U} = \frac{1}{2} (1 - \text{Re}[e^{i2\Delta_{\ell}(k)}]), \qquad (58b)$$

we can obtain the regulated elastic cross section. Defining

$$x_{\ell, \text{(un)reg}} \equiv \sigma_{\ell, \text{(un)reg}}^{\text{elas}} / \sigma_{\ell}^{U}, \qquad (59a)$$

$$y_{\ell, \text{(un)reg}}^{j} \equiv \sigma_{\ell, \text{(un)reg}}^{\text{inel}, j} / \sigma_{\ell}^{U}, \qquad (59b)$$

$$y_{\ell, \text{(un)reg}} \equiv \sum_{j} y_{\ell, \text{(un)reg}}^{j}, \qquad (59c)$$

we find, using Eqs. (55) and (58),

$$x_{\ell, \text{reg}} = \frac{x_{\ell, \text{unreg}} + (1 - x_{\ell, \text{unreg}}) y_{\ell, \text{unreg}}^2}{(1 + y_{\ell, \text{unreg}})^2}, \qquad (60a)$$

$$y_{\ell, \text{reg}}^{j} = \frac{y_{\ell, \text{unreg}}^{j}}{(1 + y_{\ell, \text{unreg}})^{2}}. \qquad (60b)$$

We reiterate that both the regulated and unregulated cross-sections include the effect of the real potential.

Equations (60) are manifestly consistent with unitarity: $x_{\ell, \text{reg}} \le 1$ and $y_{\ell, \text{reg}}^{j} \le y_{\ell, \text{reg}} \le 1/4$, and the combined bound (11) is also respected. Note that $x_{\ell, \text{unreg}} \le 1$ since it already includes resummation of $\text{Re}(K^{2\text{PI}})$. The maximum inelastic cross-section is reached for $y_{\ell, \text{unreg}} = 1$, where, consistently with (11), $x_{\ell, \text{reg}} = y_{\ell, \text{reg}} = 1/4$, irrespectively of $x_{\ell, \text{unreg}}$ (cf. Fig. 1). We plot Eqs. (60) in Fig. 4.

## 8. Phenomenological implications

The ramifications of the regularization Eqs. (60) (or of the more general Eqs. (37), (38) and (51)) are far-reaching:

* The resummation of the squared inelastic diagrams affects both the inelastic and elastic cross-sections rather significantly. It is thus expected to affect the DM density in various production mechanisms, the DM indirect detection signals, and the DM self-interaction cross sections; the latter are important for the galactic structure (see [34,35] for reviews).

* Each inelastic channel is regulated by the *inclusive* inelastic cross-section. It follows that even if a certain inelastic channel appears seemingly irrelevant for a particular purpose due to its nature, it may in fact be very important, provided that it is sufficiently strong.

<figure><img src="image_10.png" /><figcaption>Fig. 4. The regulated elastic (colored lines) and 2-to-2 inclusive inelastic (black line) partial-wave cross-sections, vs the unregulated 2-to-2 inclusive inelastic cross-section.</figcaption></figure>

This is because such a process regulates the pertinent inelastic processes.

For example, in DM freeze-out, metastable bound-state formation processes that are comparable to or faster than direct annihilation into radiation can significantly renormalize down the strength of the latter (cf. Ref. [18, Fig. 11]), even when they themselves cannot deplete DM efficiently due to their rapid dissociation by the radiation bath. This may suppress the DM depletion and modify predictions.

Similarly, for DM indirect detection, rapid inelastic processes with (nearly) unobservable products can renormalize down inelastic processes that produce observable signals, thereby relaxing constraints.

* For $y_{\ell, \text{unreg}} \ll 1$, the relative corrections to the scattering rates are $\delta x_{\ell} / x_{\ell} \approx \delta y_{\ell}^{j} / y_{\ell}^{j} \approx -2y_{\ell, \text{unreg}}$.

For DM freeze-out, this implies that the imaginary potential affects the DM density at a level similar to or larger than the experimental uncertainty if $y_{\ell, \text{unreg}} \gtrsim \mathcal{O}(10^{-2})$, i.e., even when the annihilation cross-section is far below the unitarity limit.

* At the opposite limit, $y_{\ell, \text{unreg}} \gg 1$, we find $x_{\ell, \text{reg}} \approx 1 - x_{\ell, \text{unreg}}$ and $y_{\ell, \text{reg}} \approx 1/y_{\ell, \text{unreg}}$, i.e., perhaps counter-intuitively, the regulated cross-sections decrease as the unregulated ones increase. (See below for an interpretation.) Note that $y_{\ell, \text{unreg}} \gg 1$ can occur even for (very) small couplings, e.g. at resonant points, as well as in bound-state formation processes where the potentials of the incoming scattering and outgoing bound state are different [15–18]. The small couplings ensure the validity of the instantaneous and non-relativistic approximations, and thus of the framework for the computations of this work.

With respect to the DM density, this implies that there could be two branches of solutions (mass-coupling correlation) that yield the observed value.

* If $\sigma_{\ell, \text{unreg}}^{\text{inel}}$ scales with momentum in the same way as $\sigma_{\ell}^{U}$, then this feature is retained by the regularization. If not, the regularization generates a non-trivial velocity dependence. The earlier conclusion thus remains: the unitarity limits can be attained in a continuum of non-relativistic momenta only by long-range interactions.

* The upper bounds implied by unitarity on the mass of frozen-out DM annihilating via a single partial wave $\ell$, or all partial waves in the range $0 \le \ell \le L$, are, for self-conjugate and non-self-conjugate DM,
<table><tr><td>$m_{DM}^U/(140 \text{ TeV})$</td><td>solely $\ell$</td><td>$0 \le \ell \le L$</td></tr><tr><td>non-self-conjugate DM</td><td>$\sqrt{1+2\ell}$</td><td>$(1+L)$</td></tr><tr><td>self-conjugate DM</td><td>$2\sqrt{1+2\ell}$</td><td>$2\sqrt{(1+L)(1+L/2)}$</td></tr></table>

where we used the numerical results of [4,8], assuming for definiteness that DM annihilates into plasma of the same temperature as that of the Standard Model, and the number of extra relativistic degrees of freedom is negligible. It is straightforward to adjust these assumptions; this would affect the overall mass scale (here, 140 TeV), but not the scaling with $\ell$ or $L$. Several remarks are in order:

- The difference between the figures of 140 TeV denoted above, and 82 TeV deduced from [7] (cf. Section 4) is due to the momentum scaling of $\sigma_{\ell}^U$ that has been here properly taken into account [4,8].

- For self-conjugate DM, the above results incorporate the factor $2^\delta = 2$ of Eq. (7) that has not been previously appreciated in similar estimations. They also take into account that only $\ell$ even or odd modes survive, depending on the statistics.

- For non-self-conjugate DM, the above takes into account only particle-antiparticle interactions, as is standard. However, DM may be depleted also via particle-particle inelastic scatterings (see e.g. [15] for a model), which increases the bounds on the DM mass shown in the table by a factor $\sqrt{3}$.

- We reiterate that there is no model-independent bound on the mass of thermal-relic DM due to unitarity, since which and how many partial waves contribute significantly depends on the model.

We close this discussion with an interpretation of Eqs. (60) in hydrodynamic terms. The resummation of inelastic interactions accounts for the reduction in the particle flux due to inelastic scatterings, as suggested by the continuity equation [23]. This reduction suppresses both elastic and inelastic scattering, by the factor $(1 + y_{\ell, \text{unreg}})^{-2}$.

For inelastic scatterings, this effect competes against the fact that the probability for inelastic scattering increases with $y_{\ell, \text{unreg}}$. For $y_{\ell, \text{unreg}} \le 1$, the increase of the inelastic scattering probability with $y_{\ell, \text{unreg}}$ dominates, while for $y_{\ell, \text{unreg}} > 1$, the suppression of the particle flux takes over.

Elastic scattering arises from two contributions. Particles may scatter due to the purely elastic couplings with probability $x_{\ell, \text{unreg}}$. They evade this scattering with probability $1 - x_{\ell, \text{unreg}}$, but may scatter inelastically, with the products of the inelastic processes scattering back into the original species via the inverse reactions; this introduces an additional probability factor $y_{\ell, \text{unreg}}^2$. The regenerated flux dominates elastic scattering at sufficiently large $y_{\ell, \text{unreg}}$ (except if $x_{\ell, \text{unreg}} = 1$), where the regeneration probability $\propto y_{\ell, \text{unreg}}^2$ exactly balances out the suppression $(1 + y_{\ell, \text{unreg}})^{-2}$, resulting in $x_{\ell, \text{reg}} \approx 1 - x_{\ell, \text{unreg}}$.

## 9. Conclusions

The violation of inelastic unitarity bounds in the non-relativistic regime by state-of-the art computations hampers a variety of phenomenological investigations in the frontiers of DM research. Resolving this issue is important for improving our theoretical understanding, as well as interpreting and guiding experimental searches. We have derived a simple analytical regularization scheme, given by Eqs. (60), whose input are the inelastic cross-sections as affected only by the real part of the potential, and output are the unitarized elastic and inelastic cross-sections. The scheme applies to any partial wave, and is model-independent as it makes no assumptions about the momentum dependence of the unregulated inelastic amplitudes, except for analyticity and convergence requirements. (The more general Eqs. (37) and (38) along

with Eq. (51) can be employed when these requirements do not hold.) These results can be, and must be, employed in investigations of new physics scenarios in the non-relativistic regime, with wide-ranging implications.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Data availability

No data was used for the research described in the article.

## Acknowledgements

We thank Mark Goodsell, Miguel Paulos, Yifei He and Tobias Binder for useful discussions. This work was supported by the European Union's Horizon 2020 research and innovation programme under grant agreement No 101002846, ERC CoG CosmoChart.

## References

[1] S. Weinberg, Physical processes in a convergent theory of the weak and electromagnetic interactions, Phys. Rev. Lett. 27 (1971) 1688-1691, https://doi.org/10.1103/PhysRevLett.27.1688.

[2] J.M. Cornwall, D.N. Levin, G. Tikopoulos, Derivation of gauge invariance from high-energy unitarity bounds on the s matrix, Phys. Rev. D 10 (1974) 1145, https://doi.org/10.1103/PhysRevD.10.1145, Erratum: Phys. Rev. D 11 (1975) 972.

[3] B.W. Lee, C. Quigg, H.B. Thacker, The strength of weak interactions at very high-energies and the Higgs boson mass, Phys. Rev. Lett. 38 (1977) 883-885, https://doi.org/10.1103/PhysRevLett.38.883.

[4] I. Baldes, K. Petraki, Asymmetric thermal-relic dark matter: Sommerfeld-enhanced freeze-out, annihilation signals and unitarity bounds, J. Cosmol. Astropart. Phys. 1709 (09) (2017) 028, https://doi.org/10.1088/1475-7516/2017/09/028, arXiv: 1703.00478.

[5] K.L. Kowalski, Some consequences of off-shell unitarity, Phys. Rev. 144 (1966) 1239-1240, https://doi.org/10.1103/PhysRev.144.1239.

[6] L. Hui, Unitarity bounds and the cuspy halo problem, Phys. Rev. Lett. 86 (2001) 3467-3470, https://doi.org/10.1103/PhysRevLett.86.3467, arXiv:astro-ph/0102349.

[7] K. Griest, M. Kamionkowski, Unitarity limits on the mass and radius of dark matter particles, Phys. Rev. Lett. 64 (1990) 615, https://doi.org/10.1103/PhysRevLett.64.615.

[8] B. von Harling, K. Petraki, Bound-state formation for thermal relic dark matter and unitarity, J. Cosmol. Astropart. Phys. 12 (2014) 033, https://doi.org/10.1088/1475-7516/2014/12/033, arXiv:1407.7874.

[9] J. Hisano, S. Matsumoto, M.M. Nojiri, Unitarity and higher order corrections in neutralino dark matter annihilation into two photons, Phys. Rev. D 67 (2003) 075014, https://doi.org/10.1103/PhysRevD.67.075014, arXiv:hep-ph/0212022.

[10] A.D. Sakharov, Interaction of an electron and positron in pair production, Zh. Eksp. Teor. Fiz. 18 (1948) 631-635, https://doi.org/10.1070/PU1991v034n05ABEH002492, Usp. Fiz. Nauk 161 (5) (1991) 29.

[11] A. Sommerfeld, Über die Beugung und Bremsung der Elektronen, Ann. Phys. 403 (3) (1931) 257-330.

[12] K. Petraki, M. Postma, M. Wiechers, Dark-matter bound states from Feynman diagrams, J. High Energy Phys. 1506 (2015) 128, https://doi.org/10.1007/JHEP06(2015)128, arXiv:1505.00109.

[13] K. Petraki, M. Postma, J. de Vries, Radiative bound-state-formation cross-sections for dark matter interacting via a Yukawa potential, J. High Energy Phys. 04 (2017) 077, https://doi.org/10.1007/JHEP04(2017)077, arXiv:1611.01394.

[14] S. El Hedri, A. Kaminska, M. de Vries, A Sommerfeld toolbox for colored dark sectors, Eur. Phys. J. C 77 (9) (2017) 622, https://doi.org/10.1140/epjc/s10052-017-5168-z, arXiv:1612.02825.

[15] R. Oncala, K. Petraki, Dark matter bound state formation via emission of a charged scalar, J. High Energy Phys. 02 (2020) 036, https://doi.org/10.1007/JHEP02(2020)036, arXiv:1911.02605.

[16] T. Binder, M. Garny, J. Heisig, S. Lederer, K. Urban, Excited bound states and their role in dark matter production, Phys. Rev. D 108 (9) (2023) 095030, https://doi.org/10.1103/PhysRevD.108.095030, arXiv:2308.01336.

[17] R. Oncala, K. Petraki, Bound states of WIMP dark matter in Higgs-portal models. Part II. Thermal decoupling, J. High Energy Phys. 08 (2021) 069, https://doi.org/10.1007/JHEP08(2021)069, arXiv:2101.08667.
[18] R. Oncala, K. Petraki, Bound states of WIMP dark matter in Higgs-portal models. Part I. Cross-sections and transition rates, J. High Energy Phys. 06 (2021) 124, https://doi.org/10.1007/JHEP06(2021)124, arXiv:2101.08666.

[19] P. Ko, T. Matsui, Y.-L. Tang, Dark matter bound state formation in fermionic $Z_3$ DM model with light dark photon and dark Higgs boson, J. High Energy Phys. 10 (2020) 082, https://doi.org/10.1007/JHEP10(2020)082, arXiv:1910.04311.

[20] J. Harz, K. Petraki, Radiative bound-state formation in unbroken perturbative non-Abelian theories and implications for dark matter, J. High Energy Phys. 07 (2018) 096, https://doi.org/10.1007/JHEP07(2018)096, arXiv:1805.01200.

[21] S. Cassel, Sommerfeld factor for arbitrary partial wave processes, J. Phys. G 37 (2010) 105009, https://doi.org/10.1088/0954-3899/37/10/105009, arXiv:0903.5307.

[22] V.B. Berestetskii, E.M. Lifshitz, V.B. Pitaevskii, Relativistic Quantum Theory. Pt. 1, 1971.

[23] K. Blum, R. Sato, T.R. Slatyer, Self-consistent calculation of the Sommerfeld enhancement, J. Cosmol. Astropart. Phys. 1606 (06) (2016) 021, https://doi.org/10.1088/1475-7516/2016/06/021, arXiv:1603.01383.

[24] Y. Yamaguchi, Two nucleon problem when the potential is nonlocal but separable. I, Phys. Rev. 95 (1954) 1628-1634, https://doi.org/10.1103/PhysRev.95.1628.

[25] E. Braaten, E. Johnson, H. Zhang, Zero-range effective field theory for resonant wino dark matter. Part III. Annihilation effects, J. High Energy Phys. 05 (2018) 062, https://doi.org/10.1007/JHEP05(2018)062, arXiv:1712.07142.

[26] U. Aydemir, M.M. Anber, J.F. Donoghue, Self-healing of unitarity in effective field theories and the onset of new physics, Phys. Rev. D 86 (2012) 014025, https://doi.org/10.1103/PhysRevD.86.014025, arXiv:1203.5153.

[27] A. Kamada, S. Kobayashi, T. Kuwahara, Perturbative unitarity of strongly interacting massive particle models, J. High Energy Phys. 02 (2023) 217, https://doi.org/10.1007/JHEP02(2023)217, arXiv:2210.01393.

[28] A. Messiah, Quantum Mechanics, North-Holland Pub. Co., 1962.

[29] C.J. Joachain, Quantum Collision Theory, North Holland, 1975.

[30] M. McMillan, Separable nonlocal potentials and Regge poles, Nuovo Cimento 29 (1963) 1043-1050.

[31] M. Shah-Jahan, Comments on separable potentials and Coulomb interactions, Phys. Rev. D 8 (1973) 2744-2746, https://doi.org/10.1103/PhysRevD.8.2744.

[32] H.F. Arellano, G. Blanchon, Exact scattering waves off nonlocal potentials under Coulomb interaction within Schrödinger's integro-differential equation, Phys. Lett. B 789 (2019) 256-261, https://doi.org/10.1016/j.physletb.2018.12.004, arXiv:1806.05006.

[33] R.A. Mapleton, Characteristic Coulomb Green's function and its eigenfunction expansion, J. Math. Phys. 2 (4) (1961) 478-482, https://doi.org/10.1063/1.1703733.

[34] S. Tulin, H.-B. Yu, Dark matter self-interactions and small scale structure, Phys. Rep. 730 (2018) 1-57, https://doi.org/10.1016/j.physrep.2017.11.004, arXiv:1705.02358.

[35] J.S. Bullock, M. Boylan-Kolchin, Small-scale challenges to the ΛCDM paradigm, Annu. Rev. Astron. Astrophys. 55 (2017) 343-387, https://doi.org/10.1146/annurev-astro-091916-055313, arXiv:1707.04256.