A Practical and Consistent Parametrization of Dark Matter Self-Interactions

TL;DR

This paper introduces a simple, model-independent parametrization for dark matter self-interactions that accurately captures their velocity dependence across different theoretical scenarios, aiding astrophysical analysis.

Contribution

It presents a two-parameter effective-range formalism to approximate dark matter self-interaction cross sections, applicable to various models including Yukawa, Breit-Wigner, and bound states.

Findings

The parametrization accurately models velocity-dependent cross sections.

It provides a unified interpretation across multiple dark matter interaction models.

The approach can be refined for specific regimes where it initially performs poorly.

Abstract

Self-interacting dark matter has been proposed to explain the apparent mass deficit in astrophysical small-scale halos, while observations from galaxy clusters suggest that the corresponding cross section depends on the velocity. Accounting for this is often believed to be highly model-dependent with studies mostly focusing on scenarios with light mediators. Based on the effective-range formalism, in this work we point out a model-independent approach which accurately approximates the velocity dependence of the self-interaction cross section with only two parameters. We illustrate how this parameterization can be simultaneously interpreted in various well-motivated scenarios, including self-interactions induced by Yukawa forces, Breit-Wigner resonances and bound states. We investigate the astrophysical implications and discuss how the approximation can be improved in certain special regimes where it works poorly.

Equations

$$f(k,\theta )$$

$$f_{\ell }(k)\equiv \frac{e^{2i{\delta }_{\ell }(k)}-1}{2ik}=\frac{1}{k(\cot {\delta }_{\ell }(k)-i)}.$$

$$k^{2\ell +1}\cot {\delta }_{\ell }(k)\simeq -\frac{1}{a_{\ell }^{2\ell +1}}+\frac{1}{2r_{e,\ell }^{2\ell -1}}{k}^2.$$

$${\sigma }_0=\frac{4\pi }{k^2}{\sin }^2{\delta }_0\approx \frac{4\pi a^2}{1+k^2\left(a^2-a{r}_e\right)+\frac{1}{4}{a}^2{r}_e^2{k}^4}.$$

$$V(r)=\pm \alpha \frac{e^{-m_{\varphi }r}}{r}.$$

$${\sigma }_{\ell }=\frac{4\pi (2l+1)}{mE}\frac{\Gamma (E{)}^2/4}{{\left(E-E_R\right)}^2+{\Gamma }^2(E)/4}.$$

$$a_{\ell }=-\frac{\Gamma (E_R{)}^{\frac{1}{2\ell +1}}}{2^{\frac{1}{2\ell +1}}{E}_R^{\frac{2\ell +3}{4\ell +2}}{m}^{\frac{1}{2}}},\text{leavevmode}\text{leavevmode}{r}_{e,\ell }=-\frac{2^{\frac{2}{-2\ell +1}}{E}_R^{\frac{2\ell +1}{-4\ell +2}}}{{\Gamma }^{\frac{1}{-2\ell +1}}{m}^{\frac{1}{2}}}.$$

$$\sigma (v)=4\pi a^2{\left({\left(1-\frac{1}{8}\frac{r_e}{a}(mav{)}^2\right)}^2+\frac{1}{4}(mav{)}^2\right)}^{-1},$$

$$\begin{split}1\,\mathrm{cm^{2}/g}\lesssim\frac{\sigma}{m}\Bigg{|}_{v=10\,\mathrm{km/s}}&\lesssim 10\,\mathrm{cm^{2}/g}\,,\\ \frac{\sigma}{m}\Bigg{|}_{v=2000\,\mathrm{km/s}}&\lesssim 1\,\mathrm{cm^{2}/g}\,.\end{split}$$

$$f(v,v_0)=\frac{4v^2{e}^{-v^2/v_0^2}}{\sqrt{\pi }{v}_0^3}\Theta \left(v_{\text{max}}-v\right),$$

$$V=\frac{1}{2}{m}^2{\varphi }^2+\frac{1}{4!}\lambda {\varphi }^4$$

$${\sigma }_0=\frac{{\lambda }^2}{128\pi m^2}.$$

$${\sigma }_0=4\pi a^2.$$

$$m\lesssim 8.14\mathrm{MeV}\text{leavevmode}{\lambda }^{2/3}{\left(\frac{1\text{leavevmode}{\mathrm{cm}}^2/\mathrm{g}}{{\sigma }_0/m}\right)}^{1/3}.$$

$$\tan {\delta }_0$$

$$a=-\frac{m\alpha }{m_{\varphi }^2},\text{leavevmodenobreak leavevmodenobreak andleavevmodenobreak leavevmodenobreak leavevmodenobreak }{r}_e=\frac{4}{m\alpha }.$$

$$k_{\pm }^{\text{pole}}=\frac{i}{a}\frac{2}{1\pm \sqrt{1-2r_e/a}}.$$

$$2m_{\star }ϵ={\left(\frac{1}{a}+m_{\star }ϵr_e\right)}^2.$$

$${\delta }_0$$

$${\sigma }_0$$

$$e^{2i{\delta }_0(k)}=e^{2ikR}{e}^{2i{\delta }_a(k)}=e^{2ikR}\left(\frac{-\frac{1}{\mathfrak{a}}+\frac{1}{2}{\mathfrak{r}}_{\mathfrak{e}}{k}^2+ik}{-\frac{1}{\mathfrak{a}}+\frac{1}{2}{\mathfrak{r}}_{\mathfrak{e}}{k}^2-ik}\right),$$

$$f_0(k)=\frac{e^{2i{\delta }_0(k)}-1}{2ik}=\frac{e^{2ikR}-1}{2ik}+\frac{e^{2ikR}}{-\frac{1}{\mathfrak{a}}+\frac{1}{2}{\mathfrak{r}}_{\mathfrak{e}}{k}^2-ik},$$

$$V(r)=-m{\alpha }^2\Theta ({\lambda }_{\text{cut}}-r)+\Theta (r-{\lambda }_{\text{cut}})\frac{\alpha }{r}{e}^{-m_{\varphi }r},$$

$$k\cot {\delta }_0\simeq k\cot (\text{Re}{\delta }_0)-\frac{ik\text{Im}{\delta }_0}{{\sin }^2(\text{Re}{\delta }_0)}\simeq -\frac{1}{a}+\frac{r_e}{2}{k}^2,$$

$${\sigma }_{an,0}(k)=\frac{4\pi }{k^2}\frac{1-|e^{2i{\delta }_0}{|}^2}{4}\simeq \frac{{\sigma }_0(k)}{k}\frac{|\text{Im}a|}{(\text{Re}a{)}^2}.$$

$$\frac{1}{r^{2}}\frac{d}{dr}\Big{(}r^{2}\frac{dR_{\ell,k}}{dr}\Big{)}+\Big{(}k^{2}-\frac{\ell(\ell+1)}{r^{2}}-mV(r)\Big{)}R_{\ell,k}=0\,,$$

