Dark Matter microphysics and 21 cm observations

Dark matter interactions with massless or very light Standard Model particles, as photons or neutrinos, may lead to a suppression of the matter power spectrum at small scales and of the number of low mass haloes. Bounds on the dark matter scattering cross section with light degrees of freedom in such interacting dark matter (IDM) scenarios have been obtained from e.g. early time cosmic microwave background physics and large scale structure observations. Here we scrutinize dark matter microphysics in light of the claimed 21 cm EDGES 78 MHz absorption signal. IDM is expected to delay the 21 cm absorption features due to collisional damping effects. We identify the astrophysical conditions under which the existing constraints on the dark matter scattering cross section could be largely improved due to the IDM imprint on the 21 cm signal, providing also an explicit comparison to the WDM scenario.


I. INTRODUCTION
Interacting Dark Matter (IDM) with standard model light or massless degrees of freedom, such as photons and neutrinos, gives rise to a suppression of the small-scale matter power spectrum [1][2][3][4] (see also Refs. [5][6][7][8][9][10] for interactions with new (dark) light degrees of freedom). This damping is similar to the one caused by the free streaming of warm dark matter (WDM). In IDM scenarios, in contrast, the suppression of the small-scale overdensities is due to collisional damping [11][12][13]. These two alternatives to the standard ΛCDM model fall into the category of non-cold dark matter scenarios [14] (NCDM). Such dark matter models can provide some solutions to the the ΛCDM (where dark matter is made of purely cold and collisionless dark matter particles) small scale crisis (see e.g. the review of Ref. [15]). A large number of studies in the literature have been devoted to constrain the NCDM picture by means of cosmological probes, such as Cosmic Microwave Background (CMB) fluctuations and spectral distortions, galaxy clustering and Lyman-α forest power spectrum, the number of Milky Way satellites, the reionization history, or gravitational lensing [1-4, 13, 14, 16-58].
In this regard, the 21 cm signal offers a new cosmological probe, complementary to the existing ones, that could open a new window on the early universe and can further test the imprint of NCDM (see e.g. [59,60] for early work on the subject). Here we will focus on the cosmic dawn period and in particular on the first claimed detection of an absorption feature in the sky-averaged global 21 cm signal at a redshift z ∼ 17 by the experiment to Detect the Global Epoch of Reionization Signatures (EDGES) [61]. The measured amplitude of the dip in the 21 cm global signal appears to be much deeper than that expected in standard CDM scenarios and therefore requires new physics to heat the radio background or cool the gas temperature, see also [62,63] for other possible interpretations. This has triggered a surge of interest from the dark matter community trying to relate this effect to dark matter decay and annihilation [64][65][66][67] 1 , and investigating the dark matter scenarios that could account for the signal, see e.g. [71][72][73][74][75][76].
While NCDM scenarios are unlikely to explain the large absorption amplitude, they delay structure formation and, therefore, might delay the onset of reionization and of UV and X-ray emission, see e.g. [43,44,58,59]. As a result a shift to later times in the typical features in the 21 cm sky-averaged global signal and power spectrum is observed in the context of non-cold dark matter [58,59,[77][78][79]. Consequently, it is timely to study the compatibility between the observation reported by EDGES, located at a redshift around z 17 and the IDM scenario. We follow two possible avenues. The first of them relies on exploiting the non-negligible Lymanα coupling between the gas and the spin temperature characterizing the 21 cm signal at z 20 [78,79]. The second one consists in imposing the minimum in the absorption feature to happen before z 17 [77]. For both strategies, a number of degeneracies between the details of dark matter microphysics and the astrophysical parameters will appear. We will briefly discuss their impact on the constraints on NCDM scenarios.
Universal fits to the halo mass functions from N-body simulations for the IDM scenarios have been obtained in [3,80,81]. In particular, here we will use the results of [80] derived for IDM scenarios involving dark matterphoton scatterings. A possible particle physics model related to this cosmological scenario is the case of millicharged dark matter [82,83]. IDM scenarios including dark matter-neutrino scatterings have been shown to give rise to a very similar damping in the power spectrum and also to a very similar mass function as for dark matterphoton scatterings, see Refs. [3,12] and the Appendix A. Unfortunately, in the latter case, no publicly available dedicated analysis provides the necessary fits to the associated halo mass functions necessary for our study. We therefore use the IDM scattering on photons as a toy model to evaluate the impact of the EDGES signal on the more general case of IDM with light degrees of freedom. In order to ease comparison with previous studies, we shall also study the case of thermal warm dark matter (WDM) with mass in the keV range (see Refs. [33,34,[42][43][44][45][46][47][48][49][50][51][52][53][54][55][56][57][84][85][86][87][88][89][90][91][92][93][94][95][96][97][98][99][100][101][102] and the most recent works of Refs. [77,79]).
The structure of the paper is as follows. We start in Sec. II by describing the physics of the 21 cm global signature. We account for the effect of IDM in the 21 cm global signature in Sec. III, presenting the constraints on the dark matter photon elastic cross sections arising from (i) the presence of a rich Lyman-α background at z 20 (see Sec. III.1), and (ii) the location of the EDGES minimum (see Sec. III.2). Finally, we summarize our results and conclude in Sec. IV.
II. THE 21 CM SIGNAL II.1. The differential brightness temperature The brightness of a patch of neutral hydrogen (HI) relative to the CMB at a given redshift z is expressed in terms of the differential brightness temperature, δT b . The sky-averaged δT b scales as [103][104][105][106] where ν = ν 21 /(1 + z) with ν 21 = 1420 MHz, x HI represents the fraction of neutral hydrogen and Ω b h 2 and Ω m h 2 refer to the current baryon and matter contributions to the universe's mass-energy content. The ratio of the populations of the two ground state hyperfine levels of hydrogen is quantified by the spin temperature, T S , which is determined by three competing effects [107]: 1) absorption and stimulated emission of CMB photons coupling the spin temperature to the CMB temperature T CMB in contrast with 2) atomic collisions (which are important at high redshifts z 30); and 3) resonant scattering of Lyman-α photons that couple the spin temperature to the gas kinetic temperature T k . The latter effect is the so-called Wouthuysen-Field effect [108,109] that turns on with the first sources. Assuming that the Lyman-α color temperature is T α T k [105], the spin temperature results from: with x tot = x α + x c and x c and x α are the coupling coefficients for collisions and Lyman-α scatterings. At the low redshifts of interest here, collision coupling effects can be safely neglected and therefore, x tot = x α . The Lyman-α coupling is defined as where T = hν 21 = 68.2 mK is the hyperfine energy splitting, e and m e the charge and mass of the electron, f α is the oscillator strength of the Lyman-α transition, A 10 is the spontaneous decay rate of the 21 cm transition, J α is the specific intensity of the background radiation evaluated at the Lyman-α frequency and S α is an order unity correction factor which accounts for the detailed shape of the spectrum near the resonance [110]. In particular, in the framework considered here, S α 1 and the Lyman-α flux gets two types of contributions. One results from the X-ray excitation of HI (J αX ), while the other one results from direct stellar emission of UV photons between Lyman-α and the Lyman-α limit J α , thus J α = J αX + J α [111]. We will see in Sec. II.2 that for the X-ray efficiencies considered here, J αX only represents a small contribution to the total Lyα flux. On the other hand, the direct stellar emission contribution to J α is computed assuming by default a Pop II stars spectral model. This gives rise to the emission of 9690 photons per baryon between Lyman-α and the Lyman limit, see Appendix B (see also [112]). Notice from Eq. (2.2) that, when x tot = x α = 1, δT b will be at the half of the value that it would have if T S were completely coupled to T k , which happens when x α 1. The authors of Ref. [79], following the EDGES results [61], have imposed that x α should be one or larger at redshift z 20. We shall apply this constraint in our numerical analyses of IDM scenarios, see Sec. III.1.
In order to extract the imprint of NCDM on the 21 cm signal, we profit from the publicly available tool 21cmFast. The main purpose of the code is the study of variations in the 21 cm signal due to changes in a given set of astrophysical and cosmological parameters. We make use of the output values of x α and of T S and x HI to extract the differential brightness temperature as in Eq. (2.1). In practice, we use a version of the 21cmFast code adapted to account for the IDM and the WDM as detailed in Ref. [58]. We have modified the default WDM scenario implementation modifying the definition of both the transfer function and the halo mass function according to the prescription given in the Appendix A. This halo mass function plays an important role in the evaluation of the production rate of ultraviolet (UV), X-rays, and Lyman-α radiation, responsible of the ionization, heating and Lyman-α coupling respectively. These production rates are proportional to the star formation rateρ . In 21cmFast, this quantity is evaluated in terms of the growth of the fraction of mass collapsed in haloes which are able to host star-forming galaxies, f coll (> M min vir ), defined as where ρ m,0 is the current matter density and dn/dM is the halo mass function. For the NCDM cosmologies, the halo mass function is always suppressed at small masses compared to the CDM scenario, giving rise to a smaller f coll , at fixed redshift, as illustrated in Fig. 1, see also [1,59,113]. 2 With the purple area we show the case of CDM with threshold masses M min vir between 10 6 M (upper curve) and 3 × 10 7 M (lower curve) at z = 20 or equivalently T min vir between 10 3 and 10 4 K, see eq. (2.5). The upper purple curve for CDM can be compared to case of NCDM scenarios in the form of WDM (red, blue and cyan continuous curves) and of IDM (red, blue and cyan dotted curves) for the same M min vir = 10 6 M at z = 20 (T min vir = 10 3 K). The IDM scattering cross-sections are normalized in terms of the Thompson cross section σ T = 6.65 × 10 −25 cm 2 .
Overall, in the framework considered here, the NCDM modification of f coll will result into a delayed reionization, heating and Lyman-α coupling, giving rise to an absorption feature in the 21 cm signal located at lower redshifts with respect to CDM scenarios [43,44,58,59]. Notice also that curves of a fixed color in Fig. 1 correspond to IDM and WDM models giving rise to a fixed value of the half-mode mass or breaking scale of the linear power spectrum, see Eqs. (A6), (A7) and (A9) in Appendix A. For fixed half-mode mass, the suppression of the f coll is always more severe in the WDM scenarios than in the IDM scenarios due to a relatively larger number of low mass haloes in the IDM case, see [80]. The rate at which f coll increases with redshift is also different between IDM and WDM. These features might help to discriminate between NCDM models using both the 21 cm global signal and its power spectrum [58].

II.2. The astrophysical parameters
The minimum virial mass M min vir from which haloes begin to efficiently form stars (see Eq. (2.4)) is related to the threshold temperature T min vir as [115] M min vir (z) = 10 8 1 + z 10 5) where µ is the mean molecular weight and it is equal to 1.2 (0.6) for a neutral (fully ionized) primordial gas. The value of the minimum virial temperature depends on the cooling mechanism considered. The lower threshold to make the atomic cooling channel effective is T min vir = 10 4 K [116][117][118][119][120], corresponding to M min vir = 3 × 10 7 M at z = 20 for µ = 0.6. In contrast, the molecular H 2 cooling channel can be effective down to temperatures of T min vir = 10 3 K corresponding to M min vir = 10 6 M at z = 20. The hydrogen molecules may however be destroyed by the ionizing radiation. Nevertheless, several hydrodynamical works in the literature [121,122] have shown that in the presence of a large soft UV background, molecular cooling could be highly effective, i.e. metal-enriched star formation is not restricted to atomic cooling. Furthermore, even in the absence of electrons, molecular cooling could cool down the gas in haloes associated to virial temperatures much lower than the ones required for atomic cooling. Consequently, in the following, we shall consider T min vir = 10 3 K as the minimum threshold temperature for star formation and we assume the same threshold temperature T min vir for haloes host-ing ionizing and X-ray sources. This parameter plays a crucial role in extracting the constraints from the 21 cm absorption signal in both CDM and NCDM scenarios, see also [77,78]. The impact of varying T min vir in the 10 3 K to 10 4 K range within the CDM paradigm is illustrated in Fig. 1 with the purple area. The lower value of the threshold parameter corresponds to higher values of the faction of collapsed haloes at a given redshift.
The comoving star formation rate density is described byρ = f * ρ b,0ḟcoll (> M min vir ), where f * , the fraction of baryons that are converted into stars, considered here as a constant parameter (neglecting any dependence on the halo masses or redshift). As f * controls the amplitude of the star formation rate density, it also sets the amplitude of both the ionizing and heating rates, as well as the Lyman-α flux. This parameter is quite uncertain as no observations of low mass mass haloes of mass 10 6 − 10 8 M at redshift z ∼ 20 are available. Nevertheless, several previous works, based on radiationhydrodynamic simulations of high-redshift galaxies in a neutral medium [122][123][124][125] or based on the comparison of the star formation rate density to the one derived from UV luminosity function measurements [79,126,127], have found values of f * ∼ O(0.01). In the following, we therefore consider f * = 0.01, leaving for the discussion in Sec. IV the impact in our results of slightly larger values of f * .
Finally, in order to characterize the overall normalization of the X-ray luminosity, 21cmFast makes use of the X-ray efficiency parameter ζ X , expressed in units of M −1 . This parameter is varied here within the range ζ X = [1 − 5] × 10 56 /M . One can relate this parameter to the integrated X-ray soft band emissivity (below 2 keV) per unit of star formation rate escaping the galaxy L X<2keV /SFR, varied in the range log 10 (L X<2keV /SFR) ∈ [39.5, 40.2] erg/s / (M /yr). 3 This range is similar to one extracted from observations of the hot interstellar medium, which lead to log 10 (L X<2keV /SFR) ∼ [39,40] erg/s / (M /yr) [129], and also to the one adopted in Ref. [77] (see also Ref. [130]).

III. IMPRINT OF THE IDM ON THE 21 CM SIGNAL
Non-cold dark matter scenarios are expected to delay structure formation and therefore the absorption feature in the 21 cm signal at cosmic dawn. This effect is illustrated for the IDM model under study in Fig. 2, where we show the coupling coefficient for Lyman-α scattering x α (top), and the sky-averaged 21 cm brightness temperature (bottom). The different colors correspond to different values of the dark matter-photon scattering cross-sections over the dark matter particle mass, σ IDM /σ T (GeV/m DM ). As can be noticed, large cross sections, inducing a stronger suppression at small scales, also induce a stronger suppression of the Lyman-α coupling at a given redshift (see the top panel of Fig. 2) and a larger shift towards smaller redshifts of the 21 cm features (see the bottom panel of Fig. 2). We also illustrate, with a vertical line, the position of the minimum of absorption reported by the EDGES experiment (at ν = 78 MHz corresponding to z = 17.2). The dashed lines correspond to the largest signal redshift range at the minimum of absorption within the 99% CL interval reported by Ref. [61]. Figure 2 also shows the impact of varying the X-ray heating efficiency ζ X in the 10 56 M −1 − 5 × 10 56 M −1 range with the width of the colored bands. We see from the top panel of Fig. 2 that the range of ζ X considered here, the Lyman-α flux resulting from X-ray excitation (J αX ∝ ζ X ) is usually a negligible contribution to the Lyman-α coupling at the redshift at which x α ∼ 1. On the bottom panel of Fig. 2, we see the impact of the X-ray heating parameter ζ X on the differential brightness temperature. Here the impact is more significant. The deepest absorption dip corresponds to ζ X = 10 56 M −1 , while the shallowest one is obtained for ζ X = 5 × 10 56 M −1 . Notice also that a larger value of ζ X implies an earlier minimum, i.e. the X-ray heating of the IGM occurs earlier in time. As a result, in order to extract conservative constraints on the NCDM parameters from the redshift at which the absorption minimum is located, we shall consider the value ζ X = 5 × 10 56 /M , see Sec. III.2. Figure 3 illustrates the equivalent to Fig. 2 with ζ X fixed to 5 × 10 56 M −1 and considering this time T min vir = 10 3 K (continuous curves) as well as T min vir = 10 4 K (dotted curves). Going from left to right panel of Fig. 3 we increase the value of f * by a factor of 3. Based on these plots, we discuss in the next two subsections the constraints that could be derived on NCDM scenarios following two different approaches.

III.1.
Constraints from the Lyman-α background In order to account for the EDGES results, the authors of Ref. [79] based their analyses on the assumption that a sufficiently strong Lyman-α background is present by z ∼ 20. They imposed (3.1) This limit results from the observation that the absorption signal reported by the EDGES experiment is equal to half of the maximum amplitude of absorption at z 20. We shall follow this assumption here, applying this condition to our simulations within NCDM models. Let us first focus on the top left panel of Fig. 3, which shows the Lyman-α coupling coefficient x α as a function of the redshift for f * = 0.01. The condition reported in (3.1) disfavors the cosmological scenarios associated to a prediction of x α (z) lying within the shaded area at z = 20. The astrophysical parameter T min vir has a significant impact on the Lyman-α coupling coefficient x α . Larger virial temperature shifts star formation to lower redshifts, giving rise to a lower Lyman-α background at a given redshift. In a similar way, a large value of the IDM scattering cross-section implies a longer delay in structure formation. As a result, one can deduce from Fig. 3 (left panel) that, for molecular cooling (T min vir = 10 3 K), the scattering cross-section must be below 10 −10 σ T × (m DM /GeV). If one assumes instead that the only efficient cooling mechanism is atomic cooling with T min vir ≥ 10 4 K, the limits on the IDM scattering cross section become much tighter, excluding scattering cross-section lower than 10 −11 σ T × (m DM /GeV). Our most conservative bound (assuming f * = 0.01) is therefore σ IDM 10 −10 σ T × (m DM /GeV). This constraint is stronger than the 95% CL upper limit of σ IDM < 8×10 −10 σ T ×(m DM /GeV) reported in Ref. [58], based on observations of Milky Way satellite galaxy number counts and assuming a mass for our galaxy of M MW = 2.6 × 10 12 M .
Notice that the Lyman-α coupling coefficient x α is directly proportional to the fraction of baryons converted into stars (f * ), here considered as constant. The left (right) panels of Fig. 3 have been simulated with f * = 0.01 (f * = 0.03). Notice that larger values of f * increase the Lyman-α coupling at a fixed redshift, weakening the bound on IDM scenarios resulting from Eq. (3.1). In the case of f * = 0.03, the limits quoted above are translated into σ IDM 10 −9 σ T × (m DM /GeV) and σ IDM 10 −10 σ T × (m DM /GeV) for T min vir = 10 3 K and T min vir = 10 4 K respectively.

III.2.
Constraints from the position of the absorption minimum in the 21 cm global signature Another possible avenue to constrain NCDM models using the EDGES observations is based on the location of the minimum of the absorption. Reference [77] imposed that it should be located at a redshift higher than z = 17.2 and studied the resulting bounds on a large set of NCDM models. We show in the bottom panels of Fig. 3 the effect of the IDM scenario considered here on the global sky-averaged 21 cm brightness temperature obtained by means of Eq. (2.1).
Let us focus first in the f * = 0.01 case (i.e. left panel). Considering atomic cooling (dotted curves with T min vir = 10 4 K), it appears that for scattering cross sections larger than ∼ 10 −11 σ T ×(m DM /GeV) the absorption minimum takes place at redshifts lower than z = 17.2. Such crosssections should therefore be regarded as disfavored. Considering molecular cooling softens this constraint by ∼ one order of magnitude (see the continuous curves for T min vir = 10 3 K). Notice that we have considered a conservative X-ray efficiency of ζ X = 5 × 10 56 M −1 for all curves. Considering lower values of ζ X will give rise to a later X-ray heating and thus a minimum of absorption located at lower redshifts.
If instead, we consider a larger fraction of baryons converted into stars, f * = 0.03, (see the bottom right panel), and T min vir = 10 3 (10 4 ) K, the limit on the IDM interactions is relaxed, excluding cross sections above ∼ 10 −9 (10 −10 ) σ T × (m DM /GeV). This is due to the fact that both the X-ray and the Lyman-α emission rates are directly proportional to f * . Increasing f * implies an earlier Lyman-α coupling and X-ray heating periods, displacing the minimum of the absorption in the 21cm signal to a larger redshift.

IV. DISCUSSION AND CONCLUSIONS
Interacting Dark Matter (IDM) models, in which dark matter is not made out of purely cold, pressureless particles, are an interesting alternative to the standard CDM paradigm and provide a possible avenue to alleviate the so-called small-scale crisis of the ΛCDM. The IDM could be scattering off light or massless degrees of freedom such as photons or neutrinos. Here we have considered the case of dark matter scattering on photons, characterized by the size of the scattering rate over the DM mass, σ IDM /m DM . The reason for this choice is driven by the availability of a fitting function for the halo mass function relevant for our study. Let us emphasize though that scatterings on neutrinos are expected to give rise to a similar imprint on the 21 cm signal.
Several studies have constrained IDM models based on their suppression of clustering at small scales, exploiting galaxy power spectrum, gravitational lensing, CMB, number of Milky Way dwarf galaxies and Lyman-α forest observations, among others. Here we focus on the imprint of IDM on the 21 cm signal arising from cosmic dawn. Based on a modified version of the 21cmFast code, our simulations show that IDM delays the formation of haloes capable of star formation, shifting the timing of the 21 cm signal features compared to the standard CDM scenario. A similar effect has been reported in the case of other NCDM models [43,44,58,59]. In this paper, we have considered two possible ways to test the IDM properties against the 21 cm signal. First, following Ref. [79], a significant Lyman-α coupling between the gas and the spin temperature characterizing the 21 cm signal should be present at z = 20. Secondly, as argued in Ref. [77], the location of the absorption minimum in the EDGES signal at z = 17.2 implies that any scenario with sufficiently enough delayed structure formation could be discarded.
We have first identified which are the most relevant astrophysical parameters showing large degeneracies with the IDM scattering cross-section over the mass σ IDM /m DM . Namely, the fraction of baryons into stars f * , the threshold temperature for haloes to host starforming galaxies T min vir and the X-ray efficiency, ζ X , have been shown to interfere with the eventual extraction of a non-zero σ IDM . Fortunately, the parameter ζ X only plays a significant role in extracting the location of the absorption minimum, and we have adopted the conservative value of ζ X = 5 × 10 56 M −1 (corresponding to an integrated soft band X-ray emissivity of L X<2keV /SFR = 10 40.2 erg/s / (M /yr)). We also considered a lower limit on the threshold virial temperature of T min vir = 10 3 K (corresponding to molecular cooling) as well as a constant value of f * < 0.03.
Our results are summarized on the left panel of Fig. 4 in the case of f * = 0.01. The top panels shows the Lyman-α coupling coefficient x α at z = 20 for T min vir = 10 3 (10 4 ) K with the top blue (bottom red) curve as a function of the IDM scattering cross-section. Notice that the x α (z = 20) curves saturate to a maximum value at low enough values of the scattering cross-section. This corresponds to the limiting x α (z = 20) value that one would get in the CDM scenario (indicated with the horizontal dotted lines). In addition, the shaded region refers to the parameter space in which the condition x α (z = 20) > 1 can not be satisfied and can therefore be considered as disfavored by the condition given by Eq. (3.1). As a result, for f * = 0.01 and ζ X = 5 × 10 56 M −1 , a value of σ IDM > 10 −10 σ T × (m DM /GeV) fails to satisfy the condition x α (z = 20) > 1 when considering molecular cooling (T min vir = 10 3 K). This improves the previous bound derived on such IDM model in Ref. [58]. When considering higher threshold temperatures for efficient cooling, the bounds gets tighter, disfavouring even the canonical CDM scenario if T min vir = 10 4 K. All these limits would be relaxed if the fraction of baryons converted into stars is larger, as it was illustrated in Fig. 3 for f * = 0.03. For instance, for T min vir = 10 3 K, the limit would be reduced by one order of magnitude (i.e. σ IDM < 10 −9 σ T × (m DM /GeV)).
In the bottom left panel of Fig. 4, we show the redshift at which the sky-averaged 21 cm brightness temperature exhibits its minimum of absorption at cosmic dawn, z(δT min b ), as a function of the IDM scattering cross-section. Again z(δT min b ) saturates for low enough scattering cross-section to a value corresponding to the the CDM limit denoted with a dotted horizontal line. The shaded area denotes the region in which this minimum is located below z = 17.2. More concretely, in order to derive the curves shown in the bottom panels of Fig. 4, we have fitted the brightness temperature curves obtained from our simulations to the flattened gaussian shape that the EDGES collaboration uses to model the 21 cm absorption profile [61]. This procedure allows us to extract the central frequency for each possible value of σ IDM /m DM , T min vir and f * , and also to constrain the model by imposing that this mean frequency ν 0 of our flattened gaussian fits should lie at frequencies below the lower 99% confidence interval (including estimates of systematic uncertainties) reported by the EDGES collaboration for the ν 0 = 78 MHz parameter. From Fig. 4, the position of the minimum of absorption disfavours the region σ IDM > 10 −10 σ T × (m DM /GeV) (σ IDM > 10 −11 σ T × (m DM /GeV)) for T min vir = 10 3 K (T min vir = 10 4 K). As in the case of the x α (z = 20) > 1 condition, these limits would be relaxed by one order of magnitude in the case of a value of f * = 0.03. Also, for the astrophysical parameters considered here, we see that the constraints from the Lyman-α coupling condition are typically tighter than those arising from the position of the minimum of absorption.
Summarizing, for f * 0.01, ζ X < 5 × 10 56 M −1 and T min vir > 10 3 K, the bounds on σ IDM /m DM derived in this work imply an order of magnitude improvement over the most constraining existing limits in the literature [2,58].
Larger values of f * may compromise these upgraded limits. In order to ease the comparison to other studies on NCDM scenarios, we provide the results obtained following the same methodology in a thermal warm dark matter scenario involving light dark matter particles with a mass m WDM of a few keV. The prescription considered here to describe the suppression of the halo mass function at small halo masses is given in the Appendix A. We show the obtained dependence of x α (z = 20) and z(δT min b ) in the right panel of Fig. 4. From this figure with f * = 0.01 and ζ X = 5 × 10 56 M −1 , we can infer a lower limit in the WDM mass around m WDM > 6 keV (m WDM > 12 keV) if T min vir = 10 3 K (T min vir = 10 4 K). These tight limits are similar to those derived in Refs. [77,79] for slightly different WDM implementations and astrophysical parameters.

Appendix A: Halo mass functions
The halo mass function, that counts the number of haloes per unit halo mass and volume at a given redshift, can be written as [131] where n is the halo number density, ρ m = Ω m ρ c is the average matter density in the Universe at z = 0 and σ 2 = σ 2 (M, z) is the variance of density perturbations, which is a function of the halo mass M and redshift. The first-crossing distribution, f (ν), is expected to be a universal function of ν ≡ δ 2 c /σ 2 (M, z), with δ c = 1.686, the linearly extrapolated density for collapse at z = 0. The Sheth-Tormen (ST) first-crossing distribution reads as [132][133][134]: with p = 0.3, q = 0.707 and A = 0.3222. For the standard CDM scenario, we consider this first-crossing distribution. The variance at z = 0, σ(M ) ≡ σ(M, z = 0), at a given scale R can be expressed as where P X (k) is the linear power spectrum at z = 0 for a given X = {CDM, IDM or WDM} cosmology and W is the Fourier transform of a filter function that we consider to be a top-hat (TH) function in real space (see e.g. [14,77] for the possibility of using a sharp-k window function for NCDM). The redshift dependence is driven by the linear growth function, D(z) normalized to 1 at z = 0, so that the root-mean-square (rms) density fluctuation is σ(M, z) = σ(M, z = 0) D(z).
The transfer function T X for a NCDM scenario X is defined as where P CDM (k) is the linear power spectrum in ΛCDM.
Here we use the prescription of Refs. [12,84]: where µ = 1.2 is a dimensionless exponent and α X is a breaking scale. The latter takes the form: for WDM scenarios [3,90] and IDM involving dark matter-photon scattering [12] respectively. For DM-neutrino interactions, described with the same parametrization of the transfer function, one gets a breaking scale α νDM 0.8 × α γDM for a fixed value of the scattering cross section [3].
The halo mass function defined as in Eq. (A1) is well suited for CDM but it needs to be adapted for the NCDM case. In order to match the results from N-body simulations, the WDM halo mass function can be expressed as [90] where an additional mass-dependent correction to the standard ST formalism appears. We use a = −0.6 and b = 0.5, as obtained in [3,80]. The function dn ST, WDM dM in Eq. (A8) refers to the halo mass function obtained with a ST first-crossing distribution, as defined in Eq. (A2), and a linear matter power spectrum corresponding to the WDM case. In order to describe the suppression in the linear regime, one can consider the half-mode mass M hm , defined as the mass scale for which T X /T CDM = 1/2 (i.e. P X /P CDM = 1/4). Using the general fit to the transfer function, Eq. (A5), the half-mode mass M hm can be easily derived as For what concerns the IDM models, the number of lowmass structures appears to be larger than in WDM scenarios [3] (see also [7]). In order to reproduce the IDM results for masses below the half-mode mass, an extra mass-dependent correction must be introduced to the halo mass function [80]: absorbed by the medium. A photon which reaches a Lyman-n resonance at redshift z has to be emitted at a redshift below 1 + z max,n = 1−(n+1) −2 1−n −2 (1 + z). If f rec (n) is the recycled fraction of the level n, i.e., the probability of generating a Lyman-α photon from the n level [135], the total Lyman-α flux can be written as where the emission frequency is ν n = ν n 1+z 1+z , being ν n = ν LL (1 − n −2 ) and ν LL the Lyman limit frequency. The comoving emissivity α can be written as proportional to the star formation rateρ : wheren b,0 is the comoving number density of baryons.
Assuming that only Population-II stars contribute to this emissivity, the spectral distribution ε(ν) is given by a separate power law between each Lyman-n and Lymann + 1 levels: for ν n ≤ ν ≤ ν n+1 , with N n the number of photons emitted between the n and n + 1 resonances and β n the spectral index [112]. The function above is normalized as νn+1 νn dνε(ν) = N n , with ν LL να dνε(ν) = n N n 9690 the total number of photons emitted between the Lyman-α and the Lyman limit. Although we keep this normalization as constant through our analysis, notice that changes in the number of photons which contribute to the Lyman-α flux could have a deep impact in the 21cm signal [136].