Halo Substructure Boosts to the Signatures of Dark Matter Annihilation
Shin'ichiro Ando Tomoaki Ishiyama Nagisa Hiroshima




§ INTRODUCTION


One of the most popular candidates for dark matter is weakly interacting massive particles (WIMPs) <cit.>.
They are motivated by beyond-the-standard-model physics such as supersymmetry <cit.> or universal extra-dimensions <cit.>, although the non-discovery of new physics at the TeV scale with the Large Hadron Collider puts these models to serious test <cit.>.  In addition, WIMPs can naturally
explain the relic dark matter density with thermal freezeout mechanisms, where the WIMPs following the weak-scale physics were in chemical equilibrium until freezeout—when the expansion of the Universe became faster than the annihilation rate <cit.>.
Since dark matter is often the lightest particle in an extended sector, it can self-annihilate only into the standard-model particles, which end up producing gamma rays, charged cosmic rays, and neutrinos. 
Indirect detection of dark matter annihilation is therefore a direct test of the thermal freezeout of WIMPs.

WIMPs are also a subcategory of cold dark matter (CDM), where they were nonrelativistic when structure formation started. In the CDM framework, it is known that the structures form hierarchically, from smaller to larger ones.
These virialized structures are referred to as halos and they are nearly spherically symmetric.
Typical size of the smallest structure is highly model dependent.
In the case of the supersymmetric neutralino that is one of the most popular WIMP candidates, the smallest halos tend to be of the Earth mass, 10^-6 M_⊙ but with very large range of possible values of ∼10^-12–10^-3 M_⊙ <cit.>.
Smaller halos collapse at higher redshifts when the Universe was denser, and hence they are of higher density.
A larger dark matter halo today contains lots of substructures (or subhalos) of all mass scales, which can go down to the Earth masses or even smaller and hence denser.

Since the annihilation rate depends on the dark matter density squared (and ⟨ρ^2 ⟩≥⟨ρ⟩^2), the presence of the subhalos will boost the gamma-ray signatures from dark matter annihilation.
This subhalo boost of dark matter annihilation, in relation with the smallest-scale subhalos, has been a topic of interest for very many years <cit.>.