$$R_{\ell ,k}(r)\propto \cos {\delta }_{\ell }{j}_{\ell }(kr)-\sin {\delta }_{\ell }{n}_{\ell }(kr)\approx \frac{1}{r}\sin \left(kr-\frac{l\pi }{2}+{\delta }_{\ell }\right),$$

$$t_{\ell ,k}(r)=\frac{j_{\ell }(kr)\left(\frac{R_{\ell ,k}^{\prime }(r)}{R_{\ell ,k}(r)}-\frac{\ell }{r}\right)+k{j}_{\ell +1}(kr)}{n_{\ell }(kr)\left(\frac{R_{\ell ,k}^{\prime }(r)}{R_{\ell ,k}(r)}-\frac{\ell }{r}\right)+k{n}_{\ell +1}(kr)}.$$

$$\frac{d{t}_{\ell ,k}(r)}{dr}=-km{r}^{2}V(r){\left(j_{\ell }(kr)-t_{\ell ,k}(r)n_{\ell }(kr)\right)}^2.$$

$$t_{\ell ,k}(0)=0$$

$$a_{\ell }^{2\ell +1}\equiv -\underset{k\to 0}{\lim }\frac{\tan {\delta }_{\ell }}{k^{2\ell +1}}.$$

Full text

††preprint: DESY 19-137††preprint: IPMU19-0113

1. I Introduction

2. II The effective-range formalism

3. III Astrophysical implications

4. III.1 Velocity dependence in effective range theories

5. III.2 Realistic velocity distributions

6. IV Interpreting $a$ and $r_{e}$ in terms of model parameters

7. IV.1 Contact interaction

8. IV.2 SIDM with a light mediator

9. IV.3 SIDM via bound states or virtual levels

10. IV.4 Resonant SIDM

11. IV.5 SIMPs

12. V Improving the effective-range approximation

13. V.1 Antiresonances

14. V.2 Sharp resonances

15. V.3 Inelastic scatterings

16. VI Summary and Outlook

17. A The effective range theory

18. A.1 A simple method to find the phase shift

19. A.2 The $S$-wave case

20. A.3 The Hulthén Potential

21. B Velocity-averaged cross sections

I Introduction

In contrast to the ordinary substances found on Earth, more than three quarters of the matter in the Universe is not made of protons, neutrons, or electrons. This is the so-called dark matter (DM) and identifying its particle nature is one of the chief goals of particle physics and cosmology today. According to the $\Lambda$CDM model Aghanim et al. (2018), which accurately describes the Universe at cosmological scales, DM interacts very weakly with normal matter and it was cold and collisionless during the formation of structures in the early universe. Although DM can be treated as collisionless particles at large scales, non-gravitational DM scatterings can still occur in the dense central regions of small-scale halos such as those of dwarf or low-surface-brightness galaxies. This is the self-interacting dark matter (SIDM) hypothesis, which was proposed Spergel and Steinhardt (2000) to explain the seeming discrepancies between observations of the smallest DM halos that we can currently observe and certain predictions of the $\Lambda$CDM model; see Tulin and Yu (2018); Bullock and Boylan-Kolchin (2017) for recent reviews.

The aforementioned discrepancies can be explained if DM elastically scatters with a cross section per unit of mass as large as several cm${}^{2}/$g when it moves at approximately $10$ km/s, i.e., roughly the DM velocity dispersion in small-scale objects Dave et al. (2001); Vogelsberger et al. (2012); Rocha et al. (2013); Peter et al. (2013); Elbert et al. (2015); Fry et al. (2015). Meanwhile, recent studies on halo dynamics at cluster scales provide upper bounds on the self-interaction cross section of around $0.2\,\mathrm{-}1\,\mathrm{cm^{2}/g}$ Randall et al. (2008); Kaplinghat et al. (2016); Robertson et al. (2017); Bondarenko et al. (2018); Elbert et al. (2018); Harvey et al. (2018), which are associated with typical DM velocities of the order of $10^{3}\,\mathrm{km/s}$.

A natural question then arises: how can these observations be interpreted in terms of the properties of the DM particle? One possibility is to postulate a specific DM model, and translate the previous velocity-dependent cross section in terms of masses and couplings. For instance, this has been done for scenarios where DM interacts by means of a light mediator Feng et al. (2010); Buckley and Fox (2010); Tulin et al. (2013a) and for models in which DM resonantly self-scatters Chu et al. (2019); Ibe and Yu (2010); Duch and Grzadkowski (2017). Nevertheless, from the phenomenological point of view, this is not very practical, not only because the cross sections typically depend in a complicated way on the model parameters, but also because there are a myriad of SIDM scenarios. The essence of this work is to propose a simple parametrization of the DM self-interaction cross section, which approximates with great accuracy the velocity-dependent effects, and most importantly, which interpolates the predictions of different DM scenarios allowing to establish comparisons among them. More precisely, here we advocate the use of the effective-range theory as a model-independent way to study the velocity dependence of SIDM. Notice that this approach has been adopted in concrete models of DM before March-Russell and West (2009); Braaten and Hammer (2013); Cline et al. (2014); Braaten et al. (2018); Mahbubani et al. (2019).

The effective-range approach was formulated Bethe (1949); Blatt and Jackson (1949) as an effort to explain the non-relativistic scattering of neutrons by protons. Based on simple assumptions from quantum mechanics, this approach suggests that the scattering observables can be parametrized in terms of two quantities: the scattering length $a$, and the effective range $r_{e}$. While being very predictive, the effective-range theory does not demand a precise knowledge of the underlying interactions among the colliding particles, apart from the requirement that the scattering force must vanish at sufficiently large distances. In fact, due to this, it can describe the non-relativistic scattering induced by contact interactions, light mediators, and Breit-Wigner resonances, among others.

A brief historical remark may be helpful to readers.111Hans Bethe personally described the history in a YouTube video https://www.youtube.com/watch?v=hbcQMG2XpTI. In late 1940s people proposed different models to describe nucleon-nucleon scattering cross sections at low energies. After many explicit calculations, it became clear that only two parameters are relevant, independent of details of models. Schwinger came up with a proof why that was the case in an unpublished lecture note. Blatt and Jackson Blatt and Jackson (1949) showed with more calculations that indeed only two parameters were necessary to explain the data. Then Bethe Bethe (1949) came up with a simple and elegant proof (reproduced in Appendix A) to understand this observation. The flip side of this remarkable simplicity is that we gain very little information on the detailed model from the data. Bethe wrote “practically no information could be obtained, from classical scattering experiments, on the shape of the potential.” On the other hand, for the purpose of describing the impact of self-interaction among dark matter particles in various halos, this is a boon; we need to specify only two parameters (and the mass of dark matter) in order to simulate the impact of self-interactions without the need for dealing with explicit models. This is why we propose the use of the effective range theory for the study of SIDM.

