Properties of Subhalos in the Interacting Dark Matter Scenario

One possible and natural derivation from the collisionless cold dark matter (CDM) standard cosmological framework is the assumption of the existence of interactions between dark matter (DM) and photons or neutrinos. Such a possible interacting dark matter (IDM) model would imply a suppression of small-scale structures due to a large collisional damping effect, even though the weakly-interacting massive particle (WIMP) can still be the DM candidate. Because of this, IDM models can help alleviate alleged tensions between standard CDM predictions and observations at small mass scales. In this work, we investigate the properties of the DM halo substructure or subhalos formed in a high-resolution cosmological N-body simulation specifically run within these alternative models. We also run its CDM counterpart, which allowed us to compare subhalo properties in both cosmologies. We show that, in the lower mass range covered by our simulation runs, both subhalo concentrations and abundances are systematically lower in IDM compared to the CDM scenario. Yet, as in CDM, we find that median IDM subhalo concentration values increase towards the innermost regions of their hosts for the same mass subhalos. Similarly to CDM, we find IDM subhalos to be more concentrated than field halos of the same mass. Our work has a direct application to studies aimed at the indirect detection of DM where subhalos are expected to boost the DM signal of their host halos significantly. From our results, we conclude that the role of the halo substructure in DM searches will be less important in interacting scenarios than in CDM, but is nevertheless far from being negligible.