The main difficulty is the fact that subhalos of all the scales ranging from the Earth mass (or even smaller) to larger masses (a significant fraction of their host's mass) give a substantial contribution to the annihilation rate.
Covering this very large dynamic range is challenging even with the state-of-the-art numerical simulations.
In simulations of Milky-Way-size halos (10^12M_⊙) <cit.>, one can resolve only down to 10^4–10^5 M_⊙, and there still remains more than ten orders of magnitude to reach.

We will review recent progress on the subhalo contribution to dark matter annihilation.
(See also Reference <cit.> for a review on generic processes that subhalos undergo.)
We first discuss approaches using the numerical N-body simulations and estimate of the annihilation boost factor by adoping the results and extrapolation down to very-small-mass ranges.
To complement the approach based on simulations, we then review an analytical approach.
In the CDM framework, fraction of halos that collapse is described with the Press-Schechter formalism <cit.> based on spherical or ellipsoidal collapse models.
This has been further extended to accommodate collapsed regions within larger halos (excursion set or extended Press-Schechter formalism <cit.>), which can be applied to address statistics of halo substructure.
More recent literature suggests that the annihilation boost factor, defined as the luminosity due to subhalos divided by the host luminosity, is modest, ranging from order of unity to a few tens for galaxy-size halos <cit.>.
This relatively mild amount of the annihilation boost makes the prospect of indirect dark matter searches less promising compared with earlier more optimistic predictions <cit.>.
We note that our focus is mainly on subhalo boost factors in extragalactic halos. For the subhalo boosts in the Galactic halo, on the other hand, we need to assess the spatial distribution of the subhalos too. The N-body simulations described in Section <ref> can address this issue but again are subject to resolution issues as well as the baryonic effect. See, for example, Reference <cit.> for an alternative approach adopting analytical prescription.

This review is organized as follows.
In Section <ref>, we introduce basic concepts of density profiles, mass functions, and the annihilation boost factors of the subhalos, starting with simple formulations.
Here we make some simplifying assumptions, which are to be addressed in later sections.
In Section <ref>, we summarize the progress from the numerical simulations for the subhalos and the annihilation boost factors.
Section <ref> presents more recent approaches based on realistic formulation than Section <ref>. 
In Section <ref>, we first show new analytic models that predict the subhalo mass functions well in agreement with the results from the numerical simulations for various ranges of the host masses and redshifts, and that the annihilation boost factors are on the order of unity even for cluster-size halos.
Then, we summarize other semi-analytic approaches for computing the annihilation boost factors, based on self-similarity (Section <ref>) and universal phase-space clustering (Section <ref>) of the subhalos.
We conclude the review in Section <ref>.
Finally, for convenience, we summarize fitting functions for the subhalo mass functions, and annihilation boost factors that can be applicable to nearly arbitrary masses and redshifts in Appendix <ref>.





§ FORMULATION



In this section, we introduce several important quantities such as density profiles, subhalo mass function, and the annihilation boost factors.
This section is based on a simplified analytic model, which in several aspects are unrealistic but sets the basis for the latter discussions according to numerical simulations (Section <ref>) and more sophisticated semi-analytical models (Section <ref>).



 §.§ Subhalo Boost Factor


The rate of dark matter annihilation is proportional to dark matter density squared, ρ_χ^2, where χ represents the dark matter particle.
In the presence of substructure, ρ_χ is divided into two terms:

    ρ_χ( x) = ρ_χ, sm( x)+ρ_χ, sh ( x),

where ρ_ sm and ρ_ sh represent smooth and subhalo components, respectively.
(In the following, we omit the subscript χ.)
The volume average of the density squared in a host halo characterized by its virial mass M and redshift z, which is the relevant quantity for the indirect dark matter researches, is therefore written as

    ⟨ρ^2( x)⟩_M,z=⟨ρ_ sm^2( x)⟩_M,z + ⟨ρ_ sh^2( x)⟩_M,z + 2⟨ρ_ sm( x)ρ_ sh( x)⟩_M,z.

We assume that the smooth component ρ_ sm is characterized by the following Navarro-Frenk-White (NFW) profile <cit.>:

    ρ_ NFW(r) = ρ_s/(r/r_s)(1+r/r_s)^2,

where ρ_s is a characteristic density and r_s is a scale radius.
These parameters, ρ_s and r_s, are evaluated such that the volume integral of ρ_ NFW yields the total halo mass M, and thus we have ρ_ sm(r) = (1-f_ sh) ρ_ NFW(r), where f_ sh is defined as the mass fraction in the subhalos.
The first term is then simply

    ⟨ρ_ sm^2( x)⟩_M,z = 1/V∫ d^3x ρ_ sm^2(r) =  4π(1-f_ sh)^2/3Vρ_s^2 r_s^3[1-1/(1+c_ vir)^3],

where V = 4π r_ vir^3/3 is the volume of the host out to its virial radius r_ vir, c_ vir≡ r_ vir/r_s is the concentration parameter.
The parameters characterizing the host profile – ρ_s, r_s, and c_ vir – are all functions of M and z.[We note, however, that the concentration c_ vir has a scatter, which is often characterized by a log-normal distribution, whose mean c̅_ vir is the function of M and z. We will include this in the latter sections.]

Next, we evaluate the second term of Equation (<ref>), ⟨ρ_ sh^2( x)⟩.
We characterize each subhalo i with the location of its center x_i and mass m_i.
Density due to all the subhalos at a coordinate x is written as a sum of the density profile around the seed of each subhalo, that is,

    ρ_ sh(x) = ∫ dm' ∫ d^3x' ∑_i δ_D(m'-m_i) δ_D^3(x'-x_i)m'u_ sh(x-x'|m'),

where δ_D^N is the N-dimensional Dirac delta function, and u_ sh(r|m) defines the density profile of the subhalo with mass m and is normalized to one after the volume integral.[For the sake of simplicity for analytic expressions, we assume that the suhbalo mass is the only parameter characterizing its density profile. One can introduce many more parameters to make the model more realistic.]
We define the ensemble average of the product of these delta functions as

    dn_ sh( x,m)/dm = ⟨∑_i δ_D(m-m_i) δ_D^3(x-x_i)⟩,

its volume integral over the host halo as

    dN_ sh/dm = ∫ d^3x dn_ sh( x,m)/dm,

and call both dn_ sh/dm and dN_ sh/dm the subhalo mass function.

We also obtain the mass fraction in the subhalos as

    f_ sh(M,z) = 1/M∫ dm mdN_ sh/dm.



By multiplying Equation (<ref>) by itself and taking both the ensemble and the volume averages, we have

    ⟨ρ_ sh^2( x)⟩_M,z   ≡   1/V∫ d^3 x ⟨ρ_ sh^2( x)⟩
       =   1/V∫ d^3 x∫ dm'∫ d^3 x' ∫ dm”∫ d^3 x” 
    m'u_ sh(x-x'|m')m”u_ sh(x-x”|m”)
    
          ×⟨∑_iδ_D(m'-m_i)δ_D^3( x'- x_i)
    ∑_jδ_D(m”-m_j)δ_D^3( x”- x_j)⟩
       =   1/V∫ d^3 x ∫ dm' ∫ d^3 x' dn_ sh( x',m')/dm' m'^2 u_ sh^2( x- x'|m') 
       =   4π/3V∫ dm dN_ sh/dmρ_s, sh^2r_s, sh^3[1-1/(1+c_t, sh)^3],

where at the last equality we adopted the NFW function for the subhalo density profile m u_ sh(r|m) with the scale radius r_s, sh, characteristic density ρ_s, sh, and tidal truncation radius r_t, sh≡ c_t, sh r_s, sh beyond which the subhalo density abruptly decreases to zero.
At the third equality of Equation (<ref>), we ignored the term arising from j ≠ i as we evaluate the quantity at one point x and assume that subhalos do not overlap. We note, however, that such a term becomes relevant for obtaining the two-point correlation function, or the power spectrum; see References <cit.> for more details.

We define the subhalo boost factor as the ratio of the total luminosity from dark matter annihilation in the subhalos and that from the smooth component in the case that there is no substructure.
By comparing Equations (<ref>) and (<ref>), and remembering that the luminosity is proportional to the volume integral of the density squared, the boost factor is simply written as

    B_ sh(M,z) = 1/L_ host, 0(M,z)∫ dmdN_ sh(m|M,z)/dmL_ sh(m),

where the subscript 0 shows that this is a quantity in the case of no subhalo contributions.
Equation (<ref>) is also valid for any other spherically symmetric density profiles than the NFW.

Finally, we evaluate the last cross-correlation term in Equation (<ref>). See also References <cit.>.
Following a similar procedure as in Equation (<ref>), we have

    2⟨ρ_ sm( x)ρ_ sh( x)⟩_M,z   =   2/V∫ d^3x⟨ρ_ sm( x)ρ_ sh( x)⟩
    
    
    
       =   2/V∫ d^3x ρ_ sm( x)∫ dm' ∫ d^3x' m'u_ sh( x- x'|m')
    
    dn_ sh( x',m')/dm'
       ≈   2/V∫ d^3xρ_ sm( x)∫ dm mdn_ sh( x,m)/dm,

where in the last equality, we first used the fact that the subhalo density profile is much more sharply peaked than their spatial distribution, and take dn_ sh/dm out of x' integration adopting x' ≈ x as its spatial variable.
Second, we performed volume integral for u( x- x'|m) over x' variable, which simply returns one, to reach the last expression of Equation (<ref>).
Then we assume that the spatial distribution of the subhalos is independent of their masses:

    dn_ sh( x,m)/dm = P_ sh( x) dN_ sh(m)/dm,
 
where P_ m sh( x) d^3x represents the probability of finding a subhalo in a volume element d^3x around x.
With this and Equation (<ref>), we have

    2⟨ρ_ sm( x)ρ_ sh( x)⟩_M,z = 2f_ shM/V∫ d^3xρ_ sm( x)P_ sh( x).

For simplicity, we assume that the subhalos are distributed following the smooth NFW component.
In this case, we have ρ_ sm( x) = (1-f_ sh)MP_ sh( x), and 
    2⟨ρ_ sm( x)ρ_ sh( x)⟩_M,z = 2f_ sh/1-f_ sh⟨ρ_ sm^2( x)⟩= 8π f_ sh(1-f_ sh)/3Vρ_s^2r_s^3[1-1/(1+c_ vir)^3].




The luminosity from the smooth component in the presence of the subhalos (L_ sm) is related to the host luminosity in the subhalos' absence via L_ sm = (1-f_ sh)^2L_ host,0, because the density in the smooth component gets depleted by a factor of 1-f_ sh, if there are subhalos.
Thus, the total luminosity from both the smooth component and the subhalos are given by

    L_ total   =    L_ sm + L_ sh + L_ cross
       =   [(1-f_ sh)^2+B_ sh+2f_ sh(1-f_ sh)]L_ host,0
       =   (1-f_ sh^2+B_ sh)L_ host,0

Often, L_ total/L_ host,0 is also referred to as the subhalo boost factor in the literature.
Note, however, that we have not included the effect of sub-subhalos (and beyond) yet in this formalism.
In order to accommodate it, in the right-hand side of Equation (<ref>), we need to include the sub-subhalo boost to the subhalo luminosity L_ sh. Thus, we replace L_ sh with (1-f_ ssh^2+B_ ssh)L_ sh, where the subscript “ssh” represents the contribution from the sub-subhalos. If the subhalo mass fraction f_ sh and the boost factor B_ sh depend only on the host mass, then one can assume f_ ssh(m) = f_ sh(m) and B_ ssh(m) = B_ sh(m), and repeat the calculations in an iterative manner. See, however, Section <ref> for a more realistic treatment.





 §.§ Characterization of Dark Matter Halos



We shall discuss the density profile of dark matter halos that are characterized by the virial radius r_ vir, the scale radius r_s, and the characteristic density ρ_s.
The halo is virialized when a mean density within a region reaches some critical value times the critical density of the Universe at that time: Δ_ vir(z)ρ_c(z), where ρ_c(z) = ρ_c,0 [Ω_m(1+z)^3+Ω_Λ], ρ_c,0 = 3 H_0^2/(8π G) is the present critical density, H_0 = 100 h km s^-1 Mpc^-1 is the Hubble constant,  Ω_m and Ω_Λ are the density parameters for matter and the cosmological constant, respectively.
In CDM cosmology with the cosmological constant, this critical value is given as <cit.>

    Δ_ vir(z) = 18π^2 + 82 d(z) - 39 d^2(z),

where d(z) = Ω_m(1+z)^3/[Ω_m(1+z)^3+Ω_Λ] - 1.
Given the virial mass M and the redshift z of the halo of interest, r_ vir is therefore obtained by solving

    M = 4π/3Δ_ vir(z)ρ_c(z) r_ vir^3.

Alternatively, one can define M_200 and r_200 via

    M_200 = 4π/3200ρ_c(z) r_200^3.

M_200 is often adopted to define halo masses in N-body simulations.

The concentration parameter c_ vir≡ r_ vir/r_s (or c_200≡ r_200/r_s) has been studied with numerical simulations and found to be a function of M and z.
It follows a log-normal distribution with the mean of c̅_ vir(M,z) (e.g., Reference <cit.>) and the standard deviation of σ_log c≈ 0.13 <cit.>.
The mean c̅_ vir has been calibrated at both large (galaxies, clusters) and very small (of Earth-mass size) halos, and found to decreases as a function of M and z. 
Once c_ vir is drawn from the distribution, it is used to obtain r_s = r_ vir/c_ vir.
Finally, ρ_s is obtained through the condition of having mass M within r_ vir:

    M    =   ∫_0^r_ vir dr 4π r^2 ρ(r) = 4πρ_s r_s^3 f(c_ vir),
    
        f(x)    =   ln(1+x)-x/1+x,

where the second equality of Equation (<ref>) holds in the case of the NFW profile.

In the case of the subhalos, the procedures above cannot be adopted.
This is because they are subject to tidal effects from the host, which strip masses away from the subhalos.

However, the regions well inside the scale radius r_s—because of strong self-gravity—is resilient against the tidal force and hence the annihilation rate hardly changes.
These tidal processes, therefore, make the subhalos more concentrated and hence effectively brighter compared with the field halos of the same mass.
In many analytical studies in the literature <cit.>, however, the effect of tidal stripping was ignored and the concentration-mass relation of the field halos was adopted, which resulted in underestimate of the annihilation boost factor.
This has been pointed out by Reference <cit.> and will be discussed in Section <ref> (see also References <cit.>).







§ ESTIMATES OF ANNIHILATION BOOST WITH NUMERICAL SIMULATIONS




In order to assess the annihilation boosts, one has to have reasonably good ideas on the density profiles ρ(r), the concentration-mass relation,[The concentration-mass relation is defined as the average concentration parameter as a function of halo mass.] and the subhalo mass function. 

Cosmological N-body simulations have been a powerful tool for probing all of them because
once a halo collapses from initial density fluctuations, it evolves
under a strongly nonlinear environment.  They have indeed demonstrated
that there are a large amount of surviving subhalos (see Figure <ref>) in halos and halos have cuspy density profiles.










 §.§ Subhalo Abundance


Cosmological N-body simulations predict that there are many
surviving subhalos in host halos as a consequence of hierarchical structure
formation.  <cit.> and <cit.> performed
high-resolution cosmological N-body simulations for the formation
and evolution of galaxy-scale halos.  They demonstrated that too many
subhalos existed in simulated halos in comparison with the number of
observed dwarf galaxies in the Local Group.  This discrepancy is known
as the “missing satellite problem” that has been investigated by a number of follow-up simulation studies (e.g., References <cit.>).






Even though it triggered many studies attempting to reduce small-scale structures by imposing other non-CDM candidates such as warm dark matter <cit.> and self-interacting dark matter <cit.>, it is also possible to solve it with standard baryonic physics including early reionization <cit.> and self-regulation of star formation in low-mass halos <cit.>.  Hence, it
is no longer regarded as a serious problem of the CDM model.

These studies suggest a large number of “dark satellites”
exist in halos, which do not contain optically visible components such as
gases and stars.  The population of dark satellites is more abundant
in host halos than visible satellite galaxies and could enhance
annihilation boosts significantly.  To estimate subhalo boosts to
annihilation signals accurately, understanding abundance of
subhalos as well as their density structure is crucial.

A number of studies have calculated the subhalo mass function in halos using cosmological N-body simulations (e.g., References <cit.>), indicating 
that it obeys a power law

    dN_ sh/dm∝ m^-α,

where the slope -α ranges from -2 to -1.8, although no
consensus has yet emerged. There is also a large halo-to-halo scatter
for the subhalo abundances <cit.>. 
The subhalo abundance at a fixed mass halo depends on their accretion
history. Namely, it increases with the mass of halo and decreases with the 
halo concentration (e.g., References <cit.>).

Due to the limitation of currently available computational resources,
simulations cannot resolve the full hierarchy of subhalos from the
smallest to the most massive scales, which ranges more than twenty orders of magnitude in the mass. Even in the highest resolution
simulations for galaxy-scale halos, the smallest resolved subhalo mass
is around ∼10^5 M_⊙ <cit.>, which is still more than ten orders of magnitude more massive than that of the cutoff scale.  



To study subhalo boosts to annihilation signals, a single-power-law subhalo mass function, Equation (<ref>), is traditionally extrapolated beyond the resolution. 

Another approach is to use some analytical
models (e.g., References <cit.>), which can shed light on the resolution issue. <cit.> developed  a model of the subhalo evolution calibrated with cosmological N-body simulations and found that the power-law index of the subhalo mass function is in
a rather narrow range between -2 and -1.8 with a vast range of subhalo
mass from z = 0 to 5.  
This picture is more or less consistent with the assumption of the subhalo mass function of the single power law.
More details on the analytic approach are discussed in the following section.

Note that the annihilation boost factors strongly depend on the underlying
subhalo mass function <cit.>.
Assuming that ∼10% of the halo mass is within subhalos, the difference of the boost factors in the Milky-Way-size halo could be as large as a factor of ten between the slope of -α=-2 and -1.9 <cit.>.  


More extensive simulations are needed to obtain the subhalo mass function in wide mass ranges and also to compare with analytic models <cit.>.





 §.§ Density Profile of Dark Matter Halos


By the end of the 1980's, it was already known in both
analytic <cit.> and numerical <cit.>
studies that the density profiles were described by power-law functions. Reference
<cit.> studied the density profiles of dark matter
halos using cosmological N-body simulations and argued that the
profiles were well described by a Hernquist
model <cit.>. 



Navarro et al. <cit.> simulated 
the structures of CDM halos systematically
with masses in the range of galaxy to rich cluster size.
They claimed that the radial density profile ρ(r) could be described by a simple universal profile, Equation (<ref>), the so-called NFW profile.
They also claimed that the shape of the profile was universal, independent of cosmological parameters, the primordial power spectrum, and the halo mass.  






Today, the NFW profile has been extensively used to model halos analytically for various purposes.

After the work of NFW, A number of subsequent studies (e.g., References <cit.>) performed simulations with better mass resolutions.
Whereas previous studies <cit.> used only ∼10,000 particles, they used ∼1,000,000 particles for a halo, and found that the slope was steeper than -1.  
In the original results of the NFW, the numerical two-body relaxation effects due to the small number of particles affected the structures of central regions and led to form a shallower cusp.  
Higher resolution simulations could resolve more inner structures of halos <cit.>.  
In most cases, the slope of density became shallower as the radius went inward.  
A different approach was adopted by  <cit.>, who used the triaxial model for describing the central structures.
<cit.> and <cit.> considered a more general profile,

    ρ(r) = 
    ρ_s/(r/r_s)^γ[ 1+(r/r_s)^η]^
    (β-γ)/η .

If β=3, γ=1, and η=1, the profile is the same as NFW.

More recent studies <cit.> archived one of the highest resolution dark matter only simulations for galaxy-size halos with mass resolution better than 10^4M_⊙.  
Their results are in agreement in that the density slope cannot be described by a single power law and the slope is around -1 at the radius ∼0.001r_200.  
Besides, <cit.> and <cit.> fitted the density using the Einasto profile <cit.> 

    ρ(r) = ρ_s exp{-2/α_E[(r/r_s)^α_E-1 ]},

where α_E is a free parameter.
Note that r_s and ρ_s are not the same parameters as those in Equation (<ref>).

Although we can obtain the density profile down to the radius ∼0.001r_200, the result does not converge to a single power law.
In addition, the physical origin of this flattening towards the center is not
understood at all.  However, the importance of understanding the
central structures is increasing.  In particular, if we would like to
detect signals from dark matter annihilation, the central structure of the dark matter halo is essential. 


The most important parameter to describe the halo profile is the concentration parameter.  
Assuming the universal NFW profile regardless of the halo mass, the concentration-mass relation gives the annihilation rate as a function of the halo mass.  
Combined with assumed subhalo mass functions, they enable to estimate the annihilation boost factor.  
The concentration-mass relation of halos has been widely investigated and a number of fitting functions has been suggested <cit.>.  
The concentration shows a weak dependence on the halo mass.  
The average concentration at fixed halo mass becomes smaller with increasing halo mass because the central density is tightly correlated with the cosmic density at the halo formation epoch, reflecting the hierarchical structure formation <cit.>.

Traditionally, the concentration-mass relation has been calibrated with cosmological N-body simulations for relatively massive halos (10^10M_⊙≲ M_ 200≲ 10^15 M_⊙).  
Because the mass dependence of the concentration is weak for these mass halos, it is found that a single power-law function, c ∝ M_200^-α_c, with slope α_c in the range of 0.08 to 0.13, 

gives reasonable fits <cit.>.
However, the dependence gradually becomes weaker toward less massive halos, and a clear flattening emerges <cit.>, ruling out single power-law concentration-mass relation 
for the full hierarchy of subhalos.


These fitting functions are valid for the NFW density profile.
More generally, the concentration can be defined independently
of the density profile and the subhalo mass as
(e.g.,  <cit.>)

    c_V = ρ̅ (<r_ max)/ρ_ c
    = 2 ( V_ max/H_0r_ max)^2,

where r_ max is the radius at which the circular velocity reaches its maximum value V_ max.
This definition is also used to estimate the
annihilation boost factor (e.g.,  <cit.>).



Even with the highest resolution simulations for galaxy-scale halos, the smallest resolved subhalo mass is around ∼10^5 M_⊙ <cit.>.
To estimate the annihilation boost from the full hierarchy of subhalos, we have to make some assumption of the concentration at unresolved scales, which has a significant impact on the result.  
One approach is extrapolating single power-law  fittings to the smallest scale beyond the mass range calibrated with simulations, although the literature including the above cautions the risk of such extrapolations.  
With such extrapolations, the concentration of the smallest halo can reach more than 100, substantially enhancing the annihilation boost. 
A number of studies have computed the concentration in such a manner and the resulting boost factor is a few hundreds for Milky-Way halos <cit.>, and ∼1000 for cluster-scale halos <cit.>. 

Another approach is adopting analytic models or fitting functions that can reproduce flattening of the concentration-mass relation 

(e.g., <cit.>).

In contrast to using the power-law extrapolation, the resulting boost factor is rather modest, three to a few tens <cit.> for Milky-Way halo, and less than ∼100 for cluster-scale halos <cit.>.

The density profile at fixed halo mass shows a significant halo-to-halo scatter <cit.>, possibly making a big impact on the annihilation signal.
Inferring from the cosmological Millennium simulation <cit.>, the effect of this non-universality on the annihilation flux is a factor of ∼3 <cit.>, which indicates that the uncertainty of the concentration-mass relation for low-mass halos has a more significant effect.

These discussions are based on the universal density profile and the concentration-mass relation for field halos. 

There is a concern that whether or not we can apply the universal NFW profile for the full hierarchy of halos and subhalos beyond the range that cosmological simulations have been able to tackle.
We discuss this issue in Section <ref>.
More importantly, we have to use the concentration-mass relation for subhalos, not field halos. 
We also discuss this issue in the following section.





 §.§ Density Profile of Dark Matter Subhalos



Density structures of subhalos are more challenging to be investigated
than field halos because it requires much higher resolutions.  Therefore, to evaluate the subhalo contribution to the
annihilation signals, the universal NFW profile and the
concentration-mass relation for field halos have been historically
assumed to be the same for subhalos as a first approximation, although the underlying assumption is not well studied.  Complex physical
mechanisms relevant to subhalos could change their original density
profiles, such as the tidal effect from host halos, the encounter with
other subhalos, and denser environment than the field.



Cosmological simulations have been suggesting that the density profile
of subhalos is cuspy in analogy with field halos. On the other hand,
the average concentration of subhalos tend to be higher than those of
field halos (e.g., <cit.>).  For example,
<cit.> showed that subhalos and halos tend to be more
concentrated in dense environments than in the field, and the scatter of
concentrations is larger.  This result was taken into account to
estimate the gamma-ray flux from dark matter annihilation (e.g., <cit.>).
<cit.> showed that outer regions of subhalos tend to be
tidally stripped by host halos, which gives higher concentrations.
These results suggest that both earlier formation of halos/subhalos in
dense environments and tidal effect are responsible for the increased
concentration.  <cit.> derived the concentration-mass
relation of subhalos in Milky-Way-size halos by analyzing high
resolution cosmological simulations <cit.> and showed that it depends on the location of subhalos
relative to host halos.  Subhalos have considerably large
concentrations near the center than at the edge of host halos.
<cit.> quantified the concentration of subhalos in
Milky-Way-size halos as a function of not only subhalo mass but
distance from host halo center, and found a factor 2–3 enhancement of
the boost factor compared to the estimation that relied on the
concentration-mass relation of field halos (see also Figure <ref>).





As shown in the literature in the above, higher concentrations of
subhalos than field halos could have a big impact on the annihilation
boost.  However, <cit.> argued that
subhalos even in state-of-the-art cosmological simulations 
suffer from excessive mass loss and artificial tidal disruption due to
inadequately large force softening (see also <cit.>).  If that is the case, it might be possible
that subhalos have larger concentrations than those ever considered.
These issues should be addressed by extremely high resolution
cosmological N-body simulations and analytic models
(e.g., <cit.>).






 §.§ Density Profile of Dark Matter Halos Near Cutoff Scales




In the CDM framework, smaller halos
collapse first, and then they merge into more massive halos.  
Since the smallest halos contain no subhalos, their central structures might entirely differ from that observed in more massive halos.
If the dark matter particle is the lightest supersymmetric particle such as the neutralino, the smallest halo mass is predicted to be around the Earth mass <cit.>.
Such halos are sometimes referred to as “microhalos.”



The density profiles of the microhalos have been investigated using cosmological N-body simulations <cit.>.  <cit.> simulated the formation of Earth-mass microhalos by means of cosmological N-body simulation.
They claimed that a single power law could describe the density profiles of microhalos, ρ(r) ∝ r^-γ, with a slope γ in the range of 1.5 to 2.  As a
consequence of such steep slope, most microhalos could not be completely destructed by the Galactic tide and encounters with stars, even in the Galactic center.  


<cit.> have performed N-body simulations with much higher resolution and showed that the density profile of microhalos had steeper cusps than the NFW profile.
The central density scales as ρ(r) ∝ r^-1.5, which is supported by follow-up cosmological simulations <cit.> and cold-collapse simulations <cit.>.
<cit.> has also shown that the cusp slope gradually becomes shallower with increasing halo mass.
Major merger of halos is responsible for the flattening, indicating that the process of violent relaxation plays a key role (see also <cit.>).  
Similar density structures are observed in recent simulations of ultracompact minihalos <cit.> and warm dark matter <cit.>.

The self-similar gravitational collapse models
(e.g., <cit.>)
can also give hints to understand the main physical origin of such steeper
cusps, because the smallest halos do not contain smaller density
fluctuations by definition and collapse from initially overdense
patches.


Such microhalos with steep cusps can cause a significant effect on indirect dark matter searches.
<cit.> argued that the central parts of microhalos could survive against the encounters with stars except in the very Galactic center.
The nearest microhalos could be observable via gamma rays from dark matter annihilation, with usually large proper motions of ∼0.2 deg yr^-1, which are, however, stringently constrained with the diffuse gamma-ray backgrond <cit.>.  Gravitational perturbations to the 


millisecond pulsars might be detectable with future observations by pulsar timing arrays <cit.>.
<cit.> have assumed a transition from the NFW to steeper cusps at scales corresponding to ∼100 times more massive than the cutoff and have found that such profiles can enhance moderately the annihilation boost of a Milky-Way-size halo by 5–12%.
They also have found that concentrations of microhalos are consistent with a toy model proposed by <cit.>.

<cit.> showed that the steeper inner cusps of halos in the smallest scale and near the cutoff scale could increase the annihilation rate of a Milky-Way-size halo by 12–67%, compared with estimates adopting the universal NFW profile and an empirical concentration-mass relation <cit.> 
(see Figure <ref>).
The value, however, depends strongly on the adopted subhalo mass function and concentration model.  
They have found that concentrations near the free-streaming scale show little dependence on the halo mass and corresponding conventional NFW concentrations are 60–70, consistent with the picture that the mass dependence is gradually becoming weaker toward less massive halos (e.g., <cit.>), ruling out a single power-law concentration-mass relation.





As shown in the literature above, steep density cusps of halos near the free-streaming scale have an impact on the annihilation
boost.
However, these studies rely on the density structure seen in field halos, not subhalos.
It is also important to quantify the structures of subhalos near the free-streaming scale by larger simulations.
Another concern is that the cutoff in the matter power spectrum should suppress the number of subhalos near the free-streaming scale, which should weaken the annihilation signal.
However, the shape of the mass function near the free-streaming scale is not understood well for the neutralino dark matter.
The structure of subhalos and the subhalo mass function near the free-streaming scale should be explored by larger volume cosmological N-body simulations.







§ SEMI-ANALYTIC APPROACHES







 §.§ Models Based on Structure Formation and Tidal Evolution




In an analytical approach in Section <ref>, the subhalo luminosity L_ sh is characterized with the mass of the subhalo and the redshift of interest; see, for example, Equation (<ref>).
The mass and redshift, however, are not the only quantities that fully characterize the subhalo properties.
Indeed, they depend on the accretion history and mass loss after they fall onto their host halo,
that is, two subhalos that have the identical mass could have formed with different masses and accreted at different redshifts, evolved down to z = 0 reaching the same mass.
<cit.> and <cit.> developed an analytical prescription to take these effects into account, which we follow in this section.

A subhalo is characterized with its mass and redshift when it accreted onto its host, (m_a, z_a).
The concentration parameter c_a is drawn from the log-normal distribution with mean c̅_a(m_a,z_a) <cit.> and σ_log c = 0.13 <cit.>.
Since the subhalo was a field halo when it just accreted, one can use the relations in Section <ref> to obtain r_s,a and ρ_s,a for the NFW profile.

After the accretion, the subhalos evolve by losing their mass through tidal forces.
The mass-loss rate is typically characterized by a dynamical timascale at the redshift z,

    τ_ dyn(z) =1.628h^-1  Gyr[Δ_ vir(z)/178]^1/2[H(z)/H_0]^-1,




as follows <cit.>:

    ṁ(z) = -Am(z)/τ_ dyn(z)[m(z)/M(z)]^ζ,

where H(z) = H_0 [Ω_m(1+z)^3+Ω_Λ]^1/2, m(z) and M(z) are the subhalo and host-halo masses at z, respectively.
Following <cit.>, <cit.> adopted simple Monte Carlo simulations to estimate ṁ based on the assumption that the subhalo loses all the masses beyond its tidal radius in one complete orbit at its peri-center passage.
While <cit.> found A = 0.81 and ζ = 0.04, <cit.> extended the mass and redshift ranges of applicability and found that these parameters are weakly dependent on both M and z:

    log A    =   {-0.0003log[M(z)/M_⊙]+0.02}z+0.011log[M(z)/M_⊙]-0.354, 
    ζ   =   {0.00012log[M(z)/M_⊙]+0.0033}z+0.0011log[M(z)/M_⊙]+0.026.

One can solve Equation (<ref>) to obtain the subhalo mass at a redshift of interest z, m(z), with a boundary condition of m(z_a) = m_a.
For the evolution of the host, M(z), <cit.> adopted a fitting formula given by <cit.>.

The subhalo density profile after accretion is also well described with the NFW profile with a sharp truncation at r_t:

    ρ(r) = {[ ρ_s r_s^3/[r(r+r_s)]^2,                     for                  r<r_t,;                      0,                     for                 r≥ r_t. ].

This is indeed a good approximation found in the simulations <cit.>.
In addition to r_t,  found that the internal structure changes.
If the inner profile is ∝ r^-1 just like NFW, the maximum circular velocty V_ max and its corresponding radius r_ max evolve as

    V_ max(z)/V_ max,a   =   2^0.4[m(z)/m_a]^0.3/[1+m(z)/m_a]^0.4, 
    r_ max(z)/r_ max,0   =   2^-0.3[m(z)/m_a]^0.4/[1+m(z)/m_a]^-0.3,

respectively.
After computing V_ max and r_ max at z, one can convert them to ρ_s and r_s through

    r_s    =   r_ max/2.163,
    ρ_s    =   4.625/4π G(V_ max/r_s)^2,

which are valid for the NFW profile.
Finally by solving the condition

    m(z) = ∫_0^r_tdr 4π r^2ρ(r) = 4πρ_s r_s^3 f(r_t/r_s),

the truncation radius r_t is obtained. 
<cit.> omitted subhalos with r_t < 0.77 r_s from the subsequent calculations assuming that they were tidally disrupted <cit.>. This criterion, however, might be a numerical artifact <cit.>. Either case, <cit.> checked that whether one
implements this condition or not did not have impact on the results of, for example, subhalo mass functions.

Thus, given (m_a, z_a, c_a), one can obtain all the subhalo parameters after the evolution, (m, r_s, ρ_s, r_t), in a deterministic manner.
The differential number of subhalos accreted onto a host with a mass m_a and at redshift z_a, d^2N_ sh/(dm_a dz_a), is given by the excursion set or the extended Press-Schechter formalism <cit.>.
Especially <cit.> obtained analytical formulation for the distribution that provides good fit to the numerical simulation data over a large range of m/M and z.
 their model III.


The subhalo mass function is obtained as

    dN_ sh(m|M,z)/dm = ∫ dm_a ∫ dz_a d^2N_ sh/dm_adz_a∫ dc_a P(c_a|m_a,z_a) δ(m-m(z|m_a,z_a,c_a)),

where P(c_a|m_a,z_a) is the probability distrbution for c_a given m_a and z_a, for which  the log-normal distribution with the mean c̅_a(m_a,z_a) <cit.> and the standard deviation σ_log c_a = 0.13 <cit.>.
We show the subhalo mass functions obtained with Equation (<ref>) for various values of M and z in Figure <ref>, where comparison is made with simulation results of similar host halos.


Halos and subhalos formed in these simulations were identified with ROCKSTAR phase space halo finder <cit.>. The bound mass is used as the subhalo mass, which nearly
corresponds to the tidal mass <cit.>.

For all these halos, one can see remarkable agreement between the analytic model and the corresponding simulation results in resolved regimes.
Successfully reproducing behaviors at resolved regimes, this analytic model is able to make reliable predictions of the subhalo mass functions below resolutions of the numerical simulations, without relying on extrapolating a single power-law functions, from which most of the previous studies in the literature had to suffer.
The subhalo mass fraction is then obtained as

    f_ sh(M,z)    =   1/M∫ dm mdN_ sh(m|M,z)/dm
       =   1/M∫ dm_a ∫ dz_a d^2N_ sh/dm_adz_a∫ dc_a P(c_a|m_a,z_a)m(z|m_a,z_a,c_a),

and is shown in Figure <ref> (bottom right) for various values of redshifts.
The subhalo mass fraction is found to increase as a function of M and z.
At higher redshifts, since there is shorter time for the subhalos to experience tidal mass loss, f_ sh is larger.
Again, a good agreement in f_ sh is found between the analytic model and the simulation results by <cit.>.




The annihilation boost factor is then

    B_ sh^(0)(M,z) = 1/L_ host,0(M,z)∫ dm_a ∫ dz_a d^2N_ sh/dm_adz_a∫ dc_a P(c_a|m_a,z_a) L_ sh^(0)(m_a,z_a,c_a|M,z),

which is to be compared with Equation (<ref>) that was derived with a simpler (and unrealistic) discussion.
The superscript (0) represents the quantity in the absense of sub-subhalos and beyond.
The subhalo luminosity, L_ sh^(0)(m_a,z_a,c_a|M,z), is proportional to the volume integral of density squared ρ_ sh^2(r) out to the truncation radius,

    L_ sh^(0)(m_a,z_a,c_a|M,z) ∝∫ d^3xρ_ sh^2( x) = 4π/3ρ_s, sh^2r_s, sh^3[1-1/(1+r_t, sh/r_s, sh)^3],

where ρ_s, sh, r_s, sh and r_t, sh are functions of (m_a,z_a,c_a) as well as M and z.





Then, the effect of sub^n-subhalos (for n≥ 1) can be estimated iteratively.
At nth iteraction, when a subhalo accreted onto its host at z_a with m_a, it is assigned a sub-subhalo boost factor B_ sh^(n-1)(m_a,z_a).

After the accretion, the outer region of the subhalo is stripped away by the tidal force and thus all the sub-subhalos within this stripped region will disappear, reducing the sub-subhalo boost accordingly.
<cit.> assumed that the sub-subhalos were distributed within the subhalo following n_ ssh(r) ∝ [1+(r/r_s)^2]^-3/2.
The luminosity due to sub-subhalos within a radius r is therefore proportional to their enclosed number

    N_ ssh(<r|r_s) = ∫_0^r dr' 4π r'^2 n_ ssh(r') ∝ r_s^3 [sinh^-1(r/r_s)-r/√(r^2+r_s^2)],

and it gets suppressed by a factor of N_ ssh(<r_t|r_s) / N_ ssh(<r_ vir|r_s,a) due to the tidal stripping.[We note that in estimating the effect of sub^n-subhalos in the boost factors, Reference <cit.> ignored the changes of ρ_s and r_s and hence did not include the factor of r_s^3 in Equation (<ref>) and ρ_s r_s^3 in Equation (<ref>). In addition, in Equation (<ref>), they multiplied L_ sh^(0) by a factor of 1+B_ ssh^(n) instead of 1-f_ ssh^2+B_ ssh^(n). We correct for all these effects in this review.]
The luminosity due to the smooth component also decreases as L_ sh^(0)(<r_t|ρ_s,r_s)/L_ sh^(0)(<r_ vir|ρ_s,a,r_s,a), where

    L_ sh^(0)(<r|ρ_s,r_s) ∝ρ_s^2r_s^3[1-1/(1+r/r_s)^3].

Thus the sub-subhalo boost after the nth iteration, B_ ssh^(n) is obtained by

    B_ ssh^(n)(m_a,z_a,c_a|M,z) = B_ sh^(n-1)(m_a,z_a)N_ ssh(<r_t|r_s)/N_ ssh(<r_ vir|r_s,a)/L_ sh^(0)(<r_t|ρ_s,r_s)/L_ sh^(0)(<r_ vir|ρ_s,a,r_s,a).

Similarly, the sub-subhalo mass fraction f_ ssh is obtained by

    f_ ssh(m_a,z_a,c_a|M,z) = f_ sh(m_a,z_a)N_ ssh(<r_t|r_s)/N_ ssh(<r_ vir|r_s,a)/m_ sm(<r_t|ρ_s,r_s)/m_ sm(<r_ vir|ρ_s,a,r_s,a),

where f_ sh(m_a,z_a) is obtained with Equations (<ref>) and m_ sm(<r|ρ_s,r_s) ∝ρ_sr_s^3f(r/r_s) is the enclosed mass within r of the smooth component of the subhalo.
The subhalo boost factor after nth iteration is obtained with Equation (<ref>) by replacing L_ sh^(0) with [1-f_ ssh^2 + B_ ssh^(n)]L_ sh^(0) [see discussions below Equation (<ref>)]:

    B_ sh^(n)(M,z)    =   1/L_ host,0(M,z)∫ dm_a ∫ dz_a d^2N_ sh/dm_adz_a∫ dc_a P(c_a|m_a,z_a)
        
          ×[1-f_ ssh^2(m_a,z_a,c_a|M,z) + B_ ssh^(n)(m_a,z_a,c_a|M,z)] L_ sh^(0)(m_a,z_a,c_a|M,z).

The host luminosity in the absence of the subhalos L_ host,0(M,z) is defined by marginalizing over the concentration parameter c_ vir:

    L_ host,0(M,z) ∝4π/3∫ dc_ virP(c_ vir|M,z)ρ_s^2(M,z,c_ vir)r_s^3(M,z,c_ vir)[1-1/(1+c_ vir)^3].

with the log-normal distribution P(c_ vir|M,z).


Figure <ref> shows the subhalo boost factors B_ sh as a function of the host mass M at various redshifts z (top left).
The boost factors are on the order of unity, while it can be as larger as ∼5 for cluster-size halos.
It is also noted that they are larger at higher redshifts, because the subhalos have less time to be disrupted.
The top right panel of Figure <ref> shows the effect of sub^n-subhalos, which is saturated after the second iteration. 
The contribution to the boost factors due to sub-subhalos and beyond is ≲10% for the hosts with M_ host≥10^13M_⊙.
The bottom left panel of Figure <ref> shows the luminosity ratio L_ total/L_ host,0 = 1-f_ sh^2 + B_ sh (Equation <ref>) as a function of the host masses for various values of the redshifts. 
The bottom right panel of Figure <ref> shows comparison with the results of the other work <cit.>. We note that the analytic models do not rely on the subhalo mass function prepared separately, as the models can provide them in a self-consistent manner. The resulting boost factors are, however, found to be more modest than the previous results. This is mainly because the subhalo mass function adopted in the literature is larger than the predictions of the analytic models. However, they might be larger because of halo-to-halo variance. See discrepancy between predictions of the subhalo mass function for the 1.8× 10^12M_⊙ halo by <cit.> and the result of <cit.> shown in the top left panel of Figure <ref>.




Finally, for convenience of the reader who might be interested in using the results without going into details of the formalism, we provide fitting functions for both the subhalo mass functions and the annihilation boost factors.
They are summarized in Appendix <ref>.















 §.§ Models for Self-Similar Subhalos



Assuming a self-similarity of the subhalos, <cit.> developed a fully analytic formulation for the probability distribution function of the dark matter density, P(ρ).
Then <cit.> applied the formulation to the result of cosmological N-body simulations, to obtain the fitting function for the Galactic local boost factor at Galactocentric radius r:

    B_ sh(r) = f_ sm(r) e^δ_f^2+[1-f_ sm(r)]1+α_K/1-α_K[(ρ_ max/ρ_χ(r))^1-α_K-1],

where f_ sm(r) is the volume fraction occupied by the smooth component and ρ_ max is the highest dark matter subhalo density.
Through the calibration with the numerical simulations,  found δ_f = 0.2, α_K = 0 and that the subhalo fraction was given by

    1-f_ sm(r) = κ[ρ_χ(r)/ρ_χ(100  kpc)]^-0.26,

where κ = 0.007.
<cit.> then suggested a larger value of κ = 0.15–0.2 to obtain a larger boost consistent with earlier work <cit.>.
The maximum subhalo density ρ_ max is estimated as ρ_ max = [c^3/f(c)/12]200 ρ_c(z_f), where c and z_f are the concentration parameter and collapse redshift of the smallest halos.
<cit.> adopted c = 3.5 and z_f = 40 and ρ_ max = 80 GeV cm^-3.
On the other hand, <cit.> obtained a smaller ρ_ max≈ 20 GeV cm^-3 even with a very small cutoff masses of M_ min = 10^-12 M_⊙
Within the virial radius of the Milky-Way halo, the subhalo boost factor for dark matter annihilation is found to be no greater than ∼10 <cit.>.



 §.§ Universal Clustering of Dark Matter in Phase Space



<cit.> investigated the behavior of dark matter particles that belong to the halo substructure in the phase space of distance and velocity.
Reconstructing the phase-space distribution using the Aquarius numerical simulations <cit.>, they found universality of the coarse-grained phase-space distribution ranging from dwarfs to clusters of galaxies.
They developed physically motivated models based on the stable clustering hypothesis, spherical collapses and the tidal stripping of the subhalos and applied to the obtained phase-space distribution data from the simulations to find a good agreement.
Then, they computed the nonlinear matter power spectrum based on the halo model <cit.> down to very small free-streaming cutoff scales.
Based on the power spectrum, they obtained the subhalo boost factor greater than ∼30–100 for the Milky-Way size halos, which is significantly larger than the values obtained with other analytic work <cit.>. This discrepancy might come from the treatment at very small scales, where it is very hard to calibrate the analytic models against the results of the numerical simulations.




§ CONCLUSIONS



It is established that dark matter halos are made up with lots of substructures.
Especially in the cold dark matter scenario, small structures form first and merge and accrete to create larger halos.
If the dark matter is made of weakly interacting massive particles such as the supersymmetric neutralino, the smallest halos can be as light as or even lighter than the Earth.
The rate of dark matter annihilation and hence its signatures such as gamma-ray fluxes are proportional to dark matter density squared and therefore, having small-scale “clumpy” subhalos will boost the signals.

It is, however, of an extraordinary challenge to estimate this subhalo boost factor, that is, the ratio of luminosity from dark matter annihilation in the subhalos to that in the smoothly distributed main component.
This is mainly because subhalos of all the mass scales ranging from Earth to galaxy masses can contribute to the boost factor nearly equally per decade in mass.
In this review, we cover recent progress to overcome this issue to obtain realistic and unbiased estimates on the subhalo boost factor that will impact on interpretation of the measurements on particle physics parameters such as the annihilation cross section.
While cosmological N-body simulations provide the most accurate avenue to study structures in highly nonlinear regime, it is inevitably limited by the numerical resolution.
Even the state-of-the-art N-body simulations <cit.> can resolve subhalos ranging for only several decades, which is still more than ten orders of magnitude in short to resolve all the subhalos.
Therefore, the boost estimates have to rely on extrapolation of the subhalo properties such as its mass function and concentration parameter, which are often well described with power-law functions. 
Danger of extrapolating trends found in resolved regime for other many orders of magnitude had been widely acknowledged but nevertheless, it was found that the estimates based on such extrapolations tended to give very large amount of boost factor of ∼100 (∼1000) for galaxy (cluster) size halos <cit.>.

As a complementary approach, analytic models have been investigated.
They are based on self similar propertiese of the subhalos <cit.>, universal phase-space distribution <cit.> and extended Press-Schechter formalism combined with tidal stripping modeling <cit.>.
(More recent numerical approach also adopts the concentration-mass relation calibrated for the subhalos in order to take the tidal effects into account <cit.>.)
Most importantly, these are all calibrated with the cosmological N-body simulations at resolved regimes and proven to reproduce the simulation results such as the subhalo mass functions.
For example, the most recent analytic models by Reference <cit.> predict the subhalo mass functions for various host masses and redshifts, which are found to be in good agreement with the simulation results <cit.>.
The annihilation boost factors based on these analytic models tend to be more modest, 𝒪(1) for galaxy-size halos and ≲𝒪(10) for cluster-size halos.
However, none of these models have been tested against simulations at very small host halos that are less massive than 10^6 M_⊙.
Simulations of microhalos with 10^-6M_⊙ suggest cuspier profiles towards halo centers such as r^-1.5 <cit.> and if this is the case for the subhalos too, it would boost the annihilation rate further.

It is known that including baryons in the simulations affects properties of subhalos such as spatial distribution and density profiles (e.g., References <cit.>) and hence there might be some effect on the annihilation boost factors.
This, however, remains largely unexplored and has to wait for future progress.
However, since the subhalos of all masses ranging down to about the Earth mass contribute to the boost factors and the baryons will likely affect only halos of dwarf galaxies or larger, we anticipate that it is not a very important effect for the annihilation boost factors.

The subhalo boosts directly impact the obtained upper limits on the dark matter annihilation cross section from the extragalactic halo observations.
Therefore, to obtain the most accurate estimates of the boost factor by reducing uncertainties on structure formation at small scales as well as the physics of tidal stripping is of extreme importance for the indirect searches for particle dark matter through self-annihilation with the current and near future observations of high-energy gamma-rays, neutrinos and charged cosmic rays.









SA wrote most of the review sections except for Section <ref>, which was contributed by TI including Figures <ref>–<ref>.  NH prepared Figures <ref> and <ref> and wrote Appendix <ref>. 





This work was supported by JSPS KAKENHI Grant Numbers JP17H04836 (SA), JP18H04340 (SA and NH), JP18H04578 (SA), JP15H01030 (TI), JP17H04828 (TI), JP17H01101 (TI) and JP18H04337 (TI). TI has been supported by MEXT as “Priority Issue on Post-K computer” (Elucidation of the Fundamental Laws and Evolution of the Universe) and JICFuS. Numerical computations were partially carried out on the K computer at the RIKEN Advanced Institute for Computational Science (Proposal numbers hp150226, hp160212, hp170231, hp180180), and Aterui supercomputer at Center for Computational Astrophysics, CfCA, of National Astronomical Observatory of Japan.





The authors declare no conflict of interest. 




yes 
multiple 


§ FITTING FORMULAE



In this appendix, we provide fitting functions obtained with the analytical calculation in <cit.> that covers more than twenty orders of magnitude in the mass range and the redshift up to ∼10. The mass function of the subhalo m^2dN_ sh/dm, the luminosity of the subhalo separated from the particle physics factors L_ sh, and the boost factor B_ sh=L_ sh/L_ host,0 [see Equation (<ref>)] as functions of the host mass and the redshift are provided here. We also summarize fitting functions for the boost factor B_ sh at z=0 in the literature. Note that the host mass is always measured in units of the solar mass (M_⊙) in this appendix.



 §.§ Subhalo Mass Function


The fitting formula is written in the follwing form:

    m^2dN/dm=(a+bm^α)exp[-(m/m_c)^β],

introducing a cutoff mass m_c; a, b, and m_c in Equation (<ref>) are functions of the host mass and the redshift:

    a(M_ host,z)   =   6.0×10^-4M_host(z) (2-0.08log_10[M_ host(z)])^2(log_10[M_ host(z)/10^-5])^-0.1
          ×(1+z), 
    
    b(M_ host,z)   =   8.0×10^-5M_ host(z)(2-0.08log_10[M_ host(z)])(log_10[M_ host(z)/10^-5])^-0.08z
          ×(log_10[M_ host(z)/10^-8])^-1(log_10[M_ host(z)/10^18])^2, 
    
    m_c   =   0.05(1+z)M_host(z),
    α   =   0.2+0.02z, 
    β   =   3.


In Figure <ref>, we show the comparison between the mass function obtained in analytical calculations <cit.> and Equation (<ref>). By integrating Equation (<ref>), we obtain the mass fraction shown in Figure <ref>.







 §.§ Subhalo Luminosity


The luminosity of a subhalo is written as

    L_ host(M)∝ρ_s^2r_s^3[1-1/(1+c_ vir^3)],

assuming a subhalo of NFW profile. In this section, the characteristic density ρ_s and the scale radius r_s are measured in units of g/cm^3 and cm, respectively. 
For simplicity, we do not show the integration over the distribution of the virial concentration parameter P(c_ vir|M,z) in the above expression. The constant of  proportionality is detemined by fixing the particle physics model. However, we do not include it in the following expression since it cancels in the calculation of the boost factor B_ sh by taking the ratio of the host and subhalo luminosity L_ sh/L_ host. 

The fitting formula for the luminosity takes the form of

    log_10 L_ sh=b + a log_10 m,

where

    a   =   (-0.025z+0.18)(log_10[M_host(z)/10^-5])^0.3+(0.06z+0.53),
    
    b   =   -0.95log_10[M_ host(z)]+(0.1-0.015z)log_10[M_host(z)/10^4]+0.07.




 §.§ Annihilation Boost Factor

The boost factor is sensitive to the models of the concentration-mass relation. We provide fitting functions for two different concentration-mass relation models here. One coresponds to the canonical model in Reference <cit.> which assumes the concentration-mass relation derived in Reference <cit.>. The other corresponds to the concentration-mass relation in Reference <cit.>. For both cases, the fitting function of the boost factor is written in a combination of two sigmoid functions, f(x)=(1+e^-ax)^-1,

    log_10B_ sh=X(z)/1+e^-a(z)(log_10[M_ host]-m_1(z))+c(z)(1+Y(z)/1+e^-b(z)(log_10[M_ host]-m_2(z))).

Funcitons X, Y, a, b, c, m_1, and m_2 depend on the redshift but they do not on the host mass.


  * For Correa's concentration <cit.>

    X(z)   =   2.7e^-0.2z+0.15,
    
    Y(z)   =   0.4+(-0.224z+0.56)e^-0.8z, 
    
    a(z)   =   0.10+0.095e^-0.5z, 
    
    b(z)   =   0.03z^2-0.08z-0.83, 
    
    c(z)   =   0.004z^2-0.04z-0.6, 
    
    m_1   =   -3.17z+17.4, 
    
    m_2   =   -(0.2z-1)^5-4.


  * For Okoli's concentration <cit.>

    X(z)   =   2.2e^-0.75z+0.67, 
    
    Y(z)   =   2.5e^-0.005z+0.8, 
    
    a(z)   =   0.1e^-0.5z+0.22, 
    
    b(z)   =   0.8e^-0.5(z-12)^4-0.24, 
    
    c(z)   =   -0.0005z^3-0.032z^2+0.28z-1.12, 
    
    m_1   =   -2.6z+8.2, 
    
    m_2   =   0.1e^-3z-12.


All of these formulae are applicable for hosts at arbitrary redshifts up to z∼7.



 §.§ Fitting Functions for the Boost Factor in the Literature


Several works provide the fitting function for the boost factor at z=0. We summarize functions provided in <cit.>, and <cit.> here.




    
  * <cit.> have analyzed cluster scale halo of M_ host=[5,20]×10^14h^-1M_⊙ in the Phoenix project <cit.>. Subhalos down to m∼10^6 can be resolved in their calculations.

    B_ sh=1.6×10^-3(M_ host/M_⊙)^0.39.

    
  * <cit.> derive the boost factor based on the concentration-mass relation in Reference <cit.>. The fitting function is provided for their fiducial model assuming the minimum halo mass to be M_ min=10^-6M_⊙ and the subhalo mass function dN/dm∝ m^-2. Each subhalo is assumed to be a field halo.

    log_10B_ sh(z=0)=∑^5_i=0b_i(lnM_ host/M_⊙)^i

    with

    b_0   =   -0.442
    
        b_1   =   0.0796
    
        b_2   =   -0.0025
    
        b_3   =   4.77×10^-6
    
        b_4   =   4.77×10^-6
    
        b_5   =   -9.69×10^-8

    
  * <cit.> derive the boost factor taking the dependence of the survived halo properties on the distance from the host, that is, the host potential. The form of the function is similar to that of <cit.> but adopitng different log bases:

    log_10B_ sh(z=0)=∑^5_i=0b_i(log_10M_ host/M_⊙)^i

    and parameters b_i are

    b_0   =   -0.186
    
        b_1   =   0.144
    
        b_2   =   -8.8×10^-3
    
        b_3   =   1.13×10^-3
    
        b_4   =   -3.7×10^-5
    
        b_5   =   -2×10^-7

    for α=2 and

    b_0   =   -6.8×10^-2
    
        b_1   =   9.4×10^-2
    
        b_2   =   -9.8×10^-3
    
        b_3   =   1.05×10^-3
    
        b_4   =   -3.4×10^-5
    
        b_5   =   -2×10^-7

    for α=1.9.










References
999


[Bertone et al.(2005)Bertone, Hooper, and Silk]Bertone:2004pz
Bertone, G.; Hooper, D.; Silk, J.
Particle dark matter: Evidence, candidates and constraints.
Phys. Rept. 2005, 405, 279–390,
 
 
   doi:10.1016/j.physrep.2004.08.031.

[Bringmann and Weniger(2012)]Bringmann:2012ez
Bringmann, T.; Weniger, C.
Gamma Ray Signals from Dark Matter: Concepts, Status and Prospects.
Phys. Dark Univ. 2012, 1, 194–217,


  doi:10.1016/j.dark.2012.10.005.

[Jungman et al.(1996)Jungman, Kamionkowski, and
  Griest]Jungman:1995df
Jungman, G.; Kamionkowski, M.; Griest, K.
Supersymmetric dark matter.
Phys. Rept. 1996, 267, 195–373,


  doi:10.1016/0370-1573(95)00058-5.

[Hooper and Profumo(2007)]Hooper:2007qk
Hooper, D.; Profumo, S.
Dark matter and collider phenomenology of universal extra
  dimensions.
Phys. Rept. 2007, 453, 29–115,


  doi:10.1016/j.physrep.2007.09.003.

[Boveia and Doglioni(2018)]Boveia:2018yeb
Boveia, A.; Doglioni, C.
Dark Matter Searches at Colliders.
Ann. Rev. Nucl. Part. Sci. 2018, 68, 429–459,


  doi:10.1146/annurev-nucl-101917-021008.

[Steigman et al.(2012)Steigman, Dasgupta, and
  Beacom]Steigman:2012nb
Steigman, G.; Dasgupta, B.; Beacom, J.F.
Precise Relic WIMP Abundance and its Impact on Searches for Dark
  Matter Annihilation.
Phys. Rev. 2012, D86, 023506,


  doi:10.1103/PhysRevD.86.023506.

[Hofmann et al.(2001)Hofmann, Schwarz, and
  Stoecker]Hofmann:2001bi
Hofmann, S.; Schwarz, D.J.; Stoecker, H.
Damping scales of neutralino cold dark matter.
Phys. Rev. 2001, D64, 083507,


  doi:10.1103/PhysRevD.64.083507.

[Green et al.(2004)Green, Hofmann, and Schwarz]Green:2003un
Green, A.M.; Hofmann, S.; Schwarz, D.J.
The power spectrum of SUSY—CDM on sub-galactic scales.
Mon. Not. R. Astron. Soc. 2004, 353, L23,


  doi:10.1111/j.1365-2966.2004.08232.x.

[Loeb and Zaldarriaga(2005)]Loeb:2005pm
Loeb, A.; Zaldarriaga, M.
The Small-scale power spectrum of cold dark matter.
Phys. Rev. 2005, D71, 103520,


  doi:10.1103/PhysRevD.71.103520.

[Bertschinger(2006)]Bertschinger:2006nq
Bertschinger, E.
The Effects of Cold Dark Matter Decoupling and Pair Annihilation on
  Cosmological Perturbations.
Phys. Rev. 2006, D74, 063509,


  doi:10.1103/PhysRevD.74.063509.

[Profumo et al.(2006)Profumo, Sigurdson, and
  Kamionkowski]Profumo:2006bv
Profumo, S.; Sigurdson, K.; Kamionkowski, M.
What mass are the smallest protohalos?
Phys. Rev. Lett. 2006, 97, 031301,


  doi:10.1103/PhysRevLett.97.031301.

[Bringmann(2009)]Bringmann:2009vf
Bringmann, T.
Particle Models and the Small-Scale Structure of Dark Matter.
New J. Phys. 2009, 11, 105027,


  doi:10.1088/1367-2630/11/10/105027.

[Cornell and Profumo(2012)]Cornell:2012tb
Cornell, J.M.; Profumo, S.
Earthly probes of the smallest dark matter halos.
J. Cosmol. Astropart. Phys. 2012, 1206, 011,


  doi:10.1088/1475-7516/2012/06/011.

[Diamanti et al.(2015)Diamanti, Catalan, and
  Ando]Diamanti:2015kma
Diamanti, R.; Catalan, M.E.C.; Ando, S.
Dark matter protohalos in a nine parameter MSSM and implications for
  direct and indirect detection.
Phys. Rev. 2015, D92, 065029,


  doi:10.1103/PhysRevD.92.065029.

[Silk and Stebbins(1993)]Silk:1992bh
Silk, J.; Stebbins, A.
Clumpy cold dark matter.
Astrophys. J. 1993, 411, 439–449, doi:10.1086/172846.

[Bergstrom et al.(1999)Bergstrom, Edsjo, Gondolo, and
  Ullio]Bergstrom:1998jj
Bergstrom, L.; Edsjo, J.; Gondolo, P.; Ullio, P.
Clumpy neutralino dark matter.
Phys. Rev. 1999, D59, 043506,


  doi:10.1103/PhysRevD.59.043506.

[Bergstrom et al.(1998)Bergstrom, Edsjo, and
  Ullio]Bergstrom:1998zs
Bergstrom, L.; Edsjo, J.; Ullio, P.
Possible indications of a clumpy dark matter halo.
Phys. Rev. 1998, D58, 083507,


  doi:10.1103/PhysRevD.58.083507.

[Calcaneo-Roldan and Moore(2000)]CalcaneoRoldan:2000yt
Calcaneo-Roldan, C.; Moore, B.
The Surface brightness of dark matter: Unique signatures of
  neutralino annihilation in the galactic halo.
Phys. Rev. 2000, D62, 123005,


   doi:10.1103/PhysRevD.62.123005.

[Tasitsiomi and Olinto(2002)]Tasitsiomi:2002vh
Tasitsiomi, A.; Olinto, A.V.
The Detectability of neutralino clumps via atmospheric Cherenkov
  telescopes.
Phys. Rev. 2002, D66, 083006,


  doi:10.1103/PhysRevD.66.083006.

[Stoehr et al.(2003)Stoehr, White, Springel, Tormen, and
  Yoshida]Stoehr:2003hf
Stoehr, F.; White, S.D.M.; Springel, V.; Tormen, G.; Yoshida, N.
Dark matter annihilation in the halo of the Milky Way.
Mon. Not. R. Astron. Soc. 2003, 345, 1313,


  doi:10.1046/j.1365-2966.2003.07052.x.

[Koushiappas et al.(2004)Koushiappas, Zentner, and
  Walker]Koushiappas:2003bn
Koushiappas, S.M.; Zentner, A.R.; Walker, T.P.
The observability of gamma-rays from neutralino annihilations in
  Milky Way substructure.
Phys. Rev. 2004, D69, 043501,


  doi:10.1103/PhysRevD.69.043501.

[Baltz et al.(2007)Baltz, Taylor, and Wai]Baltz:2006sv
Baltz, E.A.; Taylor, J.E.; Wai, L.L.
Can Astrophysical Gamma Ray Sources Mimic Dark Matter Annihilation
  in Galactic Satellites?
Astrophys. J. 2007, 659, L125–L128,


  doi:10.1086/517882.

[Ando(2005)]Ando:2005hr
Ando, S.
Can dark matter annihilation dominate the extragalactic gamma-ray
  background?
Phys. Rev. Lett. 2005, 94, 171303,


   doi:10.1103/PhysRevLett.94.171303.

[Oda et al.(2005)Oda, Totani, and Nagashima]Oda:2005nv
Oda, T.; Totani, T.; Nagashima, M.
Gamma-ray background from neutralino annihilation in the first
  cosmological objects.
Astrophys. J. 2005, 633, L65–L68,


  doi:10.1086/497691.

[Pieri et al.(2005)Pieri, Branchini, and Hofmann]Pieri:2005pg
Pieri, L.; Branchini, E.; Hofmann, S.
Difficulty of detecting minihalos via gamm rays from dark matter
  annihilation.
Phys. Rev. Lett. 2005, 95, 211301,


  doi:10.1103/PhysRevLett.95.211301.

[Koushiappas(2006)]Koushiappas:2006qq
Koushiappas, S.M.
Proper motion of gamma-rays from microhalo sources.
Phys. Rev. Lett. 2006, 97, 191301,


  doi:10.1103/PhysRevLett.97.191301.

[Ando et al.(2007)Ando, Komatsu, Narumoto, and
  Totani]Ando:2006cr
Ando, S.; Komatsu, E.; Narumoto, T.; Totani, T.
Dark matter annihilation or unresolved astrophysical sources?
  Anisotropy probe of the origin of cosmic gamma-ray background.
Phys. Rev. 2007, D75, 063519,


  doi:10.1103/PhysRevD.75.063519.

[Pieri et al.(2008)Pieri, Bertone, and Branchini]Pieri:2007ir
Pieri, L.; Bertone, G.; Branchini, E.
Dark Matter Annihilation in Substructures Revised.
Mon. Not. R. Astron. Soc. 2008, 384, 1627,


  doi:10.1111/j.1365-2966.2007.12828.x.

[Diemand et al.(2007)Diemand, Kuhlen, and Madau]Diemand:2006ik
Diemand, J.; Kuhlen, M.; Madau, P.
Dark matter substructure and gamma-ray annihilation in the Milky Way
  halo.
Astrophys. J. 2007, 657, 262–270,


  doi:10.1086/510736.

[Berezinsky et al.(2007)Berezinsky, Dokuchaev, and
  Eroshenko]Berezinsky:2006qm
Berezinsky, V.; Dokuchaev, V.; Eroshenko, Y.
Anisotropy of dark matter annihilation with respect to the Galactic
  plane.
J. Cosmol. Astropart. Phys. 2007, 0707, 011,


  doi:10.1088/1475-7516/2007/07/011.

[Lavalle et al.(2008)Lavalle, Yuan, Maurin, and
  Bi]Lavalle:1900wn
Lavalle, J.; Yuan, Q.; Maurin, D.; Bi, X.J.
Full Calculation of Clumpiness Boost factors for Antimatter Cosmic
  Rays in the light of Lambda-CDM N-body simulation results. Abandoning hope in
  clumpiness enhancement?
Astron. Astrophys. 2008, 479, 427–452,


  doi:10.1051/0004-6361:20078723.

[Siegal-Gaskins(2008)]SiegalGaskins:2008ge
Siegal-Gaskins, J.M.
Revealing dark matter substructure with anisotropies in the diffuse
  gamma-ray background.
J. Cosmol. Astropart. Phys. 2008, 0810, 040,


  doi:10.1088/1475-7516/2008/10/040.

[Ando et al.(2008)Ando, Kamionkowski, Lee, and
  Koushiappas]Ando:2008br
Ando, S.; Kamionkowski, M.; Lee, S.K.; Koushiappas, S.M.
Can proper motions of dark-matter subhalos be detected?
Phys. Rev. 2008, D78, 101301,


   doi:10.1103/PhysRevD.78.101301.

[Lee et al.(2009)Lee, Ando, and Kamionkowski]Lee:2008fm
Lee, S.K.; Ando, S.; Kamionkowski, M.
The Gamma-Ray-Flux Probability Distribution Function from Galactic
  Halo Substructure.
J. Cosmol. Astropart. Phys. 2009, 0907, 007,


  doi:10.1088/1475-7516/2009/07/007.

[Kamionkowski and Koushiappas(2008)]Kamionkowski:2008vw
Kamionkowski, M.; Koushiappas, S.M.
Galactic substructure and direct detection of dark matter.
Phys. Rev. 2008, D77, 103509,


  doi:10.1103/PhysRevD.77.103509.

[Springel et al.(2008)Springel, White, Frenk, Navarro,
  Jenkins, Vogelsberger, Wang, Ludlow, and Helmi]Springel:2008by
Springel, V.; White, S.D.M.; Frenk, C.S.; Navarro, J.F.; Jenkins, A.;
  Vogelsberger, M.; Wang, J.; Ludlow, A.; Helmi, A.
Prospects for detecting supersymmetric dark matter in the Galactic
  halo.
Nature 2008, 456, 73–76,


  doi:10.1038/nature07411.

[Springel et al.(2008)Springel, Wang, Vogelsberger, Ludlow,
  Jenkins, Helmi, Navarro, Frenk, and White]Springel:2008cc
Springel, V.; Wang, J.; Vogelsberger, M.; Ludlow, A.; Jenkins, A.; Helmi, A.;
  Navarro, J.F.; Frenk, C.S.; White, S.D.M.
The Aquarius Project: The subhalos of galactic halos.
Mon. Not. R. Astron. Soc. 2008, 391, 1685–1711,


   doi:10.1111/j.1365-2966.2008.14066.x.

[Ando(2009)]Ando:2009fp
Ando, S.
Gamma-ray background anisotropy from galactic dark matter
  substructure.
Phys. Rev. 2009, D80, 023520,


  doi:10.1103/PhysRevD.80.023520.

[Kamionkowski et al.(2010)Kamionkowski, Koushiappas, and
  Kuhlen]Kamionkowski:2010mi
Kamionkowski, M.; Koushiappas, S.M.; Kuhlen, M.
Galactic Substructure and Dark Matter Annihilation in the Milky Way
  Halo.
Phys. Rev. 2010, D81, 043532,


  doi:10.1103/PhysRevD.81.043532.

[Pinzke et al.(2011)Pinzke, Pfrommer, and
  Bergstrom]Pinzke:2011ek
Pinzke, A.; Pfrommer, C.; Bergstrom, L.
Prospects of detecting gamma-ray emission from galaxy clusters:
  Cosmic rays and dark matter annihilations.
Phys. Rev. 2011, D84, 123509,


  doi:10.1103/PhysRevD.84.123509.

[Gao et al.(2012)Gao, Frenk, Jenkins, Springel, and
  White]Gao:2011rf
Gao, L.; Frenk, C.S.; Jenkins, A.; Springel, V.; White, S.D.M.
Where will supersymmetric dark matter first be seen?
Mon. Not. R. Astron. Soc. 2012, 419, 1721,


  doi:10.1111/j.1365-2966.2011.19836.x.

[Fornasa et al.(2013)Fornasa, Zavala, Sanchez-Conde,
  Siegal-Gaskins, Delahaye, Prada, Vogelsberger, Zandanel, and
  Frenk]Fornasa:2012gu
Fornasa, M.; Zavala, J.; Sanchez-Conde, M.A.; Siegal-Gaskins, J.M.; Delahaye,
  T.; Prada, F.; Vogelsberger, M.; Zandanel, F.; Frenk, C.S.
Characterization of Dark-Matter-induced anisotropies in the diffuse
  gamma-ray background.
Mon. Not. R. Astron. Soc. 2013, 429, 1529–1553,


  doi:10.1093/mnras/sts444.

[Ng et al.(2014)Ng, Laha, Campbell, Horiuchi, Dasgupta, Murase,
  and Beacom]Ng:2013xha
Ng, K.C.Y.; Laha, R.; Campbell, S.; Horiuchi, S.; Dasgupta, B.; Murase, K.;
  Beacom, J.F.
Resolving small-scale dark matter structures using multisource
  indirect detection.
Phys. Rev. 2014, D89, 083001,


  doi:10.1103/PhysRevD.89.083001.

[Sánchez-Conde and Prada(2014)]Sanchez-Conde:2013yxa
Sánchez-Conde, M.A.; Prada, F.
The flattening of the concentration–mass relation towards low halo
  masses and its implications for the annihilation signal boost.
Mon. Not. R. Astron. Soc. 2014, 442, 2271–2277,


  doi:10.1093/mnras/stu1014.

[Ishiyama(2014)]Ishiyama:2014uoa
Ishiyama, T.
Hierarchical Formation of Dark Matter Halos and the Free Streaming
  Scale.
Astrophys. J. 2014, 788, 27,


  doi:10.1088/0004-637X/788/1/27.

[Bartels and Ando(2015)]Bartels:2015uba
Bartels, R.; Ando, S.
Boosting the annihilation boost: Tidal effects on dark matter
  subhalos and consistent luminosity modeling.
Phys. Rev. 2015, D92, 123508,


   doi:10.1103/PhysRevD.92.123508.

[Stref and Lavalle(2017)]Stref:2016uzb
Stref, M.; Lavalle, J.
Modeling dark matter subhalos in a constrained galaxy: Global mass
  and boosted annihilation profiles.
Phys. Rev. 2017, D95, 063003,


  doi:10.1103/PhysRevD.95.063003.

[Moliné et al.(2017)Moliné, Sánchez-Conde,
  Palomares-Ruiz, and Prada]Moline:2016pbm
Moliné, Á.; Sánchez-Conde, M.A.; Palomares-Ruiz, S.; Prada, F.
Characterization of subhalo structural properties and implications
  for dark matter annihilation signals.
Mon. Not. R. Astron. Soc. 2017, 466, 4974–4990,


  doi:10.1093/mnras/stx026.

[Zavala and Afshordi(2016)]Zavala:2015ura
Zavala, J.; Afshordi, N.
Universal clustering of dark matter in phase space.
Mon. Not. R. Astron. Soc. 2016, 457, 986–992,


  doi:10.1093/mnras/stw048.

[Hiroshima et al.(2018)Hiroshima, Ando, and
  Ishiyama]Hiroshima:2018kfv
Hiroshima, N.; Ando, S.; Ishiyama, T.
Modeling evolution of dark matter substructure and annihilation
  boost.
Phys. Rev. 2018, D97, 123002,


  doi:10.1103/PhysRevD.97.123002.

[Diemand et al.(2008)Diemand, Kuhlen, Madau, Zemp, Moore, Potter,
  and Stadel]Diemand:2008in
Diemand, J.; Kuhlen, M.; Madau, P.; Zemp, M.; Moore, B.; Potter, D.; Stadel, J.
Clumps and streams in the local dark matter distribution.
Nature 2008, 454, 735–738,


  doi:10.1038/nature07153.

[Stadel et al.(2009)Stadel, Potter, Moore, Diemand, Madau, Zemp,
  Kuhlen, and Quilis]Stadel:2008pn
Stadel, J.; Potter, D.; Moore, B.; Diemand, J.; Madau, P.; Zemp, M.; Kuhlen,
  M.; Quilis, V.
Quantifying the heart of darkness with GHALO—A multi-billion
  particle simulation of our galactic halo.
Mon. Not. R. Astron. Soc. 2009, 398, L21–L25,


  doi:10.1111/j.1745-3933.2009.00699.x.

[Berezinsky et al.(2014)Berezinsky, Dokuchaev, and
  Eroshenko]Berezinsky:2014wya
Berezinsky, V.S.; Dokuchaev, V.I.; Eroshenko, Y.N.
Small-scale clumps of dark matter.
Phys. Usp. 2014, 57, 1–36; Reprinted in  Usp. Fiz. Nauk 2014, 184, 3, doi:10.3367/UFNe.0184.201401a.0003. 





[Press and Schechter(1974)]Press:1973iz
Press, W.H.; Schechter, P.
Formation of galaxies and clusters of galaxies by selfsimilar
  gravitational condensation.
Astrophys. J. 1974, 187, 425–438, doi:10.1086/152650.

[Bond et al.(1991)Bond, Cole, Efstathiou, and Kaiser]Bond:1990iw
Bond, J.R.; Cole, S.; Efstathiou, G.; Kaiser, N.
Excursion set mass functions for hierarchical Gaussian
  fluctuations.
Astrophys. J. 1991, 379, 440, doi:10.1086/170520.

[Anderson et al.(2016)Anderson, Zimmer, Conrad, Gustafsson,
  Sánchez-Conde, and Caputo]Anderson:2015dpc
Anderson, B.; Zimmer, S.; Conrad, J.; Gustafsson, M.; Sánchez-Conde, M.;
  Caputo, R.
Search for Gamma-Ray Lines towards Galaxy Clusters with the
  Fermi-LAT.
J. Cosmol. Astropart. Phys. 2016, 1602, 026,


  doi:10.1088/1475-7516/2016/02/026.

[Navarro et al.(1996)Navarro, Frenk, and White]Navarro:1995iw
Navarro, J.F.; Frenk, C.S.; White, S.D.M.
The Structure of cold dark matter halos.
Astrophys. J. 1996, 462, 563–575,


  doi:10.1086/177173.

[Navarro et al.(1997)Navarro, Frenk, and White]Navarro:1996gj
Navarro, J.F.; Frenk, C.S.; White, S.D.M.
A Universal density profile from hierarchical clustering.
Astrophys. J. 1997, 490, 493–508,


  doi:10.1086/304888.

[Charbonnier et al.(2012)Charbonnier, Combet, and
  Maurin]2012CoPhC.183..656C
Charbonnier, A.; Combet, C.; Maurin, D.
CLUMPY: A code for -ray signals from dark
  matter structures.
Comput. Phys. Commun. 2012, 183, 656–668,  


   doi:10.1016/j.cpc.2011.10.017.

[Nezri et al.(2012)Nezri, White, Combet, Maurin, Pointecouteau,
  and Hinton]Nezri:2012tu
Nezri, E.; White, R.; Combet, C.; Maurin, D.; Pointecouteau, E.; Hinton, J.A.
gamma-rays from annihilating dark matter in galaxy clusters:
  Stacking vs single source analysis.
Mon. Not. R. Astron. Soc. 2012, 425, 477,


  doi:10.1111/j.1365-2966.2012.21484.x.

[Ando and Komatsu(2013)]Ando:2013ff
Ando, S.; Komatsu, E.
Constraints on the annihilation cross section of dark matter
  particles from anisotropies in the diffuse gamma-ray background measured with
  Fermi-LAT.
Phys. Rev. 2013, D87, 123539,


  doi:10.1103/PhysRevD.87.123539.

[Hütten et al.(2016)Hütten, Combet, Maier, and
  Maurin]Hutten:2016jko
Hütten, M.; Combet, C.; Maier, G.; Maurin, D.
Dark matter substructure modelling and sensitivity of the Cherenkov
  Telescope Array to Galactic dark halos.
J. Cosmol. Astropart. Phys. 2016, 1609, 047,


  doi:10.1088/1475-7516/2016/09/047.

[Bryan and Norman(1998)]Bryan:1997dn
Bryan, G.L.; Norman, M.L.
Statistical properties of x-ray clusters: Analytic and numerical
  comparisons.
Astrophys. J. 1998, 495, 80,


  doi:10.1086/305262.

[Correa et al.(2015)Correa, Wyithe, Schaye, and
  Duffy]Correa:2015dva
Correa, C.A.; Wyithe, J.S.B.; Schaye, J.; Duffy, A.R.
The accretion history of dark matter haloes—III. A physical
  model for the concentration–mass relation.
Mon. Not. R. Astron. Soc. 2015, 452, 1217–1232,


  doi:10.1093/mnras/stv1363.

[Ishiyama et al.(2013)Ishiyama, Makino, Portegies Zwart, Groen,
  Nitadori, Rieder, de Laat, McMillan, Hiraki, and Harfst]Ishiyama:2011af
Ishiyama, T.; Makino, J.; Portegies Zwart, S.; Groen, D.; Nitadori, K.; Rieder,
  S.; de Laat, C.; McMillan, S.; Hiraki, K.; Harfst, S.
The Cosmogrid Simulation: Statistical Properties of small Dark
  Matter Halos.
Astrophys. J. 2013, 767, 146,


  doi:10.1088/0004-637X/767/2/146.

[Anderhalden and Diemand(2013)]Anderhalden:2013wd
Anderhalden, D.; Diemand, J.
Density Profiles of CDM Microhalos and their Implications for
  Annihilation Boost Factors.
J. Cosmol. Astropart. Phys. 2013, 1304, 009, doi:10.1088/1475-7516/2013/04/009.


Erratum in  J. Cosmol. Astropart. Phys. 2013, 1308, E02,   doi:10.1088/1475-7516/2013/08/E02.  






[Ishiyama et al.(2016)Ishiyama, Sudo, Yokoi, Hasegawa,
  Tominaga, and Susa]2016ApJ...826....9I
Ishiyama, T.; Sudo, K.; Yokoi, S.; Hasegawa, K.; Tominaga, N.;
  Susa, H.
Where are the Low-mass Population III Stars?
Astrophys. J. 2016, 826, 9,
 
  doi:10.3847/0004-637X/826/1/9.

[Klypin et al.(1999)Klypin, Kravtsov, Valenzuela, and
  Prada]Klypin:1999uc
Klypin, A.A.; Kravtsov, A.V.; Valenzuela, O.; Prada, F.
Where are the missing Galactic satellites?
Astrophys. J. 1999, 522, 82–92,


  doi:10.1086/307643.

[Moore et al.(1999)Moore, Ghigna, Governato, Lake, Quinn, Stadel,
  and Tozzi]Moore:1999nt
Moore, B.; Ghigna, S.; Governato, F.; Lake, G.; Quinn, T.R.; Stadel, J.; Tozzi,
  P.
Dark matter substructure within galactic halos.
Astrophys. J. 1999, 524, L19–L22,


  doi:10.1086/312287.

[Diemand et al.(2004)Diemand, Moore, and Stadel]Diemand:2004kx
Diemand, J.; Moore, B.; Stadel, J.
Velocity and spatial biases in CDM subhalo distributions.
Mon. Not. R. Astron. Soc. 2004, 352, 535,


  doi:10.1111/j.1365-2966.2004.07940.x.

[Reed et al.(2005)Reed, Governato, Quinn, Gardner, Stadel, and
  Lake]Reed:2004vg
Reed, D.; Governato, F.; Quinn, T.R.; Gardner, J.; Stadel, J.; Lake, G.
Dark matter subhaloes in numerical simulations.
Mon. Not. R. Astron. Soc. 2005, 359, 1537–1548,


  doi:10.1111/j.1365-2966.2005.09020.x.

[Kase et al.(2007)Kase, Makino, and Funato]Kase:2006wd
Kase, H.; Makino, J.; Funato, Y.
Missing dwarf problem in galaxy clusters.
Publ. Astron. Soc. Jpn. 2007, 59, 1071,


  doi:10.1093/pasj/59.6.1071.

[Kamionkowski and Liddle(2000)]Kamionkowski:1999vp
Kamionkowski, M.; Liddle, A.R.
The Dearth of halo dwarf galaxies: Is there power on short scales?
Phys. Rev. Lett. 2000, 84, 4525–4528,


  doi:10.1103/PhysRevLett.84.4525.

[Spergel and Steinhardt(2000)]Spergel:1999mh
Spergel, D.N.; Steinhardt, P.J.
Observational evidence for selfinteracting cold dark matter.
Phys. Rev. Lett. 2000, 84, 3760–3763,


  doi:10.1103/PhysRevLett.84.3760.

[Susa and Umemura(2004)]Susa:2004ji
Susa, H.; Umemura, M.
Effects of early cosmic reionization on the substructure problem in
  galactic halo.
Astrophys. J. 2004, 610, L5–L8,


  doi:10.1086/423199.

[Benson et al.(2002)Benson, Frenk, Lacey, Baugh, and
  Cole]Benson:2001at
Benson, A.J.; Frenk, C.S.; Lacey, C.G.; Baugh, C.M.; Cole, S.
The effects of photoionization on galaxy formation. 2. Satellites in
  the local group.
Mon. Not. R. Astron. Soc. 2002, 333, 177,


  doi:10.1046/j.1365-8711.2002.05388.x.

[Stoehr et al.(2002)Stoehr, White, Tormen, and
  Springel]Stoehr:2002ht
Stoehr, F.; White, S.D.M.; Tormen, G.; Springel, V.
The Milky Way's satellite population in a lambdaCDM universe.
Mon. Not. R. Astron. Soc. 2002, 335, L84–L88,


  doi:10.1046/j.1365-8711.2002.05891.x.

[Kravtsov et al.(2004)Kravtsov, Gnedin, and
  Klypin]Kravtsov:2004cm
Kravtsov, A.V.; Gnedin, O.Y.; Klypin, A.A.
The Tumultuous lives of Galactic dwarfs and the missing satellites
  problem.
Astrophys. J. 2004, 609, 482–497,


  doi:10.1086/421322.

[Okamoto et al.(2010)Okamoto, Frenk, Jenkins, and
  Theuns]Okamoto:2009rz
Okamoto, T.; Frenk, C.S.; Jenkins, A.; Theuns, T.
The properties of satellite galaxies in simulations of galaxy
  formation.
Mon. Not. R. Astron. Soc. 2010, 406, 208–222,


  doi:10.1111/j.1365-2966.2010.16690.x.

[Ishiyama et al.(2009)Ishiyama, Fukushige, and
  Makino]Ishiyama:2008xe
Ishiyama, T.; Fukushige, T.; Makino, J.
Variation of the subhalo abundance in dark matter halos.
Astrophys. J. 2009, 696, 2115–2125,


  doi:10.1088/0004-637X/696/2/2115.

[Mao et al.(2015)Mao, Williamson, and Wechsler]Mao:2015yua
Mao, Y.Y.; Williamson, M.; Wechsler, R.H.
The Dependence of Subhalo Abundance on Halo Concentration.
Astrophys. J. 2015, 810, 21,


  doi:10.1088/0004-637X/810/1/21.

[Kravtsov et al.(2004)Kravtsov, Berlind, Wechsler, Klypin,
  Gottloeber, Allgood, and Primack]Kravtsov:2003sg
Kravtsov, A.V.; Berlind, A.A.; Wechsler, R.H.; Klypin, A.A.; Gottloeber, S.;
  Allgood, B.; Primack, J.R.
The Dark side of the halo occupation distribution.
Astrophys. J. 2004, 609, 35–49,


   doi:10.1086/420959.

[Gao et al.(2004)Gao, White, Jenkins, Stoehr, and
  Springel]Gao:2004au
Gao, L.; White, S.D.M.; Jenkins, A.; Stoehr, F.; Springel, V.
The Subhalo populations of lambda-CDM dark halos.
Mon. Not. R. Astron. Soc. 2004, 355, 819,


  doi:10.1111/j.1365-2966.2004.08360.x.

[Zentner et al.(2005)Zentner, Berlind, Bullock, Kravtsov, and
  Wechsler]Zentner:2004dq
Zentner, A.R.; Berlind, A.A.; Bullock, J.S.; Kravtsov, A.V.; Wechsler, R.H.
The Physics of galaxy clustering. 1. A Model for subhalo
  populations.
Astrophys. J. 2005, 624, 505–525,


  doi:10.1086/428898.

[van den Bosch et al.(2005)van den Bosch, Tormen, and
  Giocoli]vandenBosch:2004zs
van den Bosch, F.C.; Tormen, G.; Giocoli, C.
The Mass function and average mass loss rate of dark matter
  subhaloes.
Mon. Not. R. Astron. Soc. 2005, 359, 1029–1040,


  doi:10.1111/j.1365-2966.2005.08964.x.

[Giocoli et al.(2010)Giocoli, Tormen, Sheth, and van den
  Bosch]Giocoli:2009ie
Giocoli, C.; Tormen, G.; Sheth, R.K.; van den Bosch, F.C.
The Substructure Hierarchy in Dark Matter Haloes.
Mon. Not. R. Astron. Soc. 2010, 404, 502–517,


  doi:10.1111/j.1365-2966.2010.16311.x.

[Gao et al.(2011)Gao, Frenk, Boylan-Kolchin, Jenkins, Springel,
  and White]Gao:2010tn
Gao, L.; Frenk, C.S.; Boylan-Kolchin, M.; Jenkins, A.; Springel, V.; White,
  S.D.M.
The statistics of the subhalo abundance of dark matter haloes.
Mon. Not. R. Astron. Soc. 2011, 410, 2309,


  doi:10.1111/j.1365-2966.2010.17601.x.

[Giocoli et al.(2008)Giocoli, Pieri, and Tormen]Giocoli:2007gf
Giocoli, C.; Pieri, L.; Tormen, G.
Analytical Approach to Subhaloes Population in Dark Matter Haloes.
Mon. Not. R. Astron. Soc. 2008, 387, 689–697,


  doi:10.1111/j.1365-2966.2008.13283.x.

[Yang et al.(2011)Yang, Mo, Zhang, and Bosch]Yang:2011rf
Yang, X.; Mo, H.J.; Zhang, Y.; Bosch, F.C.v.d.
An analytical model for the accretion of dark matter subhalos.
Astrophys. J. 2011, 741, 13,


  doi:10.1088/0004-637X/741/1/13.

[Jiang and van den Bosch(2016)]Jiang:2014nsa
Jiang, F.; van den Bosch, F.C.
Statistics of dark matter substructure—I. Model and universal
  fitting functions.
Mon. Not. R. Astron. Soc. 2016, 458, 2848–2869,


  doi:10.1093/mnras/stw439.

[Gunn and Gott(1972)]Gunn:1972sv
Gunn, J.E.; Gott, III, J.R.
On the Infall of Matter into Clusters of Galaxies and Some Effects
  on Their Evolution.
Astrophys. J. 1972, 176, 1–19, doi:10.1086/151605.

[Quinn et al.(1986)Quinn, Salmon, and
  Zurek]1986Natur.322..329Q
Quinn, P.J.; Salmon, J.K.; Zurek, W.H.
Primordial density fluctuations and the structure of galactic
  haloes.
Nature 1986, 322, 329–335, doi:10.1038/322329a0.

[Dubinski and Carlberg(1991)]Dubinski:1991bm
Dubinski, J.; Carlberg, R.G.
The Structure of cold dark matter halos.
Astrophys. J. 1991, 378, 496, doi:10.1086/170451.

[Hernquist(1990)]Hernquist:1990be
Hernquist, L.
An Analytical Model for Spherical Galaxies and Bulges.
Astrophys. J. 1990, 356, 359, doi:10.1086/168845.

[Fukushige and Makino(1997)]Fukushige:1996nr
Fukushige, T.; Makino, J.
On the Origin of Cusps in Dark Matter Halos.
Astrophys. J. 1997, 477, L9,


  doi:10.1086/310516.

[Moore et al.(1999)Moore, Quinn, Governato, Stadel, and
  Lake]Moore:1999gc
Moore, B.; Quinn, T.R.; Governato, F.; Stadel, J.; Lake, G.
Cold collapse and the core catastrophe.
Mon. Not. R. Astron. Soc. 1999, 310, 1147–1152,


  doi:10.1046/j.1365-8711.1999.03039.x.

[Ghigna et al.(2000)Ghigna, Moore, Governato, Lake, Quinn, and
  Stadel]Ghigna:1999sn
Ghigna, S.; Moore, B.; Governato, F.; Lake, G.; Quinn, T.R.; Stadel, J.
Density profiles and substructure of dark matter halos. Converging
  results at ultra-high numerical resolution.
Astrophys. J. 2000, 544, 616,


  doi:10.1086/317221.

[Jing(2000)]Jing:1999km
Jing, Y.P.
The Density Profile of Equilibrium and Nonequilibrium Dark Matter
  Halos.
Astrophys. J. 2000, 535, 30–36,


  doi:10.1086/308809.

[Fukushige and Makino(2001)]Fukushige:2000ar
Fukushige, T.; Makino, J.
Structure of dark matter halos from hierarchical clustering.
Astrophys. J. 2001, 557, 533,


  doi:10.1086/321666.

[Klypin et al.(2001)Klypin, Kravtsov, Bullock, and
  Primack]Klypin:2000hk
Klypin, A.; Kravtsov, A.V.; Bullock, J.; Primack, J.
Resolving the structure of cold dark matter halos.
Astrophys. J. 2001, 554, 903–915,


  doi:10.1086/321400.

[Power et al.(2003)Power, Navarro, Jenkins, Frenk, White,
  Springel, Stadel, and Quinn]Power:2002sw
Power, C.; Navarro, J.F.; Jenkins, A.; Frenk, C.S.; White, S.D.M.; Springel,
  V.; Stadel, J.; Quinn, T.R.
The Inner structure of Lambda CDM halos. 1. A Numerical convergence
  study.
Mon. Not. R. Astron. Soc. 2003, 338, 14–34,


  doi:10.1046/j.1365-8711.2003.05925.x.

[Fukushige et al.(2004)Fukushige, Kawai, and
  Makino]Fukushige:2003xc
Fukushige, T.; Kawai, A.; Makino, J.
Structure of dark matter halos from hierarchical clustering. 3.
  Shallowing of the Inner cusp.
Astrophys. J. 2004, 606, 625–634,


  doi:10.1086/383192.

[Hayashi et al.(2004)Hayashi, Navarro, Power, Jenkins, Frenk,
  White, Springel, Stadel, and Quinn]Hayashi:2003sj
Hayashi, E.; Navarro, J.F.; Power, C.; Jenkins, A.R.; Frenk, C.S.; White,
  S.D.M.; Springel, V.; Stadel, J.; Quinn, T.R.
The Inner structure of lambda-CDM halos. 2. Halo mass profiles and
  LSB rotation curves.
Mon. Not. R. Astron. Soc. 2004, 355, 794–812,


  doi:10.1111/j.1365-2966.2004.08359.x.

[Kazantzidis et al.(2006)Kazantzidis, Zentner, and
  Kravtsov]Kazantzidis:2005su
Kazantzidis, S.; Zentner, A.R.; Kravtsov, A.V.
The robustness of dark matter density profiles in dissipationless
  mergers.
Astrophys. J. 2006, 641, 647–664,


  doi:10.1086/500579.

[Jing and Suto(2000)]Jing:1999ir
Jing, Y.P.; Suto, Y.
Density profiles of dark matter halo are not universal.
Astrophys. J. 2000, 529, L69–72,


  doi:10.1086/312463.

[Diemand et al.(2004)Diemand, Moore, and Stadel]Diemand:2004wh
Diemand, J.; Moore, B.; Stadel, J.
Convergence and scatter of cluster density profiles.
Mon. Not. R. Astron. Soc. 2004, 353, 624,


  doi:10.1111/j.1365-2966.2004.08094.x.

[Einasto(1965)]Einasto1965
Einasto, J.
On the Construction of a Composite Model for the Galaxy and on the
  Determination of the System of Galactic Parameters.
Tr. Astrofiz. Instituta Alma-Ata 1965, 5, 87–100.

[Bullock et al.(2001)Bullock, Kolatt, Sigad, Somerville, Kravtsov,
  Klypin, Primack, and Dekel]Bullock:1999he
Bullock, J.S.; Kolatt, T.S.; Sigad, Y.; Somerville, R.S.; Kravtsov, A.V.;
  Klypin, A.A.; Primack, J.R.; Dekel, A.
Profiles of dark haloes. Evolution, scatter, and environment.
Mon. Not. R. Astron. Soc. 2001, 321, 559–575,


  doi:10.1046/j.1365-8711.2001.04068.x.

[Zhao et al.(2003)Zhao, Jing, Mo, and Borner]Zhao:2003jf
Zhao, D.H.; Jing, Y.P.; Mo, H.J.; Borner, G.
Mass and redshift dependence of dark halo structure.
Astrophys. J. 2003, 597, L9–L12,


  doi:10.1086/379734.

[Maccio' et al.(2007)Maccio', Dutton, van den Bosch, Moore,
  Potter, and Stadel]Maccio:2006wpz
Maccio', A.V.; Dutton, A.A.; van den Bosch, F.C.; Moore, B.; Potter, D.;
  Stadel, J.
Concentration, Spin and Shape of Dark Matter Haloes: Scatter and the
  Dependence on Mass and Environment.
Mon. Not. R. Astron. Soc. 2007, 378, 55–71,


  doi:10.1111/j.1365-2966.2007.11720.x.

[Neto et al.(2007)Neto, Gao, Bett, Cole, Navarro, Frenk, White,
  Springel, and Jenkins]Neto:2007vq
Neto, A.F.; Gao, L.; Bett, P.; Cole, S.; Navarro, J.F.; Frenk, C.S.; White,
  S.D.M.; Springel, V.; Jenkins, A.
The statistics of lambda CDM Halo Concentrations.
Mon. Not. R. Astron. Soc. 2007, 381, 1450–1462,


  doi:10.1111/j.1365-2966.2007.12381.x.

[Maccio' et al.(2008)Maccio', Dutton, and Bosch]Maccio:2008pcd
Maccio', A.V.; Dutton, A.A.; Bosch, F.C.v.d.
Concentration, Spin and Shape of Dark Matter Haloes as a Function of
  the Cosmological Model: WMAP1, WMAP3 and WMAP5 results.
Mon. Not. R. Astron. Soc. 2008, 391, 1940–1954,


  doi:10.1111/j.1365-2966.2008.14029.x.

[Zhao et al.(2009)Zhao, Jing, Mo, and Boerner]Zhao:2008wd
Zhao, D.H.; Jing, Y.P.; Mo, H.J.; Boerner, G.
Accurate universal models for the mass accretion histories and
  concentrations of dark matter halos.
Astrophys. J. 2009, 707, 354–369,


  doi:10.1088/0004-637X/707/1/354.

[Klypin et al.(2011)Klypin, Trujillo-Gomez, and
  Primack]Klypin:2010qw
Klypin, A.; Trujillo-Gomez, S.; Primack, J.
Halos and galaxies in the standard cosmological model: Results from
  the Bolshoi simulation.
Astrophys. J. 2011, 740, 102,


  doi:10.1088/0004-637X/740/2/102.

[Prada et al.(2012)Prada, Klypin, Cuesta, Betancort-Rijo, and
  Primack]Prada:2011jf
Prada, F.; Klypin, A.A.; Cuesta, A.J.; Betancort-Rijo, J.E.; Primack, J.
Halo concentrations in the standard LCDM cosmology.
Mon. Not. R. Astron. Soc. 2012, 423, 3018–3030,


  doi:10.1111/j.1365-2966.2012.21007.x.

[Klypin et al.(2016)Klypin, Yepes, Gottlober, Prada, and
  Hess]Klypin:2014kpa
Klypin, A.; Yepes, G.; Gottlober, S.; Prada, F.; Hess, S.
MultiDark simulations: The story of dark matter halo concentrations
  and density profiles.
Mon. Not. R. Astron. Soc. 2016, 457, 4340–4359,


  doi:10.1093/mnras/stw248.

[Ludlow et al.(2016)Ludlow, Bose, Angulo, Wang, Hellwing, Navarro,
  Cole, and Frenk]Ludlow:2016ifl
Ludlow, A.D.; Bose, S.; Angulo, R.E.; Wang, L.; Hellwing, W.A.; Navarro, J.F.;
  Cole, S.; Frenk, C.S.
The mass–concentration–redshift relation of cold and warm dark
  matter haloes.
Mon. Not. R. Astron. Soc. 2016, 460, 1214–1232,


  doi:10.1093/mnras/stw1046.

[Wechsler et al.(2002)Wechsler, Bullock, Primack, Kravtsov, and
  Dekel]Wechsler:2001cs
Wechsler, R.H.; Bullock, J.S.; Primack, J.R.; Kravtsov, A.V.; Dekel, A.
Concentrations of dark halos from their assembly histories.
Astrophys. J. 2002, 568, 52–70,


  doi:10.1086/338765.

[Hellwing et al.(2016)Hellwing, Frenk, Cautun, Bose, Helly,
  Jenkins, Sawala, and Cytowski]Hellwing:2015upa
Hellwing, W.A.; Frenk, C.S.; Cautun, M.; Bose, S.; Helly, J.; Jenkins, A.;
  Sawala, T.; Cytowski, M.
The Copernicus Complexio: A high-resolution view of the small-scale
  Universe.
Mon. Not. R. Astron. Soc. 2016, 457, 3492–3509,


  doi:10.1093/mnras/stw214.

[Pilipenko et al.(2017)Pilipenko, Sánchez-Conde, Prada, and
  Yepes]Pilipenko:2017iae
Pilipenko, S.V.; Sánchez-Conde, M.A.; Prada, F.; Yepes, G.
Pushing down the low-mass halo concentration frontier with the
  Lomonosov cosmological simulations.
Mon. Not. R. Astron. Soc. 2017, 472, 4918–4927,


  doi:10.1093/mnras/stx2319.

[Diemand et al.(2007)Diemand, Kuhlen, and Madau]Diemand:2007qr
Diemand, J.; Kuhlen, M.; Madau, P.
Formation and evolution of galaxy dark matter halos and their
  substructure.
Astrophys. J. 2007, 667, 859–877,


  doi:10.1086/520573.

[Reed et al.(2011)Reed, Koushiappas, and Gao]Reed:2010gh
Reed, D.S.; Koushiappas, S.M.; Gao, L.
Non-universality of halo profiles and implications for dark matter
  experiments.
Mon. Not. R. Astron. Soc. 2011, 415, 3177–3188,


  doi:10.1111/j.1365-2966.2011.18930.x.

[Springel et al.(2005)Springel et al.]Springel:2005nw
Springel, V.; White, S.D.M.; Jenkin, A.; Frenk, C.S.; Yoshida, N.; Navarro, J.; Thacker, R.; Croton, D.; Helly, J.; Peacock, J.A; et al.
Simulating the joint evolution of quasars, galaxies and their
  large-scale distribution.
Nature 2005, 435, 629–636,


  doi:10.1038/nature03597.

[Ullio et al.(2002)Ullio, Bergstrom, Edsjo, and
  Lacey]Ullio:2002pj
Ullio, P.; Bergstrom, L.; Edsjo, J.; Lacey, C.G.
Cosmological dark matter annihilations into gamma-rays—A closer
  look.
Phys. Rev. 2002, D66, 123502,


  doi:10.1103/PhysRevD.66.123502.

[Pieri et al.(2011)Pieri, Lavalle, Bertone, and
  Branchini]Pieri:2009je
Pieri, L.; Lavalle, J.; Bertone, G.; Branchini, E.
Implications of High-Resolution Simulations on Indirect Dark Matter
  Searches.
Phys. Rev. 2011, D83, 023518,


  doi:10.1103/PhysRevD.83.023518.

[Ishiyama et al.(2015)Ishiyama, Enoki, Kobayashi, Makiya,
  Nagashima, and Oogi]Ishiyama:2014gla
Ishiyama, T.; Enoki, M.; Kobayashi, M.A.R.; Makiya, R.; Nagashima, M.; Oogi, T.
The ν^2GC simulations: Quantifying the dark side of the universe
  in the Planck cosmology.
Publ. Astron. Soc. Jpn. 2015, 67, 61,


  doi:10.1093/pasj/psv021.

[Makiya et al.(2016)Makiya, Enoki, Ishiyama, Kobayashi, Nagashima,
  Okamoto, Okoshi, Oogi, and Shirakata]Makiya:2015spa
Makiya, R.; Enoki, M.; Ishiyama, T.; Kobayashi, M.A.R.; Nagashima, M.; Okamoto,
  T.; Okoshi, K.; Oogi, T.; Shirakata, H.
The New Numerical Galaxy Catalog (ν^2GC): An updated
  semi-analytic model of galaxy and active galactic nucleus formation with
  large cosmological N-body simulations.
Publ. Astron. Soc. Jpn. 2016, 68, 25,


   doi:10.1093/pasj/psw005.

[van den Bosch and Ogiya(2018)]vandenBosch:2018tyt
van den Bosch, F.C.; Ogiya, G.
Dark Matter Substructure in Numerical Simulations: A Tale of
  Discreteness Noise, Runaway Instabilities, and Artificial Disruption.
Mon. Not. R. Astron. Soc. 2018, 475, 4066–4087,


  doi:10.1093/mnras/sty084.

[Peñarrubia et al.(2010)Penarrubia, Benson, Walker, Gilmore,
  McConnachie, and Mayer]Penarrubia:2010jk
Peñarrubia, J.; Benson, A.J.; Walker, M.G.; Gilmore, G.; McConnachie, A.;
  Mayer, L.
The impact of dark matter cusps and cores on the satellite galaxy
  population around spiral galaxies.
Mon. Not. R. Astron. Soc. 2010, 406, 1290,


  doi:10.1111/j.1365-2966.2010.16762.x.

[van den Bosch et al.(2018)van den Bosch, Ogiya, Hahn, and
  Burkert]vandenBosch:2017ynq
van den Bosch, F.C.; Ogiya, G.; Hahn, O.; Burkert, A.
Disruption of Dark Matter Substructure: Fact or Fiction?
Mon. Not. R. Astron. Soc. 2018, 474, 3043–3066,


  doi:10.1093/mnras/stx2956.

[Diemand et al.(2005)Diemand, Moore, and Stadel]Diemand:2005vz
Diemand, J.; Moore, B.; Stadel, J.
Earth-mass dark-matter haloes as the first structures in the early
  Universe.
Nature 2005, 433, 389–391,


  doi:10.1038/nature03270.

[Ishiyama et al.(2010)Ishiyama, Makino, and
  Ebisuzaki]Ishiyama:2010es
Ishiyama, T.; Makino, J.; Ebisuzaki, T.
Gamma-ray Signal from Earth-mass Dark Matter Microhalos.
Astrophys. J. 2010, 723, L195,


  doi:10.1088/2041-8205/723/2/L195.

[Angulo et al.(2017)Angulo, Hahn, Ludlow, and
  Bonoli]Angulo:2016qof
Angulo, R.E.; Hahn, O.; Ludlow, A.; Bonoli, S.
Earth-mass haloes and the emergence of NFW density profiles.
Mon. Not. R. Astron. Soc. 2017, 471, 4687–4701,


  doi:10.1093/mnras/stx1658.

[Ogiya and Hahn(2018)]Ogiya:2017hbr
Ogiya, G.; Hahn, O.
What sets the central structure of dark matter haloes?
Mon. Not. R. Astron. Soc. 2018, 473, 4339–4359,


  doi:10.1093/mnras/stx2639.

[Ogiya et al.(2016)Ogiya, Nagai, and Ishiyama]Ogiya:2016hyo
Ogiya, G.; Nagai, D.; Ishiyama, T.
Dynamical evolution of primordial dark matter haloes through
  mergers.
Mon. Not. R. Astron. Soc. 2016, 461, 3385–3396,


  doi:10.1093/mnras/stw1551.

[Gosenca et al.(2017)Gosenca, Adamek, Byrnes, and
  Hotchkiss]Gosenca:2017ybi
Gosenca, M.; Adamek, J.; Byrnes, C.T.; Hotchkiss, S.
3D simulations with boosted primordial power spectra and
  ultracompact minihalos.
Phys. Rev. 2017, D96, 123519,


  doi:10.1103/PhysRevD.96.123519.

[Delos et al.(2018a)Delos, Erickcek, Bailey, and
  Alvarez]Delos:2017thv
Delos, M.S.; Erickcek, A.L.; Bailey, A.P.; Alvarez, M.A.
Are ultracompact minihalos really ultracompact?
Phys. Rev. 2018, D97, 041303,


  doi:10.1103/PhysRevD.97.041303.

[Delos et al.(2018b)Delos, Erickcek, Bailey, and
  Alvarez]Delos:2018ueo
Delos, M.S.; Erickcek, A.L.; Bailey, A.P.; Alvarez, M.A.
Density profiles of ultracompact minihalos: Implications for
  constraining the primordial power spectrum.
Phys. Rev. 2018, D98, 063527,


  doi:10.1103/PhysRevD.98.063527.

[Polisensky and Ricotti(2015)]Polisensky:2015eya
Polisensky, E.; Ricotti, M.
Fingerprints of the initial conditions on the density profiles of
  cold and warm dark matter haloes.
Mon. Not. R. Astron. Soc. 2015, 450, 2172–2184,


  doi:10.1093/mnras/stv736.

[Gott(1975)]1975ApJ...201..296G
Gott, III, J.R.
On the Formation of Elliptical Galaxies.
Astrophys. J. 1975, 201, 296–310, doi:10.1086/153887.

[Fillmore and Goldreich(1984)]Fillmore:1984wk
Fillmore, J.A.; Goldreich, P.
Self-similiar gravitational collapse in an expanding universe.
Astrophys. J. 1984, 281, 1–8, doi:10.1086/162070.

[Bertschinger(1985)]Bertschinger:1985pd
Bertschinger, E.
Self-similar secondary infall and accretion in an Einstein-de
  Sitter universe.
Astrophys. J. Suppl. 1985, 58, 39, doi:10.1086/191028.

[Gurevich and Zybin(1995)]1995PhyU...38..687G
Gurevich, A.V.; Zybin, K.P.
Large-scale structure of the Universe. Analytic theory.
Phys. Uspekhi 1995, 38, 687–722, doi:10.1070/PU1995v038n07ABEH000094.

[Baghram et al.(2011)Baghram, Afshordi, and Zurek]Baghram:2011is
Baghram, S.; Afshordi, N.; Zurek, K.M.
Prospects for Detecting Dark Matter Halo Substructure with Pulsar
  Timing.
Phys. Rev. 2011, D84, 043511,


  doi:10.1103/PhysRevD.84.043511.

[Kashiyama and Oguri(2018)]Kashiyama:2018gsh
Kashiyama, K.; Oguri, M.
Detectability of Small-Scale Dark Matter Clumps with Pulsar Timing
  Arrays. arXiv 2018,  arXiv:1801.07847.



[Dror et al.(2019)Dror, Ramani, Trickle, and Zurek]Dror:2019twh
Dror, J.A.; Ramani, H.; Trickle, T.; Zurek, K.M.
Pulsar Timing Probes of Primordial Black Holes and Subhalos. arXiv 2019,  	arXiv:1901.04490.



[Correa et al.(2015)Correa, Wyithe, Schaye, and
  Duffy]Correa:2014xma
Correa, C.A.; Wyithe, J.S.B.; Schaye, J.; Duffy, A.R.
The accretion history of dark matter haloes—I. The physical
  origin of the universal function.
Mon. Not. R. Astron. Soc. 2015, 450, 1514–1520,


  doi:10.1093/mnras/stv689.

[Hayashi et al.(2003)Hayashi, Navarro, Taylor, Stadel, and
  Quinn]Hayashi:2002qv
Hayashi, E.; Navarro, J.F.; Taylor, J.E.; Stadel, J.; Quinn, T.R.
The Structural evolution of substructure.
Astrophys. J. 2003, 584, 541–558,


  doi:10.1086/345788.

[Behroozi et al.(2013)Behroozi, Wechsler, and Wu]Behroozi:2011ju
Behroozi, P.S.; Wechsler, R.H.; Wu, H.Y.
The Rockstar Phase-Space Temporal Halo Finder and the Velocity
  Offsets of Cluster Cores.
Astrophys. J. 2013, 762, 109,


  doi:10.1088/0004-637X/762/2/109.

[Diemand et al.(2006)Diemand, Kuhlen, and Madau]Diemand:2006ey
Diemand, J.; Kuhlen, M.; Madau, P.
Early supersymmetric cold dark matter substructure.
Astrophys. J. 2006, 649, 1–13,


  doi:10.1086/506377.

[Seljak(2000)]Seljak:2000gq
Seljak, U.
Analytic model for galaxy and dark matter clustering.
Mon. Not. R. Astron. Soc. 2000, 318, 203,


  doi:10.1046/j.1365-8711.2000.03715.x.

[Zhu et al.(2016)Zhu, Marinacci, Maji, Li, Springel, and
  Hernquist]Zhu:2015jwa
Zhu, Q.; Marinacci, F.; Maji, M.; Li, Y.; Springel, V.; Hernquist, L.
Baryonic impact on the dark matter distribution in Milky Way-sized
  galaxies and their satellites.
Mon. Not. R. Astron. Soc. 2016, 458, 1559–1580,


  doi:10.1093/mnras/stw374.

[Errani et al.(2017)Errani, Peñarrubia, Laporte, and
  Gómez]2017MNRAS.465L..59E
Errani, R.; Peñarrubia, J.; Laporte, C.F.P.; Gómez, F.A.
The effect of a disc on the population of cuspy and cored dark
  matter substructures in Milky Way-like galaxies.
Mon. Not. R. Astron. Soc. 2017, 465, L59–L63,


  doi:10.1093/mnrasl/slw211.

[Chua et al.(2017)Chua, Pillepich, Rodriguez-Gomez,
  Vogelsberger, Bird, and Hernquist]2017MNRAS.472.4343C
Chua, K.T.E.; Pillepich, A.; Rodriguez-Gomez, V.; Vogelsberger, M.;
  Bird, S.; Hernquist, L.
Subhalo demographics in the Illustris simulation: Effects of baryons
  and halo-to-halo variation.
Mon. Not. R. Astron. Soc. 2017, 472, 4343–4360,


  doi:10.1093/mnras/stx2238.

[Kelley et al.(2018)Kelley, Bullock, Garrison-Kimmel,
  Boylan-Kolchin, Pawlowski, and Graus]Kelley:2018pdy
Kelley, T.; Bullock, J.S.; Garrison-Kimmel, S.; Boylan-Kolchin, M.; Pawlowski,
  M.S.; Graus, A.S.
Phat ELVIS: The inevitable effect of the Milky Way's disk on its
  dark matter subhaloes.  arXiv 2018, arXiv:1811.12413.



[Okoli and Afshordi(2016)]Okoli:2015dta
Okoli, C.; Afshordi, N.
Concentration, Ellipsoidal Collapse, and the Densest Dark Matter
  haloes.
Mon. Not. R. Astron. Soc. 2016, 456, 3068–3078,


  doi:10.1093/mnras/stv2905.

[Gao et al.(2012)Gao, Navarro, Frenk, Jenkins, Springel, and
  White]Gao:2012tc
Gao, L.; Navarro, J.F.; Frenk, C.S.; Jenkins, A.; Springel, V.; White, S.D.M.
The Phoenix Project: The Dark Side of Rich Galaxy Clusters.
Mon. Not. R. Astron. Soc. 2012, 425, 2169,


  doi:10.1111/j.1365-2966.2012.21564.x.


