This paper is organized as follows. In Section II, we review the effective-range approach, and explain how it can be useful to describe DM self-scattering. In Section III, we discuss the implications of this approach for describing the velocity-dependence of DM self-interactions in astrophysical halos in a model-independent way. In Section IV, we investigate how to relate the scattering length and effective range to parameters of concrete DM models. In Section V, we propose a concrete method on how the effective-range approach can be extended or improved in some cases where it fails. Section VI provides the final conclusions and future prospects.

II The effective-range formalism

Before focusing on SIDM, we will first introduce the effective-range approximation. In any collision the differential cross section $d\sigma/d\Omega$ determines the scattering rate and describes the velocity dependence of the process. Up to a possible symmetry factor, it is given by $d\sigma=|f(k,\theta)|^{2}d\Omega$, where $f(k,\theta)$ is the scattering amplitude with $k$ being the incoming momentum.222Here and below, we separate the center-of-mass motion and hence $k$ is the relative momentum of two particles that scatter. The orbital angular momentum $\ell$ below is defined in the center-of-mass frame. Note also that we are dealing with low velocities, less than about $10^{-2}c$, which are typical in halos, so the use of non-relativistic quantum mechanics is justified. For collisions with definite orbital angular momentum, $\ell$, the amplitude is proportional to the Legendre polynomial, $P_{\ell}(\cos\theta)$, with the corresponding coefficient defining the partial-wave amplitude, $f_{\ell}(k)$. More precisely

[TABLE]

The second relation defines the phase shift, $\delta_{\ell}(k)$ for the $\ell$ partial wave. While the precise value of $\delta_{\ell}(k)$ must be obtained by solving the Schrödinger equation describing the scattering process, this phase shift always satisfies some general requirements. For instance, it must be real if inelastic processes are absent. In this work, we are concerned with elastic scatterings in astrophysical halos, as a result, unless stated otherwise, we will assume that inelastic processes are relatively weaker and take $\delta_{\ell}$ real (see Section V.3 for how to include the inelastic processes).

Another requirement on the phase shift is that, for finite-range interactions, the function $k^{2\ell+1}\cot\delta_{\ell}(k)$ must be analytic at $k=0$ (for more details see Appendix A). The effective-range approximation Bethe (1949); Blatt and Jackson (1949) consists in neglecting the high-order terms in the corresponding expansion in $k^{2}$, so that

[TABLE]

The quantities $a_{\ell}$ and $r_{e,\ell}$ thus defined are known as the scattering length and the effective range, respectively. This approximation describes the phase shift with good accuracy at sufficiently low energies. Consequently, if one partial wave dominates the scattering process, the velocity dependence of the cross section is determined by only two parameters.

Let us focus on the $S$-wave case, which dominates the low-energy scattering rate in many situations of interest. In this case, the cross section is given by 333From now on, for simplicity we will omit the subscript $\ell=0$ for the scattering length and the effective range in the $S$-wave.

[TABLE]

Note that the unitarity bound, $4\pi/k^{2}$, is saturated for $|a|\to\infty$. One important example of this kind is the case of particles of mass $m$ interacting via the Yukawa potential

[TABLE]

Fig. 1 compares the numerically-evaluated cross section and the approximation based on Eq. (4) for a particular region of the parameter space.444Note that it is essential to go beyond the lowest-order perturbation theory (first Born approximation) $\sigma=4\pi(\alpha m/m_{\phi}^{2})^{2}$, which does not depend on the sign of $\alpha$ nor produces spikes in the full calculation. Effective range theory reproduces both correctly. Despite its simplicity, the effective range approximation works very well, i.e. high-order terms of $k^{2}$ in Eq. (3) can be neglected. In particular, it is able to reproduce the peak structure of the cross section. As is well known, such peaks are related to zero-energy bound states induced by the attractive potential.

A similar example is the non-relativistic scattering of two nucleons. In the case of proton-neutron collision, for the spin-one channel and kinetic energies up to a few MeV, Eq. (4) accurately describes the velocity dependence of the corresponding cross section with $a=5.42\,\mathrm{fm}$ and $r_{e}=1.75\,\mathrm{fm}$ (see e.g. Noyes (1972); Stoks et al. (1993)). This is also related to a bound state: the deuteron.

Likewise, the collision of two neutrons –for which the total spin is zero– can be characterized by $a=-18.9\,\mathrm{fm}$ and $r_{e}=2.75\,\mathrm{fm}$. In contrast, in this case no real bound state exists. Instead, the scattering is induced by a virtual level,555A bound state is a pole in the scattering amplitude $f_{\ell}(k)$ along the positive imaginary axis on the complex $k$ plane, which corresponds to exponential damping of the radial wave function $\propto e^{ikr}$. When parameters of the potential are varied, a pole may move to the negative imaginary axis, which no longer describes a bound state because the wave function grows exponentially. However, the existence of a pole in the scattering amplitude can produce a pronounced enhancement in the cross section and is hence important. In this case, the pole is called a virtual level. It should not be confused with virtual particles or virtual states that refer to intermediate particles or states (propagators) in perturbation theory. commonly known as the dineutron (see Fig. 2).

Equally interesting is the fact that Breit-Wigner resonances can also be described using the effective-range approximation. For simplicity let us suppose that the colliding particles are scalars with the same mass $m$, and that the resonance has spin $\ell$. Thus, if the energy, $E=k^{2}/m$, is sufficiently close to the resonance $E_{R}$, the cross section is dominated by the partial wave $\ell$ so that

[TABLE]

The width in general varies with the energy in such way that $\Gamma(E)\propto E^{(2\ell+1)/2}$ (see e.g. Chu et al. (2019)). Using the effective-range approximation to the phase shift (Eq. (3)), we find that the cross section, $\sigma_{\ell}=4\pi(2\ell+1)\sin^{2}\delta_{\ell}/k^{2}$, exactly matches the previous formula with

[TABLE]

Far from the resonance, some deviations are expected. In fact, as shown in right panel of Fig. 1, the effective-range approximation also fails close to the antiresonances, i.e., where the cross section vanishes. We will elaborate more on these cases in Section V. Likewise, when the range of the Yukawa potential, $m_{\phi}^{-1}$, is close to or larger than the de Broglie wavelength of the incoming particles, $k^{-1}$, the approximation fails. This region, usually referred to as classical regime, corresponds to $m_{\phi}\lesssim 5\,\mathrm{MeV}$ for the parameter region of the right panel of Fig. 1. In fact, this is true for any potential, for momenta larger than the inverse of the force range, not only must one include higher-order terms in Eq. (3), but also the differential cross section receives contributions from high partial waves. In this case, the exact values of more phase shifts $\delta_{\ell}$ are needed to obtain the total scattering cross section. As a result, the effective-range approach can not be applied for long-range forces. For more details, see e.g. F. O’Malley et al. (1961); Berger and Spruch (1965); van Haeringen and Kok (1982).