Introduction
The current standard model of cosmology, ΛCDM, is based on a cosmological constant to explain the late-time accelerated expansion of the Universe and a cold dark matter (CDM) component to account for the required additional gravitational attraction to form and support the galaxies and larger structures we observe today [1]. In this framework, the structure of the Universe is formed via a hierarchical, bottom-up scenario (see, e.g., [2]) with small primordial density perturbations growing to the point where they collapse into the filaments, walls and eventually dark matter (DM) halos that form the underlying large-scale-structure filamentary web of the Universe. The galaxies are embedded in these massive, extended DM halos teeming with self-bound substructure. Any viable cosmological model has to successfully predict both the abundance and internal properties of these structures and their substructures, and match the observational data on a wide range of scales. ΛCDM achieves this challenging feat well on the largest scales [3][4][5][6][7]. Yet, on small scales tensions have been reported between its predictions and observations in our local cosmological neighbourhood. The abundance of DM substructures predicted by numerical simulations of structure formation exceeds significantly the number of satellite galaxies observed around the Milky Way and neighbouring Andromeda galaxy (see e.g., [8,9]). Various explanation attempts for this and similar discrepancies such as the "too big to fail", "cusp vs. core" and "satellite alignment" problems [10,11] were brought forward, with some of them attributed to feedback mechanisms in the baryonic sector that suppressed star formation in such small halos (see e.g., [12]), thus leaving them without any observable tracers in the observational surveys [13], or alter the DM profiles within the halos [14][15][16][17][18][19]. Others turned to alternative models for the DM to account for the lower amount of small subhalos (see, e.g., [20,21])) or deviations of their expected properties [22][23][24]. The latter pathway is not only well motivated, as the properties of DM has yet remained largely a mystery, but in return also allows us to use the study of galaxies and their structural properties as effective probes into the very nature of the elusive nature of the DM particle.
One natural derivation from the collisionless CDM in the standard model is the assumption of the existence of interactions between DM and the standard model (SM) particles we know about, in particular, photons or neutrinos [25][26][27]. This does not only affect, as we show in this article, the formation of DM structures on small scales, but also provides an explanation for the exact relic abundance of DM, Ω cdm h 2 = 0.12011, found in the Universe today [1]. With such interactions, DM was in full thermal equilibrium with SM particles at sufficiently early times and then annihilated into SM particles until the DM decoupled from the standard sector as the Universe expanded and cooled down. The cross section needed to retain the observed abundance of DM is surprisingly close to the one expected from the interaction via the weak force in the SM, thus coining the name "weakly interacting massive particles (WIMP) miracle". Beyond-SM theories provide a variety of WIMP DM candidates such as the minimal SUSY standard model with the neutralino and sneutrino and their electroweak scale interactions [28], or the minimal Universal extra dimension model of the Kaluza-Klein (KK) theory with the first excitation mode of the gauge field as the lightest KK-particle [29]. When it comes to the interaction partner, the usefulness of baryons is limited due to their relatively low abundance in the Universe at any time and the existing constrains on the cross section with DM from direct detection experiments. On the other hand, relativistic neutrinos and photons can be found in high abundance in radiation-dominated era of the early Universe and particle-physics experiments, e.g. particle colliders, provide only very few constraints on their potential interaction with DM.
In our work, we do not pick a specific model, but simply work within an effective theory, i.e. an effective interaction term between some unspecified, otherwise sterile DM particles and our SM particles of choice, photons and neutrinos in the Lagrangian. We will refer to this model as interacting dark matter (IDM). Depending on the actual type/mass of the mediator in our "black box", this can lead to a momentum/velocity-dependence of our effective cross-sections but, for simplicity, we mainly focus in the following on velocity-independent scenarios. For any given cross section, the DM remains coupled to the radiation in the early Universe until the latter is diluted enough as the Universe expands for the DM to become decoupled. As a result of this coupling, primordial perturbations and, thus, the seeds of late-time structures, are suppressed within the DM below a certain scale. This is visible as a cut-off in the linear matter power spectrum. For a DM-radiation scattering cross section of σ/σ Th = 2 × 10 −9 (m dm /GeV) with σ Th the Thomson cross section and m dm the DM mass, this characteristic scale is ∼100 kpc [30] and increases or decreases with the cross section [26,[31][32][33][34][35][36].
Returning to the premise of using the halo and subhalo population as a probe into the nature of DM, we can use this suppression and its consequences for the structure formation to find bounds on the interaction cross section. Unfortunately, as previously mentioned, a more direct study of the halo population is difficult as the distribution of its visible tracers i.e. stars and gas is also subject to not fully quantified astrophysical processes. Strong lensing may provide a way to determine the DM profile of larger halos [37], but the halos around the cut-off scale are orders of magnitude smaller. Indirect methods on the other hand, namely, the detection of the annihilation or decay products of DM particles, are highly dependent on the statistical and structural properties of the halo and subhalo population. For instance, the extragalactic γ-ray and neutrino signals due DM annihilations, when estimated via the so-called halo model [38][39][40], depend mainly on the DM halo and subhalo structural properties as well as their abundances (see e.g., [41][42][43][44][45]). Clearly, the considered cosmological model is crucial for such DM searches as different predictions for structure formation on small scales imply different gamma-ray or neutrino signal estimations. Ultimately, this may translate into different constraints on the DM annihilation cross section when compared to those obtained assuming the standard ΛCDM scenario. In [46], the isotropic extragalactic signals expected from DM annihilations into γ-rays and neutrinos were investigated for both IDM and ΛCDM models using only main halo properties as extracted from DM-only simulations. In this work, we study the properties of the halo substructure in the same IDM scenario of ref. [46], for which we now use a set of N-body, DM-only cosmological simulations with higher particle resolution.
The work is organized as follows. We briefly summarize the theory behind IDM in section 2 followed by a description of our simulations in section 3. For both IDM and ΛCDM models, in section 4 we present our results for subhalo properties such as concentrations, abundances and subhalo radial distributions within the host halos. We finally discuss these results and draw our conclusions in section 5.

Interacting dark matter
In our effective theory of IDM, the interactions between DM and photons (or alternatively neutrinos) result in additional terms in the linearized Euler equations governing the evolution of the cosmic componentsθ where ψ is the gravitational potential, H is the conformal Hubble rate, c s is the baryon sound speed and δ, θ and σ are the density, velocity divergence and anisotropic stress potential respectively, associated with the baryon (b), photon (γ) and DM fluid. For the electromagnetic interactions (EM) in the SM, the first two equations include terms with the Thomson scattering rateκ ≡ aσ Th cn e , where c the speed of light and n e the density of free electrons (the scale factor a, appears since the derivative is taken with respect to conformal time). The ratio of the baryon to photon density, R ≡ (3/4)(ρ b /ρ γ ), is a pre-factor to ensure momentum conservation. C DM−γ and C γ−DM = −S −1 C DM−γ are the new interactions terms that have to be added to include interactions between DM and the cosmic photon background with S ≡ (3/4)(ρ DM /ρ γ ) as the scaling of the counter term in the momentum and ρ DM is the dark matter energy density. Analogous to the EM interaction, depends on the new interaction rateμ ≡ aσ DM−γ cn DM . Here σ DM−γ is the elastic scattering cross-section between DM and photons while n DM = ρ DM /m DM is the DM number density. For the DM-neutrino interactions, similar modifications can be added. In [34] an implementation of these modified Euler equation for the CLASS Boltzmann solver was presented. We are using this work to calculate the linear evolution of the Universe up to the point (in this work at redshift z = 127) where we switch to simulations to also cover the full non-linear evolution and resulting structure formation accurately (for more details see also [47]).

Simulations
For this work, we calculate the non-linear evolution of the matter distribution using a suite of cosmological DM-only simulations. This includes both simulations of single-resolution periodic volumes of 100 Mpc as well as zoom-in simulations which focus on representative sub-volumes to improve the maximum resolution for a subset of the obtained DM structure samples. We perform these simulations with the parallel Tree-Particle Mesh N-body code, P-Gadget3 [48] for both a standard, collision-less CDM and a γCDM model with a cross section σ/σ Th = 2 × 10 −9 (m DM /GeV). This value is (roughly) the upper bound obtained in previous works from satellite number counts of Milky-Way-size halos [30,49]. In [50] a more conservative constraint is claimed using measurements of the ionization history of the Universe at several redshifts, results from N-body simulations and recent estimates of the number of Milky Way satellite galaxies. However, the approach implemented can generate large uncertainties since the presence of low-mass subhalos in galactic halos which simulations can not resolve and extrapolations are necessary to obtain the results. Note that whereas larger cross sections would erase most of the observed substructure, smaller cross sections would imply results in between CDM and IDM. The simulations begin at a redshift of z = 127 (the DM-radiation interaction rate is negligible at all times afterwards). For the initial conditions we use the same cosmology (WMAP7), random phases and second-order LPT method [51] as the APOSTLE project [52] and our previous studies of the impact of IDM on galactic substructures [47]. After having performed the full-volume run for both standard CDM and γCDM with a particle mass m Part = 1.96 × 10 8 M /h and a comoving softening length l soft = 2.7 kpc, we identify the DM structures within using the Rockstar halo finder [53]. All halo properties are determined for spherically overdense regions with a density of 200 times the critical density. With these results, a cubic sub-volume is chosen at z = 0 with a side length of 14 Mpc/h that reproduces the overall halo mass function on the mass scales covered by it. A 1 Mpc wide margin is added and the resulting volume traced back to the initial redshift. We checked that the sub-volume thus constructed is still convex in these Lagrangian coordinates. This ensures that the progenitors of the structures within the targeted region evolve well within the high-res region, when the resulting volume is re-run using a zooming technique [54] with m Part = 4.85 × 10 5 M /h and l soft = 860 pc in the targeted region. Throughout this work, we use the term Box to refer to the full-volume simulation (100 Mpc) at z = 0 for each cosmology. The zoom re-simulations model four Local Groups (LGs hereinafter). We filter the results to pick only those halos that are well within the higher resolution region, namely inside a ∼ 2.1 Mpc/h radius at z = 0. This is done in order to avoid boundary affects, such has halos that consist partly of higher-mass particles, which are ignored here. The total number of halos and subhalos found in both Box and LGs simulations is given in Table 1, together with the most relevant parameters of these simulations.

Results
As mentioned, IDM exhibits a linear matter spectrum different to the one of CDM [26,[31][32][33][34][35][36]. The IDM matter power spectrum features a cut-off around a smooth scale of ∼ 100 kpc for the cross section that we are considering in this work ( σ/σ Th = 2 × 10 −9 (m DM /GeV)). Therefore a suppression of the number of halos below the scale of those hosting dwarf galaxies is expected (i.e. for halo masses below ∼ 10 10 M /h). In addition, such linear matter power spectrum impacts the structural halo properties, such as shape, spin, density profile and halo concentrations [30,46,47]. In this section, we show the results we found for halo and subhalo concentrations in our simulations, as well as subhalo abundances.

Halo concentrations
We consider two different definitions for the concentration parameter. The first and more standard definition is c ∆ ≡ R vir /r −2 , i.e. the ratio between the halo virial radius, R vir , and the radius r −2 at which the logarithmic slope of the DM density profile The other definition has the advantage of being independent of the adopted DM density profile and of the particular definition used for the virial radius since is a function of the peak circular velocity, V max , and the radius at which this velocity is attained, R max : [55][56][57] with H 0 the Hubble constant. Assuming an NFW profile [58,59], the relation between c V and c ∆ is given by [55] where ρ c is the critical density of the Universe at present, ∆ is the overdensity factor that defines the halos and r ∆ is its virial radius. Using our set of simulations, both Box and LGs for IDM and CDM models, we obtain the medians of c V and c ∆ . The latter was found by applying the c V -c ∆ relation of Eq. (5) to the c V (V max ) values found for every halo in the simulations. We adopt ∆ = 200 as the value for the overdensity to define the halos. For Box, we applied a restriction on halo maximum circular velocity such that only halos with V max > 60 km/s are included; in the case of the LG data set this restriction is set at V max > 10 km/s. Both criteria are adopted in order to avoid resolution issues in the determination of c V at the smallest scales resolved by the simulations. We have grouped halos in bins of V max and have obtained the medians of c V . For both the LGs and Box simulations, similar bin sizes were chosen to cover the entire V max range, ∼ 10 km/s < V max < 10 3 km/s. For each cosmology, we consider 5 bins in LGs and 9 bins for Box simulations.
In Fig. 1   order to compare IDM and CDM subhalo concentrations one-to-one, we also include in Fig. 1 the corresponding CDM concentrations. First, it is worth noting that Fig. 1

shows an excellent agreement between the concentrations values found in both Box and
LGs at the scale where the simulations overlap. Also, as expected, both IDM and CDM yield similar results at large halo masses, while we derive significantly lower median concentration values below halo masses ∼ 10 11 M /h in the case of IDM compared to CDM. Interestingly, this decrease of concentration values is similar to that found in WDM simulations, an effect that has been explained as being due to the delayed formation time of low-mass halos [60]. In addition, similar analysis for c 200 was performed in [30] and [46] where also the dependence with redshift was presented. Our results are in good agreement with such previous ones at z = 0. As we explained above, such results for the concentration-mass relation, c 200 (M 200 ), were obtained from c V (V max ) (see Eq. 5). In this way we double check previous results for IDM halo concentrations where a NFW profile was assumed. At late times, interacting DM models become (effectively) non-collisional for the cross section studied here, in the same way that the free-streaming in WDM models becomes negligible at low redshifts. Therefore, the observed lower IDM concentration values at small halo masses also originate from the later collapse of DM halos in these models.

Subhalo concentrations
The same analysis in V max and subhalo mass, m 200 , bins was performed for c V and c 200 subhalo concentrations respectively. In this case for Box, 8 bins are considered to cover de V max range and 5 for m 200 . We applied a restriction on subhalo maximum circular velocity such that only subhalos with V max > 60 km/s are included; for the LGs, this restriction is set at V max > 10 km/s considering just 3 bins for both V max and m 200 in order to obtain the median concentration values with a good subhalo statistics. From the results of Box and LGs simulations together, the V max range covered is 10 < V max < 500 km/s in each cosmology.
In the left panel of Fig. 2, we depict median c V (V max ) values and corresponding 1σ errors as found in Box (blue) and the four LGs (red). The right panel shows the results for c 200 (m 200 ). As in Fig. 1, we also include the corresponding CDM concentrations. As it can be seen, the medians of c V (c 200 ) in both cosmologies are similar for V max > 60 km/s (m 200 > 2 × 10 9 M /h), while there is a significant departure between them at lower V max (m 200 ) values. Unfortunately, the simulations have a limited mass resolution and subhalo statistics in that range, which translates into large 1σ errors and, as a consequence, our results are not conclusive. Yet, they provide a consistent picture of the subhalos' concentration behaviour at small V max (m 200 ) values, IDM subhalos exhibiting lower concentrations than CDM subhalos in the mentioned V max < 60 km/s range. Assuming a CDM framework, previous works have shown that the subhalo concentration depends not only on the mass of the subhalo but also on the distance to the center of its host halo [45,56,61]. In order to know if the same behaviour is found for IDM subhalos, Fig. 3 depicts, for the LGs, the medians and 1σ errors of c V (left panel) and c 200 (right panel) as a function of the distance from the host halo center in units of R 200 . As before, we also include in the figure our results for the CDM case. Median IDM subhalo concentration increases towards the center of the host halo more significantly than in the CDM case. Yet, for each considered radial bin, IDM concentrations are significantly and consistently lower than CDM ones. Again, large error bars prevent us from extracting firm conclusions and, thus, we will not propose any parametric fits to the data in this paper. The latter will be left for future work instead, when higher resolution IDM simulations might become available. However, this is an interesting quantitative result that points to a significantly different distribution of subhalo concentrations inside the host halo in the IDM scenario compared to CDM.
In the standard CDM cosmological framework, it is well established from simulations that subhalos are more concentrated than field halos of the same mass [9,41,45,56,[62][63][64][65][66]. It might not be the case in the IDM model, indeed the mean subhalo concentration values (see Fig. 2) fall within the values of halos concentrations studied in previous works for CDM. However, from Fig. 1 we see that the IDM halos exhibit lower concentrations compared with the halo concentrations in CDM of the same mass and then differences are expected between the concentrations of subhalos and their hosts in the interacting models. In Fig. 4 we shape such differences between halos and subhalos in the IDM scenario by comparing their median c V (c 200 ) values and 1σ errors as a function of V max (m 200 ) as found in our set of simulations. Analogously to what occurs in CDM, we obtained that also in IDM models subhalos with mass m 200 < 10 11 M /h tend to be more concentrated than their host halos. As in previous cases above, a more quantitative statement about the observed trend is nevertheless not possible for the moment, given the relatively large uncertainties involved in our study.

Subhalo abundances
As mentioned, DM interactions lead to a matter power spectrum different from the one in CDM. This matter power spectrum features a cut-off around a smooth scale of ∼ 100 kpc, and therefore a suppression of the number of halos in the lower mass range. The impact of such IDM initial matter power spectrum on the abundance of halos was studied in [30,46], where a comparison with the standard CDM result was also presented. A suppression of the number of low-mass halos with masses below M 200 ∼ 10 11 h −1 M was found, which became particularly significant at the smallest considered halo masses. In this section, we will complement these previous studies by using our set of IDM simulations to obtain the first results for subhalo abundances. We will do so in a broad subhalo mass range, i.e., [2 × 10 6 , 10 12 ] M /h.
In Fig. 5, we show the cumulative number of subhalos, N(> m 200 ), as a function of subhalo mass, m 200 , for both IDM and CDM scenarios and for both Box and LGs. Then, we consider all subhalos residing in halos with M h > 3 × 10 13 M /h for Box, and 3 × 10 11 M /h < M h < 1.4 × 10 12 M /h for LGs. These ranges allow us to have more than 30 subhalos per host in both cosmologies and both simulation sets. For each halo, we calculate the cumulative number of subhalos by adopting 100 subhalo mass bins and by finding the mean for each subhalo mass bin over all the main halos in the corresponding simulation. In the same Fig. 5, we also show in solid lines the result of fitting the data with the following parametric expression: This fitting function follows previous works that calculated the cumulative subhalo mass function from N-body cosmological simulations, and where the subhalo mass function was found to obey a power law dN/dm ∝ m −α 200 . [67]. Both the normalization factor, β, and the slopes γ = −α + 1, will depend on the adopted cosmological model. In Tab. 2, we report the best-fit values we found in our simulations for α and β, both for CDM and IDM scenarios. As it can be seen in Fig. 5 and in Tab. 2, in the case of the LGs the normalization of the cumulative subhalo mass function in the IDM case is significantly lower than that of CDM subhalos. More precisely, we find that mean N(> m 200 ) values for IDM subhalos are almost a factor ∼ 10 lower than those of CDM for subhalos in the range 10 7 M /h < m 200 < 10 8 M /h, this factor decreasing towards large subhalo masses. In Box, which covers comparatively larger halo masses, the differences among the two considered cosmologies are not statistically significant anymore. Indeed, all these results are as expected. As discussed above, the particular differences between the IDM and CDM initial matter power spectra lead to a suppression of smaller structures in the former case with respect to the latter, an effect that must become more evident in the LGs compared to Box, as the former simulations resolve smaller subhalo masses. Finally, we also studied the radial dependence of the number of subhalos in the IDM case, and compared it to the more standard CDM subhalo radial distribution. We did so only for the LGs, since high resolution simulations are necessary to perform this kind of analysis. Indeed, we checked that the statistics in the Box simulation is not sufficient to properly perform the work. LGs (red) simulations at z = 0. We also show the corresponding fits using Eq. 7 with the best-fit parameters reported in Tab. 2 (solid colored lines).
Number Density can be seen, the radial number density of IDM subhalos increases towards the center of the host halo as in the CDM case but is significantly lower than the latter at all host radii.

Discussion and conclusions
We have investigated DM subhalo properties in models where the linear matter power spectrum is suppressed at small scales due to DM interactions with radiation (photons or neutrinos). We do so by making use of N-body cosmological simulations, which are known to be a crucial tool to study the properties of DM structures. More precisely, we use data from our own set of simulations, described in Sec. 3. The runs are performed in both the standard CDM paradigm and in the IDM scenario, where the latter assumes interactions of DM with photons. 1 This allows us to compare DM halo and subhalo properties as found in both cosmologies. Since the main impact of the DM-photon interactions on structure formation occurs mainly at small scales, we use data not only from a large simulation box (100 Mpc) but also high-resolution zoom-in simulations of four Local Groups.
First, in sections 4.1 and 4.2, we studied, respectively, halo and subhalo concentrations as a function of halo/subhalo mass (and, alternatively, V max ). Both for halos and subhalos we observed a significant reduction of the concentrations in the lower mass range (or, alternatively, small V max values). Our result for halos confirm the findings of previous works, e.g. [30,46], while this is the first time that the concentration of IDM subhalos was studied. This decrease of concentration values is expected and originates from the later collapse of low-mass DM halos and subhalos in IDM cosmologies, similarly to that observed in WDM simulations [60].
Also in section 4.2, we studied subhalo concentrations as a function of the subhalo distance to the host halo center. As in the CDM framework, we found that the median subhalo concentration values increase towards the innermost regions of the host for subhalos of the same mass. Yet, we obtained significantly lower median concentrations in the IDM case with respect to CDM at all radii (see Fig. 3). Limitations in the number of subhalos prevent us from quantifying this effect more in detail, thus it seems robust in clearly present in our data. New N-body cosmological simulations with improved resolution will be needed in order to perform a more exhaustive analysis in this direction.
In addition, when comparing our results for IDM halos and subhalos of the same mass, we conclude that in these IDM models the subhalos are more concentrated than field halos (see Fig. 4), similarly to what found for CDM, e.g. [45].
Finally, we also presented in section 4.3 our results for subhalos abundances as a function of distance to host halo center and subhalo mass. Our results are in agreement with expectations for IDM models, namely we find a significantly smaller number of subhalos in IDM with respect to that observed in our CDM simulations. But not only the normalization of the cumulative subhalo mass function decreases (up to a factor ∼10 at the smallest resolved subhalo scales); also its slope is substantially lower in IDM (α = 1.7 versus α = 1.87 for CDM in the approximated range 10 7 M /h < m 200 < 10 9 M /h; see Fig. 5 and Tab. 2). As expected from theory, these differences among both cosmologies are not observed in the larger Box simulation. The radial distribution of subhalos within host halos exhibit a similar trend: there are fewer subhalos in IDM compared to CDM. Yet, we do not find appreciable differences in behaviour, i.e., the functional form of both radial distributions are similar.
In addition to the obvious interest for structure formation and study of halo and subhalo properties, we note that our work has a direct application on studies aimed at the indirect detection of DM, namely, the detection of the annihilation or decay products of DM particles. For instance, the extragalactic γ-ray and neutrino emission due DM annihilations depends mainly of the DM halos and subhalo properties (see e.g., [41,42,44,45]). Another example is the so-called subhalo boost: subhalos are expected to boost the DM signal of their host halos significantly, e.g. [43,45]. This subhalo boost 1 We do not include the case of DM-neutrino interactions, yet the results are expected to be similar to those presented in this work; see discussions e.g. in [30,47].
is very sensitive to the details of both subhalo concentration and subhalo abundance. Overall, from our results we conclude that the role of halo substructure in DM searches will be less important in IDM scenarios than in CDM, given the fact that both the subhalo concentrations and abundances are lower in the former compared to the latter. Yet, it will not be not negligible, as we also find in our IDM simulations larger concentrations for subhalos with respect to field halos of the same mass. Although this work represents an important step on addressing this and related issues, a quantitative study of the precise role of IDM subhalos for DM searches is left for future work: the IDM cosmological model mainly impacts low mass structures, thus it will be necessary to have higher resolution simulations than those used in this work in order to do so. Likewise, for a full analysis of IDM halo and subhalo properties it will be also necessary to run IDM simulations adopting other values of the cross section of DM interactions.