In summary, the effective-range approximation properly describes many types of low-energy scattering. These cases differ in the magnitude and sign of the effective range parameters, as will be explained in Section IV and as sketched in Fig. 2.

III Astrophysical implications

The main hypothesis of SIDM paradigm is that small-scale DM halos such as those of dwarf galaxies do not develop a high central density because its DM particles self-scatter with a cross section per unit mass in the range $1-10\,\mathrm{cm^{2}/g}$ Vogelsberger et al. (2012); Rocha et al. (2013); Peter et al. (2013). On the other hand, observations of clusters of galaxies indicate that $\sigma/m\lesssim 0.2\,\mathrm{-}1\,\mathrm{cm^{2}/g}$ Robertson et al. (2017); Randall et al. (2008); Kaplinghat et al. (2016); Bondarenko et al. (2018); Elbert et al. (2018); Harvey et al. (2018). Since in the former objects the average DM relative velocity is typically of the order of $10$ km/s, whereas in clusters of galaxies it is around $2000$ km/s, a velocity-dependent cross section is required in order to accommodate both.

Before discussing this in detail, let us note that we use $m$ for the DM mass, $v$ for the relative velocity between two initial DM particles in the centre-of-mass (CM) frame, $m_{\star}=m/2$ for the reduced mass and $k=m_{\star}v=mv/2$ for each incoming DM momentum.

III.1 Velocity dependence in effective range theories

In the effective-range framework, the velocity-dependent cross section is given by 666 The transfer cross section, $\sigma_{T}=\int d\sigma(1-\cos\theta)$, is typically used as a proxy for the scattering effects in DM halos. This is because, on the one hand, $\sigma_{T}$ takes into account that perpendicular scattering is most efficient for thermalizing the DM halo and affecting structure observables. On the other hand, SIDM studies often discuss scatterings induced by the exchange of a light mediator, which exhibits a divergence in the forward direction regularized by the transfer cross section. For $S$-wave scattering, $\sigma_{T}$ and $\sigma_{0}$ coincide and are therefore interchangeable.

[TABLE]

where the signs of the scattering length and the effective range only enter in the equation via their ratio, and can not be separately constrained by studying the velocity dependence of the DM scattering.

The velocity-dependence of the scattering cross section is also shown in Fig. 3. At very low velocities the cross section is roughly constant and equal to $4\pi a^{2}$. If $r_{e}/a<1$, the cross section monotonously decreases, most appreciably for high DM velocities, $v\gtrsim(m|a|)^{-1}$. In contrast, if $r_{e}/a>1$, the cross section increases with $v$ until it reaches the maximum $4\pi r_{e}^{2}a/(2r_{e}-a)$ at $v_{\text{peak}}=2(m|r_{e}|)^{-1}\sqrt{2(r_{e}/a-1)}$ and then decreases. If $|a|\ll 1/m$, the corresponding cross section can be considered as a constant in all realistic DM halos.

To numerically illustrate this for halos of various sizes, Fig. 4 shows the contours of the self-interaction cross section per unit mass at $v=2000\,\mathrm{km/s}$ for $\sigma/m|_{v\to 0}$ equal to $1\,\mathrm{cm^{2}/g}$ (left) and $10\,\mathrm{cm^{2}/g}$ (right). We would like to note that $\sigma/m|_{v\to 0}$ approximates the corresponding values in dwarf scales at 1% level in the parameter space shown in the figure. From Fig. 4, we conclude that GeV SIDM is associated with scattering lengths of several fm and that sub-GeV SIDM is either excluded by cluster observations or requires a cross section of around $1\,\mathrm{cm^{2}/g}$ throughout all scales of interest. While some of these conclusions have been obtained in specific SIDM scenarios such as those involving a light mediator Tulin et al. (2013b) or resonant SIDM Chu et al. (2019), we would like to emphasize that these conclusions apply to any model where the effective-range approach applies.

Fig. 4 also suggests that the ratio $r_{e}/a$ is poorly constrained by the velocity dependence of the cross section. In the light of this and in order to consider a wider range of parameter space, in Fig. 5 we show three possibilities for the effective range compared to the scattering length. Concretely, for each fixed value of $r_{e}/a$, we illustrate the parameter space simultaneously satisfying

[TABLE]

The former takes place within the colored region while the latter, as in Fig. 4, corresponds to the region not excluded by the gray area. Recent studies have claimed stronger constraints of $0.2\,\mathrm{cm^{2}/g}$ from observations of galaxy clusters Bondarenko et al. (2018); Elbert et al. (2018); Harvey et al. (2018), which are indicated as a gray dashed line in each panel.

In each panel there are regions where all the constraints are simultaneously satisfied. This confirms our previous remark that the ratio $r_{e}/a$ is largely unconstrained. Likewise, scenarios with a cross section of $1\,\mathrm{cm^{2}/g}$ at dwarf and cluster scales are those for which the borders of the gray and the colored regions lie on top of each other and –as mentioned above– they correspond to DM masses below a few GeV.

III.2 Realistic velocity distributions

So far we have assumed a monochromatic velocity distribution for all DM particles in each halo. Below, we take into account the realistic distribution of DM velocities and then consider the corresponding average cross section for individual DM halos. The former is typically achieved by assuming a Maxwell-Boltzmann distribution with a cut-off scale:

[TABLE]

where $v_{\text{max}}$ is the escape velocity and $v_{0}$ is a parameter determining the typical velocities in the DM halo. For $v_{\text{max}}\gg v_{0}$, the average velocity is $\langle v\rangle=2v_{0}/\sqrt{\pi}$. More concretely, we take the average cross section as $\langle\sigma v\rangle/\langle v\rangle=({\int^{\infty}_{0}f(v,v_{0})\,\sigma v\,dv})/({\int^{\infty}_{0}f(v,v_{0})\,v\,dv})$. Using a Maxwell-Boltzmann distribution would at most modify the curves of Fig. 5 mildly, so the conclusion of the previous subsection remains unchanged. For a detailed calculation and a numerical comparison of $\langle\sigma v\rangle/\langle v\rangle$ and $\sigma(\langle v\rangle)$, see Appendix. B.

On the observational side, extracting cross sections from experimental data is challenging and generally requires delicate N-body simulations at present. An intermediate method is given by the semi-analytical method proposed in Kaplinghat et al. (2016), which allows to infer the velocity-averaged cross section per unit mass, $\langle\sigma v\rangle/m$, for a given DM halo. This method was applied to five clusters from Newman et al. (2013), seven low-surface-brightness (LSB) spiral galaxies in Kuzio de Naray et al. (2008) and six dwarf galaxies of the THINGS sample Oh et al. (2011) (see also Valli and Yu (2018)). Fig. 6 shows their resulting values in green, blue and red, respectively.

Using the velocity-averaged scattering cross section, we can fit $|a|$, $r_{e}/a$ and $m$ to the aforementioned semi-analytical results and thus constrain the parameters of the effective-range theory. The best-fit point is shown in Fig. 6 and corresponds to the benchmark S1.777As mentioned before, the sign of $a$ or $r_{e}$ can not be fit by studying the velocity dependence of the scattering. Moreover, since S1 has very small value of $r_{e}$, flipping the sign of either of them, in practice, gives the same cross section as a function of the velocity. As expected, it fulfills the condition stated in Eqs. (9).

We would like to emphasize that the points shown in Fig. 6 should be taken with caution, as subtle effects, such as tidal stripping, still need to be further studied to understand the (sub-)halo dynamics (see e.g. Sokolenko et al. (2018); Kummer et al. (2019); Kaplinghat et al. (2019); Sameie et al. (2019); Kahlhoefer et al. (2019) for recent discussions). In fact, the need of a sizeable DM self-interaction at cluster scales is under debate as this relies on the assumption that one can robustly infer the existence of cores in clusters of galaxies. This motivates us to also consider other possibilities which are not necessarily fitting the green points but in agreement with conservative bounds at cluster scales, $\sigma/m\lesssim 1\,\mathrm{cm^{2}/g}$. These are the benchmarks S2, S3 and S4 labeled in Fig. 6.

The benchmark S2 fits the dwarf and LSB data points fairly well with a relatively high cross section but is too low to fully accommodate the cluster points. This is because its peak is very pronounced with $r_{e}/a\gg 1$ and $\langle\sigma v\rangle$ decreases very rapidly for $\langle v\rangle$ greater than the peak velocity. This can be achieved with a narrow resonance Chu et al. (2019). In contrast, another benchmark S3 describes an almost constant self-interaction cross section. Interestingly, this is the benchmark that gives the lowest DM mass. As mentioned below, this is the sort of points expected in particle models with contact interaction or heavy mediators. Finally, benchmark S4 describes a velocity-averaged cross section whose peak velocity is around 50 km/s and avoids potentially stringent bounds on self-interaction cross section from massive galaxy/cluster observations.

IV Interpreting $a$ and $r_{e}$ in terms of model parameters

In this section, we discuss how the scattering length and the effective range are related to the model parameters of several SIDM scenarios, including those with a light mediator Spergel and Steinhardt (2000); Feng et al. (2010), resonant SIDM Chu et al. (2019), as well as Strongly Interactive Massive Particles (SIMP) Dolgov (1980); Carlson et al. (1992); Hochberg et al. (2014).

IV.1 Contact interaction

The simplest model discussed for SIDM is

[TABLE]

where $m$ is the DM mass and $\lambda$ the coupling. It leads to a constant cross section (within the Born approximation)

[TABLE]

In the effective-range framework, this happens for $|r_{e}|\ll|a|\ll k^{-1}$ so that

[TABLE]

For the Born amplitude to be trusted, we need $\lambda\lesssim 1$, and hence

[TABLE]

If we believe in the upper limit from the clusters, such contact interaction provides a poor fit to the data. Note that the benchmark point S3 corresponds to a large coupling $\lambda\sim 10^{5}$ for 15 GeV DM, where we can no longer trust the Born approximation. Typically, the $\phi^{4}$ theory requires a UV completion of a strongly-coupled dynamics among dark matter particles.

IV.2 SIDM with a light mediator

Models in which non-relativistic DM is coupled to a boson of mass $m_{\phi}$ predict DM self-interactions mediated by the Yukawa potential of Eq. (5). Using the numerical method discussed in Appendix A, we calculate the $S$-wave scattering length and effective range together with the corresponding exact and approximated cross sections, for both repulsive and attractive cases. The results are shown in Fig. 7, Fig. 8 and Fig. 9, respectively. Notice that fixing $\alpha m/m_{\phi}$ and $v/\alpha$ determines all quantities in units of $m_{\phi}^{-1}$. The figure also gives the ($a$, $r_{e}$) parameters of the Hulthén potential $V(r)=\pm\alpha\delta e^{-\delta r}/(1-e^{-\delta r})$, which has been used to approximate the Yukawa potential by setting $\delta=\sqrt{2\zeta(3)}m_{\phi}$ Tulin et al. (2013b). Both potentials give similar effective-range parameters, and thus similar self-interacting cross sections.

The Born regime:

In this case $\alpha\,m\ll m_{\phi}$ and the phase shift can be found by solving the Schrödinger equation perturbatively. In this way, according to Eq. (36), we have

[TABLE]

which implies

[TABLE]

Therefore, in the limit of very small $\alpha$, the scattering length is negligible, the effective range $r_{e}$ is large and they have opposite signs. This behaviour is clearly shown in Fig. 7. As can be seen from the $k^{2}$ expansion of Eq. (15), even in this Born regime, the effective range formula $k\cot\delta_{0}=-1/a+r_{e}k^{2}/2$ can only approximate the $S$-wave phase shift for $k\ll m_{\phi}$. The opposite case is the classical regime mentioned before, where higher partial waves have to be taken into account.

The resonant regime:

Now we turn to the parameter regime satisfying $\alpha m\gtrsim m_{\phi}$, where non-perturbative effects play an important role. In particular, the attractive case exhibits a very rich phenomenology. This is in sharp contrast to the repulsive case shown in Fig. 8, where the scattering cross section simply increase with larger couplings.

For attractive interactions, Fig. 7 shows that as $m\alpha/m_{\phi}$ gradually increases, a critical value is reached at which the scattering length goes to negative infinity. This corresponds to a phase shift approaching $\pi/2$ from below and a cross section of $\sigma_{0}=4\pi/k^{2}$. Notice that this is the maximum value allowed by unitarity. Immediately after $m\alpha/m_{\phi}$ exceeds such a critical value, $a$ becomes positively infinite.

Then, with even larger $m\alpha/m_{\phi}$, the scattering length starts to decrease, until it reaches zero, corresponding to a phase shift of $\pi$. This is the so-called antiresonance, where the cross section takes its minimum value. Further increasing $m\alpha/m_{\phi}$ leads to negative values for $a$, which eventually approaches negative infinity again. The same cycle repeats itself indefinitely. Notice that this behavior of the cross section is responsible for the peak structure observed in Figs. 1 and  9. As mentioned in the previous section, there is a close connection between those peaks, where $|a|\to\infty$, and the bound states that are formed due to the Yukawa potential. Below, we discuss this and how they are related to the poles of the scattering amplitude.

IV.3 SIDM via bound states or virtual levels

Eq. (4) has the following poles

[TABLE]

Even though the poles are in general complex, they can influence the low-energy scattering if they are sufficiently close to the incoming particle momentum. In fact, a close inspection of the Schrödinger equation allows us to interpret them in terms of physical states.888For a textbook review of these topics, see Landau & Lifshitz Landau and Lifshits (1991). In Fig. 10 we plot the real and the imaginary parts of $k^{\text{pole}}_{+}$, which is the closer pole to the real axis.

For simplicity, let us consider first the case of a pure imaginary $k^{\text{pole}}$. The corresponding energy $E=(k^{\text{pole}})^{2}/(2m_{\star})$ is negative, indicating the existence of a bound state with binding energy $\epsilon=-E$. Eq. (17) then leads to

[TABLE]

This formula is remarkable. Take as an example the case of proton-neutron system in the spin-1 configuration. The values quoted above ($a=5.42\,\mathrm{fm}$ and $r_{e}=1.75\,\mathrm{fm}$), which characterize the velocity dependence of the cross section $\sigma_{pn}$, can be used to solve for the binding energy of the deuteron. The result is in perfect agreement with the observed value of $\epsilon=2.2\,\mathrm{MeV}$.

Nonetheless, not every pole is related to a bound state. The latter are only associated with $k^{\text{pole}}=i|k|$ (i.e., ${\rm Im}\ k>0$). Poles with a negative imaginary part correspond to either virtual levels ($k^{\text{pole}}=-i|k|$) or resonances ($k^{\text{pole}}=\kappa_{d}-i|\gamma_{d}|$).999For $\kappa_{d}\neq 0$, the imaginary part of $k^{\text{pole}}$ cannot be positive to conserve the total probability, see e.g. Sitenko and Shepherd (1971). An example of the former is given by the collision of neutrons ($a=-18.9\,\mathrm{fm}$ and $r_{e}=2.75\,\mathrm{fm}$). No bound state of two neutrons exists in nature. In fact, the state inducing such scattering is a virtual level.

The relevance of this for SIDM is that if DM forms a bound state, as predicted in many well-motivated scenarios, the corresponding binding energy would be related to the parameters that determine velocity dependence of the self-interaction cross section by means of Eq. (18). This is particularly true for the Yukawa potential. At the peaks of the cross section (see e.g. Fig. 1), we found that $|a|\to\infty$. Eqs. (17-18) in turn suggest $|a|\sim 1/\sqrt{2m\epsilon}$ with $\epsilon\to 0$. The peaks in the cross section are thus related to the existence of nearly zero-energy bound states.

Even though the regime associated with such bound states is usually referred to as “resonant regime”, we would like to emphasize that there are no intermediate particles produced on shell, i.e., particle resonances. Instead, there are parametric resonances in the sense that for certain parameter combinations the cross section saturates the unitarity limit, where $\epsilon$ approaches zero. All this explains why the presence of poles with very small $\epsilon$ affects the scattering cross section dramatically.

IV.4 Resonant SIDM

In the case of a particle resonance mediating the self-scattering, it is straightforward to see that the kinetic energy $E=(k^{\text{pole}})^{2}/(2m_{\star})$ is complex, with its real and imaginary parts corresponding to the energy above the threshold $E_{R}$, and the decay width of the resonant state $\Gamma(E)$, respectively. More precisely, for the $\ell=0$ case, $E=E_{R}-i\Gamma(E)/2$, which together with Eq. (7) leads to the well-known formulas

[TABLE]

To conclude, when the scattering is induced by a bound state, a virtual level or a resonance, this shows up as momentum poles in the complex $k$ plane. Depending on the sign of the scattering length and the ratio $r_{e}/a$, the effective range theory allows to predict which one actually takes place. In fact, one can elaborate further on the nature of the intermediate state in the scattering process using $r_{e}$ and $a$. For instance, in the context of the deuteron, Weinberg showed that one can infer whether the intermediate state is composite or not from the sign of the effective range Weinberg (1965). Discussing these interesting topics lies beyond the scope of this work.

IV.5 SIMPs

The Strongly Interacting Massive Particle (SIMP) is a proposal where the thermal freeze-out occurs by a $3\rightarrow 2$ transition, which is important when the dynamics is strongly coupled, hence the name Hochberg et al. (2014). It can be naturally realized in QCD-like gauge theories where pions interact via the Wess–Zumino–Witten term Hochberg et al. (2015). Many variations and mediation mechanisms are discussed in the literature Yamanaka et al. (2014); Lee and Seo (2015); Hansen et al. (2015); Bernal and Chu (2016); Bernal et al. (2015, 2016); Choi and Lee (2015); Hochberg et al. (2016); Kuflik et al. (2016); Choi and Lee (2016); Choi et al. (2016); Pappadopulo et al. (2016); Farina et al. (2016); Dey et al. (2016); Cline et al. (2017); Choi et al. (2017a, b); Chu and Garcia-Cely (2017); Choi et al. (2018); Berlin et al. (2018); Choi et al. (2019); Bhattacharya et al. (2019).

The SIMP mechanism prefers dark matter mass in the range from 100 MeV to GeV, and is in marginal conflict with the cluster data as seen in Fig. 5. On the other hand, the strong dynamics often leads to existence of real resonances, bound states, and/or virtual levels, which can improve the agreement by suppressing the cross section at high velocities. In fact, such resonances in QCD-like models of SIMPs are possible McGehee et al. (2019).

V Improving the effective-range approximation

Although Figs. 8 and 9 show that the effective-range approximation works remarkably well in large portions of the parameter space of the Yukawa potential, they make clear that the approximation fails close to the parameter points where the cross section vanishes, i.e., at the anti-resonances. In fact, for realistic $S$-wave Breit-Wigner resonances, it may not work for all possible values of the momentum. Likewise, so far we have used the effective range approach to discuss the non-relativistic DM scattering induced by short-range interactions while inelastic scatterings have been neglected to make sure that the potential, as well as the phase shift, is always real. In this section, we demonstrate that all these effects can be properly described by extending the effective range formalism.

V.1 Antiresonances

$S$-wave antiresonances are probably the simplest example where the effective range formalism fails. In contrast to the prediction of Eq. (4), the scattering amplitude and the cross section vanish at a particular value of the momentum, but not everywhere. One possible way to account for this is to decompose the total phase shift into two pieces –one of them satisfying the effective range approximation– in such a way that they interfere destructively. More precisely,

[TABLE]

which leads to the scattering amplitude

[TABLE]

vanishing at certain values of the momentum, as required. Note that here ${\mathfrak{a}}$ and $\mathfrak{r_{e}}$ are not the standard scattering length and effective range.

At small $k$, this is equivalent to $f_{0}(k)=R+(-1/\mathfrak{a}+\mathfrak{r_{e}}k^{2}/2-ik)^{-1}$, as suggested in $\S$ 134 of Ref. Landau and Lifshits (1991). However, the latter expression does not respect unitarity because $|e^{2i\delta_{0}(k)}|\neq 1$ as follows from Eq. (2). In contrast, our parametrization of Eq. (21) respects unitarity manifestly.

Fig. 11 illustrates our parametrization for the first antiresonance of the attractive Yukawa potential. The difference between the numerical result and the approximation based on Eq. (21) is imperceptible.

V.2 Sharp resonances

Eq. (21) can also describe sharp resonances (second term) accompanied with a continuum piece (first term), while the standard effective range approximation can only describe one of the two. To illustrate this fact, let us consider the potential

[TABLE]

which we depict in Fig. 12 for a particular parameter choice.

This potential gives rise to positive-energy bound states that decay through quantum tunneling. These are real resonances in contrast to the peaks associated with the Yukawa potential in Fig. 1, as explained above. In fact, as shown in the left panel of Fig. 13, the scattering cross section exhibits a resonant enhancement for certain values of the momentum. These correspond to the formation of unstable bound-states. Interestingly, while the standard effective range formula describes the continuum part of cross section fairly well, it fails to describe scattering cross section at the peak. In spite of this, the improved formula approximates the exact result very well. Notice that, as explained above, such an improvement is in practice adding a continuum piece to cross section (except for negative and positive interference just below and above the resonance), which has already been considered in phenomenological studies of resonant SIDM Chu et al. (2019).

In the same fashion, Eq. (21) can simultaneously describe a resonance and an antiresonace as it is the case for potential in Eq. (23) for certain points of the parameter space. This is shown in the right panel of Fig. 13 for another parameter choice. In this case, the antiresonance is induced by the destructive inferences between the continuum and resonance parts.

V.3 Inelastic scatterings

Even though we have assumed that inelastic scatterings –such as DM annihilation or radiative capture March-Russell and West (2009); Braaten and Hammer (2013)– are subleading, in principle they can play a role. For instance, in the vanilla light-mediator model, DM typically annihilates into the mediators, leading to a complex potential Blum et al. (2016). In this case the corresponding phase shift can be decomposed as $\delta_{\ell}=\text{Re}\,\delta_{\ell}+i\text{Im}\,\delta_{\ell}$. Assuming $|\delta_{\ell}|\gg\text{Im}\delta_{\ell}\geq 0$, one can adopt the $S$-wave effective range approximation in the following way

[TABLE]

where $a$ and $r_{e}$ contain a subleading imaginary component Mitroy (2002). If we further neglect the imaginary component of $r_{e}$, the $S$-wave annihilation and scattering cross sections are related via

[TABLE]

Note that such expression does not violate the unitarity limit as its last factor vanishes at $|a|\to\infty$. Moreover, it shows that $\sigma_{an,\,0}(k)\,k$ becomes constant at $k\ll 1/|a|$ as long as the effective range approximation applies. For more discussions on parametrizing the relation between elastic and inelastic cross sections, see e.g. Braaten and Hammer (2013); Blum et al. (2016).

VI Summary and Outlook

In this work, we have studied the effective range approach as a model-independent way to parametrize DM scattering cross sections in astrophysical halos. While it only contains two parameters besides the DM mass, it provides a good description of the self-scattering in most appealing SIDM scenarios, including SIMPs, SIDM with a light mediator and resonant SIDM models.

Starting with a brief introduction to the effective-range approach, we have studied the astrophysical implications. In general, there exists a velocity scale, $(m|a|)^{-1}$, below which the scattering cross section can be treated as a constant. For velocities well above this scale, the cross section quickly decreases. Taking bounds derived from current cluster observations, we have reached the conclusion that DM masses below several GeV are excluded for $\sigma/m\sim 10$ cm 2/g in dwarf-sized halos. See Figs. 4 and  5. The tentative non-vanishing values of $\sigma/m$ –extracted from observational data at various scales– can also be fit in terms of the effective-range parameters. See Fig. 6. Our results suggest that more precise measurements and better extractions are needed to identify or constrain such parameters. Nevertheless, since the value of $r_{e}/a$ only affects the cross section around $v\sim 1/(m|a|)$, obtaining its value from future observations might be challenging.

In addition, we have further investigated the correspondence of the scattering length $a$ and the effective range $r_{e}$ to the model parameters of several popular SIDM scenarios. In general, the scattering cross sections calculated from $a$ and $r_{e}$ agrees well with the exact values as long as the range of the interaction is sufficiently short. Moreover, the effective-range approach demonstrates that significant enhancements in the self-scattering cross section are induced by the poles in the complex plane of the DM momentum. In analogy to nuclear physics, such poles can be interpreted as intermediate physical states, such as a bound state, a virtual level or a resonance.

In the end, we have briefly commented on possible extensions of the effective-range approach, especially for the cases that contain anti-resonances or sharp resonances. Besides, we have also shown that it is possible to study subleading inelastic processes such as DM annihilations using the same framework.

We believe the effective-range approach provides a simple, yet very useful, parametrization to consistently take into account the velocity dependence of DM self-interactions in cosmological simulations involving different astrophysical scales. For instance, such velocity-dependence may play an important role in better understanding the evolution of the sub-halos that move inside the Milky Way halo. This is left for future work.

Acknowledgements

We thank Kai Schmidt-Hoberg for discussions. X.C. is supported by the ‘New Frontiers’ program of the Austrian Academy of Sciences. C.G.C. is supported by the ERC Starting Grant NewAve (638528). X.C. and C.G.C. thank the Erwin Schrödinger International Institute for hospitality while this work was completed. H.M. thanks the Alexander von Humboldt Foundation for support while this work was completed. H.M. was supported by the NSF grant PHY-1638509, by the U.S. DOE Contract DE-AC02-05CH11231, by the JSPS Grant-in-Aid for Scientific Research (C) (17K05409), MEXT Grant-in-Aid for Scientific Research on Innovative Areas (15H05887, 15K21733), by WPI, MEXT, Japan, by the Binational Science Foundation (grant No. 2016153), and by Hamamatsu Photonics K.K.

Appendix A The effective range theory

The phase shifts associated with the self-scattering of DM particles are obtained by solving the Schrödinger equation for the radial wavefunction $R_{\ell,k}(r)$ of the reduced DM two-particle system. This is given by

[TABLE]

together with a boundary condition demanding that $rR_{l,k}$ must vanish at $r=0$. In fact, close to the origin it is expected that the angular-momentum term dominates for sufficiently well-behaved potentials, in which case $R_{l,k}\propto r^{l}$. At large distances from the origin, the potential vanishes and the wave function must be that of a free particle, i.e., a superposition of two spherical waves. The phase shift, $\delta_{\ell}$, parametrizes such a superposition. More precisely, at $r\to\infty$ the asymptotic behavior $R_{\ell,k}(r)$ is given by

[TABLE]

where $j_{\ell}$ and $n_{\ell}$ are respectively the spherical Bessel functions of first and second order.

A.1 A simple method to find the phase shift

In the SIDM context, Ref. Tulin et al. (2013b) presented a systematic method for solving Eq. (26). Here we would like to point out a simpler possibility that will not only provide a powerful method to solve for the phase shift but will also allow us to define the scattering length and the effective range. Let us first define

[TABLE]

Simple algebra shows that

[TABLE]

The fact that $R_{\ell,k}\propto r^{\ell}$ and Eq. (27) fix the boundary conditions of this differential equation to

[TABLE]

Notice that $j_{\ell}(kr)\propto k^{\ell}$ and $n_{\ell}(kr)\propto k^{-(\ell+1)}$ in the limit $k\to 0$, which together with Eq. (29) imply that $\tan\delta_{\ell}\propto k^{2\ell+1}$ for small momenta. The corresponding coefficient of proportionality defines scattering length $a_{\ell}$. More precisely,

[TABLE]

The function $k^{2\ell+1}\cot\delta_{\ell}$ is thus analytic at $k=0$. The next-to-leading term determines the effective range, $r_{e,\ell}$, by means of

[TABLE]

As a by-product we have found a powerful method to solve for the phase shift.101010Simple changes of variable on Eq. (29) allow to simplify the method further. For instance, in the presence of resonances –for which the angle $\delta_{\ell}$ goes beyond $\pi/2$ and its tangent takes values in different branches– it is more convenient to use

$\frac{d\delta_{\ell,k}(r)}{dr}=-k\,m\,r^{2}V(r)\left(\cos\delta_{\ell,k}(r)j_{\ell}(kr)-\sin\delta_{\ell,k}(r)n_{\ell}(kr)\right)^{2}\,,$

(33)

with

$\displaystyle\delta_{\ell,k}(0)=0$ and $\displaystyle\delta_{\ell,k}(r)\to\delta_{\ell}\quad\text{at}\quad r\to\infty\,.$

(34)

Alternatively, using the spherical Hankel function of first kind, $h^{(1)}_{\ell}$, one finds the even simpler formula

$\frac{d\delta_{\ell,k}(r)}{dr}=-k\,m\,r^{2}V(r){\text{ Re}}\left[e^{i\delta_{\ell,k}(r)}h^{(1)}_{\ell}(kr)\right]^{2}\,,$

(35)

which must be solved together with Eq. (34). Due to the fact that $h^{(1)}_{\ell}(r)\propto r^{-(\ell+1)}$, this algorithm is numerically unstable for very large $\ell$. Since those values are only relevant in the classical regime, in that case it might be more efficient to use the analytical formula for $\delta_{\ell}$ from classical physics Landau and Lifshits (1991).

In fact, it is numerically much more efficient to integrate Eq. (29) than to integrate Eq. (26), not only because the former is of first order but also because solving Eq. (29) does not require matching the solution to a plane wave at infinity in order to find the phase shift.

Eq. (29) can be solved by expanding on the potential, with the first term determining the Born regime, which is given by

[TABLE]

A.2 The $S$-wave case

Let us take $\ell=0$ and introduce $u_{k}(r)=rR_{k\,,0}(r)$. Then, Eqs. (26) and (27) read

[TABLE]

and

[TABLE]

Here we have chosen a convenient normalization factor for $u_{k}$. In the following we will find it useful to employ the previous definition of $\psi_{k}(r)$ for any positive value of $r$. Simple algebra proves that for any potential

[TABLE]

Moreover, using the fact that $\psi_{k}(r)$ is the solution of the Schrödinger equation for $V(r)=0$, we find that

[TABLE]

Notice that $\psi_{0}(r)=1-r/a_{0}$, where $a_{0}$ is the scattering length. Subtracting Eq. (39) from Eq. (40), taking $r\to\infty$ and using the fact that $u_{k}$ and $\psi_{k}$ approach to each other in that limit, we find that

[TABLE]

where

[TABLE]

This is the original expression found by Bethe Bethe (1949), who argued that the expansion in Eq. (41) approximates the phase shift with a great accuracy because $\psi_{k}$ and $u_{k}$ differ only where the potential is non-negligible. Note that this conclusion is based on the assumption that in this region both wave functions depend very weakly on $k$, which is generally true, when the potential energy is much larger than kinetic energy and $kr$ is small.

To qualitatively understand the effective range, one can consider the following upper bound, which is valid for potentials that effectively vanish at distances greater than certain range $R$ Phillips and Cohen (1997), as it is the case of the Yukawa potential. Then, $\psi_{0}$ and $u_{0}$ in Eq. (42) coincide for $r\gtrsim R$, which implies that

[TABLE]

Moreover, for shallow attractive potentials, $u_{0}$ behaves like a slowly-varying sine function, where mostly $u_{0}\lesssim\psi_{0}$, resulting in a positive $r_{e,0}$. See also $\S$. 133 of Landau and Lifshits (1991).

A.3 The Hulthén Potential

In the main text, it has been mentioned that the Hulthén potential

[TABLE]

approximates well the Yukawa potential if one sets $\delta=\sqrt{2\zeta(3)}m_{\phi}$, where $\alpha$ gives the coupling and $m_{\phi}$ is the mediator mass of the Yukawa potential Cassel (2010). The advantage of employing the Hulthén potential is that its corresponding Schrödinger equation is analytically solvable, and yield the $S$-wave phase shift Cassel (2010); Tulin et al. (2013b)

[TABLE]

Here, the dimensionless function $\lambda_{\pm}$ is given by $1+imv/(2\delta)\pm\sqrt{\alpha m/\delta-m^{2}v^{2}/(2\delta)^{2}}$.

Using Eqs. (41) and (42), one can obtain the analytical expressions of the $S$-wave effective-range parameters as

[TABLE]

where $\eta=\sqrt{\alpha m/\delta}$, $\psi^{(n)}(z)$ are the polygamma functions of order $n$ and $\gamma\simeq 0.5772$ is the Euler-Mascheroni constant.

Appendix B Velocity-averaged cross sections

The averaged cross section $\langle\sigma v\rangle$ can be calculated in terms of

[TABLE]

with

[TABLE]

We also show the ratio of $\langle\sigma v\rangle$ and $\sigma(\langle v\rangle)\langle v\rangle$ as a function of $am\langle v\rangle$ and $r_{e}/a$ in Fig. 14, which shows that within the effective-range approach both coincide at one percent level, except for the large-velocity regime $\langle v\rangle\gg(|a|m)^{-1}$, where $\sigma$ can be sensitive to $v$.