3D simulations with boosted primordial power spectra
and ultracompact minihalos

Mateja Gosenca,¹,* Julian Adamek,²,† Christian T. Byrnes,¹,‡ and Shaun Hotchkiss¹,³,§

¹Astronomy Centre, School of Mathematical and Physical Sciences, University of Sussex,
Brighton BN1 9QH, United Kingdom

²Laboratoire Univers et Théories, Observatoire de Paris—PSL Research University—CNRS—Université Paris Diderot—Sorbonne Paris Cité, 5 place Jules Janssen, 92195 Meudon CEDEX, France

³Department of Physics, University of Auckland, Private Bag 92019, Auckland, New Zealand
(Received 24 October 2017; published 18 December 2017)

We perform three-dimensional simulations of structure formation in the early Universe, when boosting the primordial power spectrum on ~kpc scales. We demonstrate that our simulations are capable of producing power-law profiles close to the steep $\rho \propto r^{-9/4}$ halo profiles that are commonly assumed to be a good approximation to ultracompact minihalos (UCMHs). However, we show that for more realistic initial conditions in which halos are neither perfectly symmetric nor isolated the steep power-law profile is disrupted, and we find that the Navarro-Frenk-White profile is a better fit to most halos. In the presence of background fluctuations, even extreme, nearly spherical initial conditions do not remain exceptional. Nonetheless, boosting the amplitude of initial fluctuations causes all structures to form earlier and thus at larger densities. With a sufficiently large amplitude of fluctuations, we find that values for the concentration of typical halos in our simulations can become very large. However, despite the signal coming from dark matter annihilation inside the cores of these halos being enhanced, it is still orders of magnitude smaller compared to the usually assumed UCMH profile. The upper bound on the primordial power spectrum from the nonobservation of UCMHs should therefore be reevaluated.

DOI: 10.1103/PhysRevD.96.123519

# I. INTRODUCTION

The cosmic microwave background (CMB) has yielded precise constraints on the primordial power spectrum of density perturbations that gave rise to all of the gravitationally bound structures observed in the Universe today [1]. These constraints are especially tight over about 2 orders of magnitude in scales, with the constraints on the largest scales being weakened by cosmic variance limits, which cannot be overcome. Hence, the search for new information about the primordial power spectrum must focus on the smaller-scale perturbations for which there is potentially enormously more information to gain. There are numerous challenges to measuring the power spectrum on scales $k \gtrsim 0.3 \text{ hMpc}^{-1}$), e.g. Silk damping of the CMB, radiation pressure during structure formation, and the nonlinear dynamics of matter in the late Universe.

A promising constraint comes from the nondetection of primordial black holes (PBHs), which puts upper limits on the primordial power spectrum over a vast range of scales spanning 20 orders of magnitude [2,3]. However, the constraints are weak, around $P^\zeta \lesssim 10^{-3} - 10^{-2}$, 6 orders of magnitude above the amplitude observed on CMB

scales. While these upper limits do provide useful constraints on certain models of inflation [4,5], they are very weak because the critical threshold to collapse to a PBH is extremely large and only a curvature perturbation of order unity will do so (with an exception if there was an early matter-dominated era, during which time the critical collapse threshold was greatly reduced by the absence of pressure [6,7]).

If the initial perturbations had a higher amplitude than expected from extrapolating the CMB measurement, then another possible signature would be ultracompact minihalos (UCMHs), which are dense dark matter halos with a very steep density profile and large central density. It has been claimed that the nondetection of annihilating dark matter signatures from such compact objects provides relatively stringent upper limits on the primordial power spectrum, around $P^\zeta \lesssim 10^{-7} - 10^{-6}$ over the scales corresponding to $k$ from 10 to $10^7 \text{ Mpc}^{-1}$ [8]. Because of the exponential dependence of the number density of UCMHs and PBHs on the amplitude of the power spectrum, confirmation of the UCMH constraints would also completely rule out the existence of any PBHs forming on these scales (corresponding to a PBH mass greater than about a solar mass). This would rule out the observed supermassive black holes having a primordial origin [9]. The relation between PBH and UCMH formation has also been studied in Refs. [10–12].

<footnote>*M.Gosenca@sussex.ac.uk</footnote>

<footnote>†Julian.Adamek@obspm.fr</footnote>

<footnote>‡C.Byrnes@sussex.ac.uk</footnote>

<footnote>§S.Hotchkiss@auckland.ac.nz</footnote>
Assuming that weakly interacting massive particles (WIMPs) make up the majority of dark matter, the most stringent observational constraint on the existence of UCMHs comes from the expected annihilation signal from the cores of UCMHs, where the density is extremely large and hence the probability of WIMPs annihilating is massively boosted by their large number density. The constraints only weakly depend on the assumed cross section and mass of the WIMPs but can be totally evaded if dark matter does not annihilate. There are, nonetheless, other possible signatures of UCMHs, such as from the Shapiro time delay of millisecond pulsars (although note that no actual constraints exist yet, but there are forecasts [13,14]) or microlensing [15–17].

It is in any case, interesting to ask how the dark matter substructure would change if the initial small-scale perturbations were larger than expected. Other papers that constrain the power spectrum on small scales include Refs. [18–35]. In particular, we note that CMB spectral distortion constraints provide an upper bound on the primordial power spectrum on about the same scales as UCMHs, which are about an order of magnitude weaker than the currently claimed constraints but independent of the dark matter (DM) model [20].

To date, the theoretical forecasts on UCMH formation have been made using an idealized analytical model of their formation. Josan and Green [36] assumed that any density contrast $\delta = \delta\rho/\bar{\rho}$ satisfying $\delta_c > 10^{-3}$ at horizon entry would form an UCMH (following Ref. [15]) with a steep power-law profile $\rho \propto r^{-9/4}$. This is motivated by the analytical calculation of Bertschinger 1985 [37], who calculated the profile of a spherically symmetric perturbation accreting from a homogeneous background. The analytical calculation of the requisite density threshold was refined by Bringmann et al. [8], who included the effects of radiation in their calculation, which is relevant because UCMHs are assumed to form shortly after matter-radiation equality. They found comparably tight constraints. We show that in practice no $\delta_c$ actually exists, because the formation and final profile of a small halo strongly depend on its environment in a way that cannot be captured even approximately by a single number, unlike the case for PBHs.

In this paper, we revisit and dramatically revise the calculation of UCMH formation and for the first time perform realistic 3D simulations of this complex process. The analytical calculations assumed that UCMHs are sufficiently rare high-density peaks that form in isolation, surrounded by an unperturbed background until redshift ~10. They also neglected the effects of angular momentum, thereby using spherical symmetry to reduce the problem to a one-dimensional one. This gives rise to the steep $\rho \propto r^{-9/4}$ density profile, with a correspondingly boosted dark matter annihilation signal. We verify these results with

an N-body simulation of an isolated spherical overdensity, obtaining a UCMH-like density profile.

In practice, however, all UCMHs will be formed with particles with nonzero initial peculiar velocities, and they are not really isolated from nearby perturbations, even if the nearest perturbations have a smaller amplitude. We analytically estimate that even an exceptionally large-amplitude perturbation will be relatively close (in comparison to the size of the object) to another perturbation with the same scale and at least half its amplitude, even in the case of non-Gaussian initial perturbations. Moreover, even if an overdense region appears spherical when smoothed on a particular scale, it will still have fluctuations on smaller scales. These fluctuations will themselves grow, making the initial sphericity unstable [38].

We then perform N-body simulations using two classes of initial conditions: the first with a 1 kpc/h sized proto-halo of large initial density surrounded by successively larger surrounding fluctuations and the second starting with the standard featureless power-law spectrum of the curvature perturbation but boosted around the kpc scale, with boost factors of 10, 100, and 1000 and without the boost [corresponding to the standard $\Lambda$ cold dark matter (CDM) cosmology].

We find no clear evidence that UCMHs are formed even in the case in which we boost the initial power spectrum by a factor of 1000. Furthermore, a UCMH that formed when simulated in isolation becomes unexceptional compared to other halos in the simulation box when we add background perturbations with a typical amplitude one-fifth of that of the protohalo. This means that the analytical estimates of UCMH formation are a poor match to realistic initial conditions. Instead, we find that the emergence of Navarro-Frenk-White (NFW) profiles is more generic [39–41]. NFW profiles have a central core with a much gentler $\rho \propto 1/r$ slope and a correspondingly much weaker WIMP-annihilation signal.

The NFW minihalos formed when the initial power spectrum is boosted are nonetheless of much greater density than would exist in the standard $\Lambda$CDM Universe. Because every structure that forms earlier forms at a greater density, these NFW minihalos are also much more abundant than proposed UCMHs, which would only have formed at the most extreme locations of the primordial density field. Therefore, it remains interesting to study to which observational signatures they may give rise. We provide estimates of how sharply the existing constraints have to be weakened when taking our results into account by calculating the WIMP-annihilation signal from an NFW rather than the $r^{-9/4}$ profile.

This paper is structured as follows. In the next section, we introduce the power-law density profile usually assumed for the UCMHs and the universal NFW profile and explain the fitting procedure we used to find the parameters of these models. We also describe the initial
conditions of large-amplitude fluctuations. In Sec. III, we discuss the simulation results from an initial spherical overdensity, of which the evolution into an UCMH in a homogeneous background we disrupt by adding random perturbations around it, with successively larger amplitude. We then boost the power spectrum on the kpc scale in Sec. IV by up to a factor of 1000, and we describe how this gives rise to many NFW-like halos with very large concentration parameters. In Sec. V, we briefly review the theory of WIMP annihilation in the dense center of halos and calculate the signal from the halos we simulated, before concluding in Sec. VI. Technical details of how we performed the simulations and convergence testing are left to the Appendixes.

## II. HALO PROFILES AND PROPERTIES

### A. UCMH profiles

Ultracompact minihalos are dark matter halos, expected to form around matter-radiation equality, featuring a very steep power-law density profile [15],

$$ \rho(r) \propto r^{-9/4}. \qquad (1) $$

If a UCMH forms, its extreme compactness might allow it to retain its shape until the present time, making it and the extreme density at its center potentially observable.

This $r^{-9/4}$ power-law profile is the late-time form for the density of an initial spherical overdensity accreting from a homogeneous background [37].<fnref n="1" /> This is true irrespective of the density profile of the initial seed overdensity, so long as its size is finite, it is spherical, and it eventually starts accreting from a homogeneous background. However, note that even the smallest amount of asphericity in initial conditions is unstable and will result in a triaxial profile [38].

For fitting to the power-law profiles, we use the ansatz

$$ \rho(r) = Cr^{-\alpha}. \qquad (2) $$

If we know the mass within a certain radius, then C can be calculated via

$$ M(r_{\text{max}}) = \int_{0}^{r_{\text{max}}} 4\pi r^2 dr \rho(r). \qquad (3) $$

The choice of $r_{\text{max}}$ is somewhat ambiguous because halos do not have a clear “edge”; their density instead asymptotes to the average density of the Universe, and for $\alpha \le 3$, the above integral diverges with $r_{\text{max}}$. There are several

possible choices for $r_{\text{max}}$ in the literature; we use the virial radius $r_{\text{vir}}$ as determined by the ROCKSTAR halo finder [42]. It is defined as the radius of a sphere inside which the average density contrast is $18\pi^2 \approx 178$ [43]. Correspondingly, the mass of a halo up to the virial radius $M_{\text{vir}} = M(r_{\text{vir}})$ is the virial mass.

For the above power-law profile, C can therefore be expressed in terms of the virial mass and virial radius as

$$ C = \frac{M_{\text{vir}}}{4\pi} (3 - \alpha) r_{\text{vir}}^{(\alpha-3)}. \qquad (4) $$

### B. Navarro-Frenk-White profile

Most of the halos in the Universe appear to exhibit a density profile close to the NFW [44,45] profile,

$$ \rho(r) = \frac{\rho_0}{\frac{r}{r_s}(1 + \frac{r}{r_s})^2}, \qquad (5) $$

characterized by two parameters: $r_s$, the “scale radius,” and $\rho_0$. For small radii, $r \ll r_s$, the profile’s radial dependence is $\rho(r) \propto 1/r$, and for large $r$, it is $\rho(r) \propto 1/r^3$. The scale radius $r_s$ therefore determines the radius at which the profile changes from one power law to the other. The characteristic density $\rho_0$ corresponds to $\rho_0 = 4\rho(r_s)$. Integrating Eq. (3) for $\rho(r)$, we obtain the mass within radius $r_{\text{max}}$:

$$ M(r_{\text{max}}) = 4\pi\rho_0 r_s^3 \left( \ln\left(\frac{r_s + r_{\text{max}}}{r_s}\right) - \frac{r_{\text{max}}}{r_s + r_{\text{max}}} \right). \qquad (6) $$

Defining the concentration parameter [46] as

$$ c = \frac{r_{\text{vir}}}{r_s}, \qquad (7) $$

$\rho_0$ can be expressed as

$$ \rho_0 = \frac{M_{\text{vir}}}{4\pi r_s^3 (\ln(1+c) - c/(1+c))}. \qquad (8) $$

### C. Fitting the profiles

We can use the relation between $M_{\text{vir}}$, $r_{\text{vir}}$ and the two profiles’ parameters [i.e., Eqs. (4) and (8)] to eliminate one free parameter from each profile. In logarithmic space, the NFW profile can be expressed as

$$ \log(\rho(r)) = \log\left(\frac{M_{\text{vir}}}{4\pi r_{\text{vir}}^3}\right) + \log\left(\frac{c^2}{\ln(1+c) - c/(1+c)}\right) - \log(r/r_{\text{vir}}) - 2\log(1 + cr/r_{\text{vir}}), \qquad (9) $$

and the power law can be expressed as

<footnote n="1">¹Note that this is only strictly true when the profile is sufficiently coarse grained. The true profile will be made up of a series of separate caustics arising from shells that have passed through the center of the overdensity a different number of times. See Fig. 8 of Ref. [37] for the shape of the full non-coarse-grained profile.</footnote>
$$ \log(\rho(r)) = \log\left(\frac{M_{\text{vir}}}{4\pi r_{\text{vir}}^3}\right) + \log(3-\alpha) - \alpha \log(r/r_{\text{vir}}). \quad (10) $$

For $M_{\text{vir}}$ and $r_{\text{vir}}$, we use the values calculated with ROCKSTAR.

To perform the fit, we then logarithmically weight every point in the profile. Specifically, the function we minimize when fitting to a profile is

$$ S = \sum_i w_i (\log(\rho_{\text{model}}(r_i)) - \log(\rho_{\text{data}}(r_i)))^2, \quad (11) $$

with the weights $w_i = \log(r_i/r_{i-1})$ chosen so that each point is weighted according to the logarithm of the size of the bin it represents. This would be equivalent to having bins of equal length in $\log(r)$ with equal weighting. This fitting was chosen because we fit over more than an order of magnitude in $r$ and the density within this range also changes by more than an order of magnitude.

The fitting is performed over the range from $r = 0.004 \text{ kpc/h}$ to $r = r_{\text{vir}}$, independently for each halo. The lower bound is forced upon us by the resolution of the profiles (see Sec. A 4), and the upper bound is where we assume the “edge” of the halo to be.

For the NFW fits, we can also plot rescaled density $\rho/\rho_0$ against rescaled radius $r/r_s$. In these units, all halos with an NFW profile “collapse” onto a single curve.

One selection criterion we apply is that we only consider halos that are not in the process of merging. To determine this, we compare two measures of the halo mass. The first is the combined mass of the particles that ROCKSTAR associates to that halo given their six-dimensional phase-space information. The second one includes all particles, even the ones of which the momenta suggest that they are not part of the halo. We only keep halos for which these two mass parameters agree within 10%.

## D. Extreme fluctuations

To date, the literature on UCMHs has focused on treating these halos as isolated objects, forming in an otherwise homogeneous and uniform background. This allows analytical calculations to be made, but the validity of this approximation has not been verified.

To answer how extreme the most dense region of a certain size in some larger volume is, and assuming the density perturbation $\delta$ at horizon entry to be distributed according to the Gaussian distribution with $\bar{\delta} = 0$, we can express the fraction of volume where $\delta > x$ as

$$ \begin{aligned} \frac{V_{\delta>x}}{V} &= \frac{1}{\sigma\sqrt{2\pi}} \int_x^\infty e^{\left(-\frac{\delta^2}{2\sigma^2}\right)} d\delta \\ &= \frac{1}{2} \text{Erfc}\left(\frac{x}{\sqrt{2}\sigma}\right). \end{aligned} \quad (12) $$

<figure><img src="image_1.png" /><figcaption>FIG. 1. A histogram showing how many standard deviations the maximum $\delta_{1 \text{ kpc/h}}$ is away from the average in a 32 kpc/h-sized volume (taken from 1000 samples). To derive this histogram, we have made the approximation that each 1 kpc/h region in each 32 kpc/h volume is independent of all the other 1 kpc/h regions in that volume.</figcaption></figure>

If we further consider a region of the size of $(1 \text{ kpc/h})^3$ in a $(32 \text{ kpc/h})^3$ volume, we find that the most extreme overdensity is a $4.0-\sigma$ fluctuation. We show the distribution of the most extreme values in Fig. 1, which shows that a $4-\sigma$ fluctuation is the most likely value but that values as extreme as $5.5-\sigma$ are possible if one draws from 1000 samples.

We can also reverse the question and estimate the fraction of volume in which the density is greater than some threshold value. For two and three times the typical density contrast, we get this fraction to be 0.023 and 0.001, respectively. In a $(32 \text{ kpc/h})^3$ volume, this corresponds to 745 and 44 regions of $(1 \text{ kpc/h})^3$ size. Assuming that all peaks are randomly placed independently of each other, an estimate of the distance between the most overdense region and the one with overdensity greater than or equal to twice what is typical is given by $32/(745)^{1/3} \text{ kpc/h} \approx 3.5 \text{ kpc/h}$, i.e., rather close compared to the scale of the halo itself. Because of the clustering of density peaks, which arises from superimposing multiple scales, the distances between large overdensities will in fact be even smaller than the estimate made above [47].

The spherical infall model, which leads to a steep power-law profile also assumes spherical symmetry. While it is true that the most extreme peaks are expected to be close to spherical initially [48], the unsmoothed density field will still have substructure. This substructure will grow with time and pull any matter on a purely radial trajectory off this trajectory. Moreover, gravitational tidal forces from the nearby surrounding halos will further break the spherical symmetry.

In the presence of large primordial non-Gaussianities (which are observationally unconstrained on the small scales we are interested in), one might expect the existence
of isolated large-amplitude density peaks to become much more probable. For some probability density functions, extremely large-amplitude peaks do become exponentially more likely, but so do other comparably large-amplitude peaks, implying that the extremely large peak will still not be isolated. For example, if the density perturbation is drawn from a chi-squared probability density function with 1 degree of freedom, $\delta = \delta_G^2 - \langle \delta_G^2 \rangle$, where $\delta_G$ has a Gaussian distribution (which follows in the asymptotic limit of very large local $f_{NL}$ [49]), then the typical $(32 \text{ kpc/h})^3$ volume will have a $10-\sigma$ fluctuation, and the most extreme fluctuation in 1000 realizations of the initial conditions will have a $20-\sigma$ fluctuation, twice the amplitude of the second largest perturbation expected in the same volume. Furthermore, the mode coupling between different scales would in practice mean that the largest overdensity is likely to be situated close to the second-largest perturbation, and so it still appears to be very hard to come up with a situation in which a large-amplitude perturbation is likely to form in isolation from its nearest neighbours. We will show in Sec. III B that a UCMH-like halo can form if we place a smooth, spherical overdensity surrounded by perturbations with a typical amplitude 15 times smaller. It would be interesting to find a model of the early Universe capable of generating comparable initial conditions.

### III. GAUSSIAN-PROFILE SPHERICAL OVERDENSITY

#### A. Completely isolated halo

We first simulate an isolated, spherically symmetric overdensity to test the theoretical prediction of a power-law profile (1). For the initial profile of the overdensity, we chose a three-dimensional function with a Gaussian profile

$$ \delta_{\text{in}}(r) = A \exp\left(\frac{-r^2}{2\sigma_w^2}\right), \qquad (13) $$

where A represents the initial amplitude of the seed overdensity in the center and $\sigma_w$ regulates the size of the overdensity. We place the overdensity at the center of the box, and r denotes the distance from this center. We have set $\sigma_w = 0.5 \text{ kpc/h}$ so that the overdensity is well contained inside the simulation volume. We set $A = 0.3$ at initial redshift $z_{in} = 10000$ such that the perturbation starts off marginally within the linear regime but quickly triggers the spherical collapse process shortly after matter-radiation equality. After some time, the collapse dynamics give rise to a self-similar steady state, and we determine the slope of the density profile by fitting

$$ \rho(r) = Cr^{-\alpha}, \qquad (14) $$

which in log-log space is just a straight line with the slope $-\alpha$. In Fig. 2, we show that we indeed obtain a $\rho(r) \sim r^{-\alpha}$

<figure><img src="image_2.png" /><figcaption>FIG. 2. The profile of the Gaussian-seed halo in the absence of all other perturbations. We fit the power law and thus obtain the parameter $\alpha$. The number of particles used for this simulation was $512^3$. By $z=3$, the halo mass has grown to $1.8 \times 10^4 M_{\odot}/h$, which is still a negligible fraction of the total mass inside the box (meaning that the finite box size is not slowing down the halo accretion).</figcaption></figure>

profile. It is not precisely $\alpha = -9/4$, but this is not entirely unexpected because of numerical instabilities [38,50] and the fact that $\alpha = -9/4$ is only the coarse-grained form of this type of structure growth.

Collisionless particles that can freely move through the center of the overdensity form caustics, which can be seen as additional structure at the edges of profiles. This effect is physical, and we have tested that by changing the resolution of the simulation, as shown in Fig. 20.

#### B. Peak-to-background ratio 15

Having seen that we can form a UCMH starting with a spherically symmetric overdensity surrounded by a homogeneous background, we will now investigate how the evolution of the central halo changes when we drop the assumption of a homogeneous background. Following the discussion in Sec. II D, we would not expect the amplitude of the central halo to be more than about five times larger than the typical amplitude of the other perturbations. We study that case in Sec. III C. Here, we study as an intermediate case the situation in which the amplitude of the typical fluctuations in the box smoothed on the same kpc/h scale was 15 times smaller.

Given our choice of $A = 0.3$, we can make the typical amplitude of the perturbations 15 times smaller by boosting the standard $\Lambda$CDM power spectrum by a factor of 16 on all scales relevant to our simulation, meaning that the linearized root-mean-square perturbation (excluding the large-amplitude central halo) is 4.0 times larger than it would be in $\Lambda$CDM on all scales. Starting from such initial conditions, we simulate and study the formation and evolution
<figure><img src="image_3.png" /><figcaption>FIG. 3. Slice through the simulation box at redshift $z = 30$ for the peak-to-background ratio 15 simulation. The halo that was seeded by the Gaussian overdensity peak is slightly to the left and to the top of the center, highlighted in light blue.</figcaption></figure>

of halos. A snapshot at redshift $z = 30$ is illustrated in Fig. 3.

To determine $\alpha$, the slope of the profile in log-log space, we fit the linear function (14) to the central halo.<fnref n="2" />

In Fig. 4, we show a scatter plot of the slope parameter $\alpha$ against the virial mass for every halo in this simulation (at $z = 100$ and $z = 30$). The red circle in each panel represents the seed halo that started with amplitude 15 times the background. The seed halo clearly has an exceptional slope compared to its background, even down to $z = 30$.

In Fig. 5, we show the density profiles of halos at multiple redshifts. For $z = 100$, we show a random subsample of halos with a mass $M_{vir} > 10^2 M_{\odot}/h$, and for the other three redshifts, the mass cut is $M_{vir} > 10^3 M_{\odot}/h$. We also include a histogram of $\alpha$ values for all of the halos that pass this mass cut. It is again clear that among the most massive halos the seed halo's slope remains exceptional, even down to $z = 10$ (although less and less so as the background structure develops).

Finally, in Fig. 6, we collapse the halo profiles on to the NFW function in the manner described at the end of Sec. II C (we apply the same mass cuts as for Fig. 5). We

also show a histogram of the difference in the value of our fitting statistic $S$, Eq. (11) for the best-fitting power-law and NFW profiles (i.e., $S_{PL} - S_{NFW}$). It is clear from the profiles that the NFW profile is not an excellent fit for the special seed halo. There is a noticeable upturn at small $r/r_s$. Also note how far to the right the special halo's profile extends in each panel. This is a consequence of the NFW fit requiring a small $r_s$ value, exceptionally smaller than any other halo that satisfies the mass cuts.

We cannot resolve the profile of any halo on scales below the resolution limit of our simulation. However, because the special halo remains exceptional compared to its background, even down to $z = 10$, it seems reasonable to expect that the profile would remain close to a power law even beyond the limits of our resolution.

### C. Peak-to-background ratio 5

We now increase the amplitude of the surrounding perturbations even more and repeat the analysis. In the results of this subsection, the central overdensity is only five times larger than a typical one. A snapshot at redshift $z = 30$ is shown in Fig. 7.

Figures 8, 9, and 10 are analogous to Figs. 4, 5, and 6, however, for the simulations with even larger background perturbations and showing only $z = 100$ and $z = 30$.

The most striking observation to make is that the special halo is no longer at all exceptional compared to its surroundings. Its slope compared to its mass is large compared to other halos in the simulation, but not exceptionally so. Moreover, its value of $S_{PL} - S_{NFW}$, although still favoring a power law, is not at all special compared to other halos, some of which also favor the power law.

We cannot say what happens below the limits of our resolution; however it is clear that when an initial, spherically symmetric perturbation is embedded in a background five times smaller it is no longer exceptional at least by $z = 100$. Therefore, it does not behave as if it were in a homogeneous background, and it is reasonable to assume that it has not formed a UCMH-like profile on smaller scales. Although we do not show it here, this holds even at $z = 300$.

### D. Summary of special halo simulations

In this section, we have shown that as the background of an initial spherically symmetric perturbation is increased the halo descending from this seed perturbation becomes less and less exceptional. This might not appear too surprising; however, it is not trivial that an initial $5-\sigma$ fluctuation at $z = 10000$ will be entirely unexceptional by $z = 100$.

The conclusion we take from this section is that unless a fluctuation has an initial amplitude at least between 5 and 15 times larger than the typical background fluctuations its halo will soon become comparable in slope, mass, and density to many other nearby halos. Unless the perturbation

<footnote n="2">²Since halos move around the simulation volume, it is not always obvious which halo formed from the Gaussian overdensity. We identified that halo by tracking one of the particles that was initially closest to the center of the simulation volume. At each redshift, we looked for the halo that contained the tracked particle at its core.</footnote>
<figure><img src="image_4.png" /><figcaption>FIG. 4. The best-fit power-law exponent $\alpha$ for all halos in the peak-to-background 15 simulation at $z = 100$ and 30. The special halo has an exceptionally steep slope compared to all other halos with a comparable mass.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 5. The profiles of the halos in the peak-to-background ratio 15 simulation at $z = 100, 30, 15,$ and 10. The profile of the halo descending from the seed halo is highlighted in red so it can be distinguished from the other halos. For each snapshot, we fit the profiles to Eq. (14) and present the resulting $\alpha$ in a histogram. The range we fit over is $r$ between 0.004 kpc/h and $r_{vir}$ for each halo individually. To plot the profiles, for $z = 100$, we only kept halos with $M_{vir} > 10^2 M_{\odot}/h$, and for $z = 30, 15, 10$, we kept the halos with $M_{vir} > 10^3 M_{\odot}/h$. The power-law fit to the seed halo's profile is plotted as a violet dashed line.</figcaption></figure>
<figure><img src="image_6.png" /></figure>

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

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

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

FIG. 6. The rescaled profiles of halos in the peak-to-background ratio 15 simulation. The NFW analytical prediction is shown in green, and the descendant of the seed halo is highlighted in red. Histograms show the difference of a measure for the goodness of fit between the power-law and NFW profile. A negative value of this quantity points toward the power law being a better fit than NFW. With the exception of $z = 100$, the special halo seems to favor the power-law profile over NFW. The exceptionality is particularly noticeable near the center of the halo. Fitting the special halo with an NFW profile gives a very small $r_s$. Here, this is manifested by the special halo's profile being shifted toward larger radii than any other halo.

is more extreme initially than one would expect to find given Gaussian initial conditions, see Sec. II D, it will not grow into a UCMH-like halo.

We stress that this is not the result of lowering the density of the initial seed fluctuation. The initial profile was identical in all three subsections; the only thing that was changed was the amplitude of background fluctuations. Therefore, when examining the evolution of objects on the very small scales in the early Universe, a single density contrast $\delta_c$ is not enough to describe the subsequent evolution. The evolution depends sensitively on the environment, and therefore any constraints must, too.

## IV. BOOSTING THE POWER SPECTRUM AROUND THE 1 kpc/h SCALE

In this section, we no longer consider the evolution of a specific overdensity peak that was planted *by hand*. Instead,

we increase the probability of extreme *random* fluctuations on a similar scale by boosting the variance of the primordial perturbation modes over some range of scales. In principle, this boost can take many forms (e.g., a step or a bump in Fourier space), and the precise nature of structure formation will depend on the form of the boost. However, UCMHs are claimed to form with a sufficiently large abundance to be observed under *any* boost of the power spectrum with sufficiently large amplitude.

We see no halo with convincing UCMH properties; however, we are somewhat limited by the resolution. We do, however, see many compact structures. With hindsight, this is not surprising either. If we boost the primordial spectrum, then structures on the boosted scales form earlier. This means these structures form when the Universe is more dense. Therefore, they also reach a larger virial density than they would had they formed from a non-boosted initial power spectrum.
<figure><img src="image_10.png" /><figcaption>FIG. 7. Slice through the simulation box at redshift $z = 30$ for the peak-to-background ratio 5 simulation. The halo that was seeded by the Gaussian overdensity peak is slightly to the left and to the top of the center, highlighted in light blue.</figcaption></figure>

Note, however, this is not only true for the most extreme, and therefore rare, structures forming under a boosted power spectrum. In fact every structure that forms on the boosted scales forms earlier and is therefore more compact. As a consequence, even though we do not find structures as compact as a hypothetical UCMH, the structures we find may, due to their increased abundance, have their own unique cosmological signals. It would be useful and interesting future work to examine how the new type of

compact halos we describe depends on the specific nature of the boost to the power spectrum.

The form of our boost in the power spectrum is as follows. First, we take the following unboosted power-law power spectrum of the primordial curvature perturbation $\zeta$, $\mathcal{P}_0^\zeta(k) = A_s(k/k_{\text{pivot}})^{n_s-1}$ with $A_s = 2.26 \times 10^{-9}$, $k_{\text{pivot}} = 0.05$ Mpc$^{-1}$, and $n_s = 0.96$.

We then boost this power spectrum to form the primordial power spectrum used for our simulations, $\mathcal{P}^\zeta(k)$. Specifically,

$$ \mathcal{P}^\zeta(k) = \mathcal{P}_0^\zeta(k)(1 + B(k)), \qquad (15) $$

where

$$ B(k) = A_b \exp(-2.77(\ln k - \ln k_*)^2). \qquad (16) $$

The scale at which we boost the simulation's power spectrum is $k_* = 1$ h kpc$^{-1}$.

The value of 2.77 in the definition of $B(k)$ is chosen such that the full width at half maximum of the boost is 1 (in units of $\ln k$). This means that our power spectrum is boosted over a width corresponding to approximately one e-folding, which is a natural length scale during inflation.

Before starting each simulation, we use CAMB [51] to calculate the linear transfer functions and hence the matter power spectrum at the simulation's starting redshift. The shape of such a boosted power spectrum is plotted in Fig. 11 for $A_b = 10^3$, $k_* = 1$ h kpc$^{-1}$ and at redshift $z = 100000$. Note that this is the power spectrum of dark matter only. The power spectrum of baryons is significantly different due to baryons still being tightly coupled to photons at this high redshift.

We also show $\sigma_R$ (i.e., the rms of the density contrast smoothed by a spherical top hat on a scale $R$) in Fig. 12 for

<figure><img src="image_11.png" /><figcaption>FIG. 8. Figures analogous to Fig. 4, but for the peak-to-background ratio 5. The power-law parameter $\alpha$ of the halo, formed from the special seed, is no longer very different from the other halos.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>FIG. 9. The profiles of the halos in the peak-to-background ratio 5 simulation at $z = 100$ and 30. We have again only included halos with $M_{vir} > 10^2 M_{\odot}/h$ for the $z = 100$ profiles and $M_{vir} > 10^3 M_{\odot}/h$ for $z = 30$. The profile of the descendant of the seed halo is again highlighted in red. Unlike the situation in the peak-to-background ratio 15 simulation, the slope of the special halo is no longer more extreme than the other halos in the box.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>FIG. 10. The rescaled profiles of halos in the peak-to-background ratio 5 simulation. The analytical prediction is shown in green, and the special halo is highlighted in red. In contrast to the peak-to-background ratio 15, the special halo does not appear distinguishable from the rest. The $S_{PL} - S_{NFW}$ measure also demonstrates the profile of the special halo to be a better fit to NFW than to the power law (note a smaller value of $S$ indicates a better fit; see Eq. (11) and the discussion around it).</figcaption></figure>

the same input power spectrum with four different boost amplitudes at the same redshift.

We ran simulations with boosts of $A_b = 0, 10, 100$, and 1000. Each simulation was given the same seed. At $z = 100 000$, the largest fluctuation in the $A_b = 1000$ simulation, when the density was smoothed with a spherical top hat filter of radius 1 kpc/h, was $\delta_R = 0.23$. With this boost, $\sigma_R$ at this radius and redshift is 0.054, and therefore this corresponds to a $4.3-\sigma$ fluctuation.

We start each simulation at $z = 5 \times 10^6$ and run them to $z = 15$. We do not go to lower redshifts than this because at $z = 10$ the scale of our box, 32 kpc/h, becomes nonlinear. Snapshots at redshift $z = 30$ are illustrated in Fig. 13.

A. Results of boosted simulations

As might be expected, given our results from Sec. III, we find that most of the halos fit NFW profiles better than the power law (see Figs. 14 and 15). Moreover, there are no obvious exceptional halos, even when we set $A_b = 10^3$. There are a few halos that fit the power law better than NFW; however, none stands out as much as the special halo in the peak-to-background ratio 15 simulation in Sec. III B (see, for example, the red line in Fig. 6 for $z = 30, 15$, and 10). We expect this better fit arises here simply because our resolution is not good enough to resolve the $r_s$ of the NFW profile for these *lower mass* halos. A similar effect was seen for the nonspecial halos in the previous section.
<figure><img src="image_14.png" /><figcaption>FIG. 11. Boosted matter power spectrum used in our simulations. Plot shown has boost parameters $A_b = 10^3$ and $k_* = 1 \text{ hpc}^{-1}$ from Eq. (16). The matter power spectrum has been calculated using CAMB and is output at $z = 100000$. Note that this is only the dark matter power spectrum.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>FIG. 12. The root-mean-square fluctuations of the density perturbation when smoothed with a top hat filter of radius $R$ (i.e., $\sigma_R$). This is for the same input power spectrum as that shown in Fig. 11 with four different boosts and at the same redshift.</figcaption></figure>

Histograms of the concentration parameter $c$ are shown in Fig. 16. At redshift $z = 30$, concentration parameters of $c \gtrsim 100$ can occur for the $A_b = 10^3$ simulation. Given that we expect this to grow with time, the concentration today would be much larger. This shows that, although there may be no UCMH candidates in our simulations, the structures that do form are still much more compact than those in a simulation without a boosted power spectrum. It is worth stressing again that when the input power spectrum is boosted all the structures that form are much more compact. It is not just the rare, extreme fluctuations that experience this effect. This is because the boosted initial conditions formed structures earlier, when the whole Universe was more dense. Therefore, these structures also

virialize at much larger densities and remain much more dense at later times.

A similar, but less pronounced, effect is also seen for the smaller boosts.

## V. UCMH OBSERVABILITY AND WIMP ANNIHILATION

Numerous observational tests have been either proposed or enacted to constrain the number density of UCMHs, and these observational limits translate into upper bounds on the primordial power spectrum. These constraints assumed a dense $r^{-9/4}$ inner profile of the halos, which we have shown to be unlikely to form. In this section, we briefly calculate how the observability of these halos changes when using an NFW halo with a large concentration parameter, of the form observed in our simulations. While calculating the constraint on the power spectrum goes beyond the scope of this work, we use our simulation results in this section to show how the expected signature of WIMP annihilation in the dense center of the halos reduces when using a realistic halo density profile.

The expected gamma-ray flux from WIMP annihilation within a halo at distance $d$ from the Earth is given by [36]

$$ \Phi_{\gamma} = \frac{\Phi_{\text{astro}} \Phi_{\text{particle}}}{2d^2}, \qquad (17) $$

where $\Phi_{\text{particle}}$ depends on the WIMP particle mass and annihilation cross section, which is assumed to be independent of the halo properties, e.g., the velocity distribution of the particles (see Ref. [52] for a justification).

The astrophysical part is given by an integral of the density squared

$$ \Phi_{\text{astro}} = \int_{0}^{r_{h}} \rho(r)^{2} r^{2} dr, \qquad (18) $$

where $r_h$ is the radius of the halo, which we always take to be $r_{vir}$. In practice, although the majority of the NFW halo mass is within the outer part of the profile where $\rho \propto 1/r^3$, the annihilation signal is dominated by the innermost part with the greatest density, which satisfies $r < r_s$, so it does not matter much where we cut off the integral. Note that the physical radius and density must be used.

For halos with an $r^{-9/4}$ profile, the signal is dominated by the central density, which must be cut off at some maximum value. Bringmann et al. estimate a maximum possible density of the UCMHs today of

$$ \rho_{\text{max}} \approx \frac{m_{\chi}}{\langle \sigma v \rangle (t_0 - t_i)} = K \rho_{c,m}, \qquad (19) $$

where

$$ K \approx 5 \times 10^{16} \qquad (20) $$
<figure><img src="image_16.png" /><figcaption>FIG. 13. Slices through the simulation volume at redshift $z = 30$ for the simulations with boosted power spectrum. The amplitudes of the boosts are $A_b = 0, 10, 100, 1000$, clockwise, starting from the top left.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>FIG. 14. The evolution of halos' profiles with the increasing amount of boost. The NFW analytical prediction is shown in green. Histograms show the difference of a measure for goodness of fit [Eq. (11)] between the power law and NFW. Almost all halos seem to be a better fit to NFW than to the power law.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>FIG. 15. The rescaled profiles for $A_b = 10^3$ at $z = 100, 30$, and 15. From the positions of the lines on top of the NFW reference line, it is clear that $r_s$ is getting smaller as time increases. As with Fig. 14, the NFW profile is a better fit than the power law for almost all of the profiles.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>FIG. 16. Concentration parameter $c$ for the halos in simulations in which we boosted the power spectrum. The halos in higher-boost simulations are much more compact. Concentration also grows with time. We only retain halos with $M_{vir} > 10^2 M_{\odot}/h$. Note that at $z = 15$ there already exist halos in the $A_b = 10^3$ simulation with $c > 200$.</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>FIG. 17. The total $\Phi_{astro}$ part of the WIMP-annihilation signal from the most massive halos plotted against $z$ for the power spectrum boosted by four different amounts. Notice how the strength of the signal typically decreases slowly with redshift, which we caution may be a numerical artifact; see the text after Eq. (26). To calculate this signal, we take into account all the halos with $M_{vir} \gtrsim 3 M_{\odot}$. For the two smallest boosts at $z = 100$, no such halos are identified in the simulation.</figcaption></figure>

calculated assuming fiducial values of the WIMP mass $m_{\chi} = 1$ TeV, thermally averaged cross section $\langle \sigma v \rangle = 3 \times 10^{-26} \text{ cm}^3 \text{ s}^{-1}$, and the age of the Universe $t_0 = 13.7$ Gyr. The UCMH formation time $t_i$ is irrelevant provided that $t_i \ll t_0$ [52]. The critical density of the Universe today is $\rho_c = 415 M_{\odot}h^2/\text{kpc}^3$, and $\rho_{c,m} \equiv \Omega_m\rho_c$. The UCMH profile is given by

$$ \rho_{\text{UCMH}} = \begin{cases} \rho_{\text{max}} & \text{if } r < r_{\text{cut}}, \\ \rho_{\text{max}} \left(\frac{r}{r_{\text{cut}}}\right)^{-\frac{9}{4}} & \text{if } r > r_{\text{cut}}. \end{cases} \qquad (21) $$

Using the definition of the virial mass in terms of the virial radius (which is independent of the density profile),

$$ M_{\text{virial}} = \frac{4\pi}{3} 178 \rho_{c,m} r_{\text{vir}}^3, \qquad (22) $$

where $\rho_{c,m}$ is the critical density of matter (of both CDM and baryons, because our simulations treat them equally). Calculating the density contrast centered on a UCMH to a radius $r_{vir}$,
<figure><img src="image_21.png" /></figure>

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

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

<figure><img src="image_24.png" /><figcaption>FIG. 18. Scatter plots for the WIMP-annihilation signal plotted against $M_{vir}$ for three different redshifts (top and left panels) and against $c$ for $z = 30$ (bottom right). Different amplitudes of the boost are represented with different colors. It is curious that the WIMP-annihilation signal appears to depend on concentration in a way that does not depend on the size of the boost. However, there is a lot of scatter in this relationship.</figcaption></figure>

$$1 + \delta_{\text{UCMH}}(r_h) = \frac{1}{\rho_{c,m}} \frac{3}{r_{\text{vir}}^3} \int_0^{r_{\text{vir}}} \rho_{\text{UCMH}}(r) r^2 dr \\ \approx 4K \left(\frac{r_{\text{cut}}}{r_{\text{vir}}}\right)^{9/4}, \qquad (23)$$

leads to

$$\frac{r_{\text{cut}}}{r_{\text{vir}}} \approx \left(\frac{179}{4K}\right)^{4/9} \ll 1. \qquad (24)$$

Using the above results, we find the WIMP-annihilation signal is

$$\Phi_{\text{astro,UCMH}} \approx r_{\text{cut}}^3 \rho_{\text{max}}^2 = K^{2/3} \left(\frac{179}{4}\right)^{4/3} r_{\text{vir}}^3 \rho_{c,m}^2. \qquad (25)$$

For a NFW profile, we can similarly evaluate (18) to calculate the WIMP-annihilation signal from such a halo

and then compare it to the UCMH result for a halo with the same mass (and hence the same $r_{vir}$),<fnref n="3" />

$$\begin{aligned} \Phi_{\text{astro,NFW}} &= r_s^3 \rho_0^2 \frac{c(3+c(3+c))}{3(1+c)^3} \\ &= 1.2 \times 10^3 r_{\text{vir}}^3 \rho_{c,m}^2 \frac{c^4(3+c(3+c))}{(1+c)^3} \\ & \quad \times \frac{1}{(\log(1+c) - c/(1+c))^2}, \end{aligned} \qquad (26)$$

where we note that the result strongly depends on the concentration parameter, $c = r_{vir}/r_s$, with an approximate $c^3$ dependence in the limit of $c \gg 1$. For this reason, the large values of $c$ generated by boosting the power spectrum

<footnote n="3">³Although the density of the NFW profile also diverges at the center, the impact of including the maximum density restriction given by (19) is negligible.</footnote>
<figure><img src="image_25.png" /></figure>

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

FIG. 19. The $\Phi_{astro}$ part of the WIMP-annihilation signal from all halos with $M_{vir} > 2 \times 10^2 M_{\odot}/h$. For the peak-to-background ratio 15 simulations, the value of $\Phi_{astro}$ is significantly boosted compared to all other halos with the same mass from the same simulation, due to the steeper profile density, as shown by Fig. 4 for the same two redshifts. However, for the peak-to-background 5 simulation, we see that the special seed no longer has a larger value than average but that its mass is considerably larger than in the peak-to-background 15 simulation, especially at $z = 30$, meaning that the special seed has undergone more merging.

do give rise to much larger WIMP-annihilation signals than would be the case with an unboosted power spectrum for halos of the same mass and also with an NFW profile. We show the mild redshift evolution of the WIMP-annihilation signal in Fig. 17 for all four levels of boost. Increasing $A_b$ by an order of magnitude has a much larger effect on $\Phi_{astro,tot}$, the total value of $\Phi_{astro}$ added up for all halos, than the redshift evolution of any given boost with redshift. Although the WIMP-annihilation signal does initially increase with redshift for the unboosted simulation, we expect this is a numerical artifact of many protohalos having not yet reached sufficient size to be detected by ROCKSTAR. Similarly, we caution that the slight decrease in $\Phi_{astro,tot}$ with redshift may be due to the decreasing minimum physical length we can resolve in an expanding background, meaning that we resolve the cores of the halos less well at late times.

We calculated the WIMP-annihilation signal by integrating (18), which assumes spherical symmetry, and compared the result to (26) using the best-fit values of $c$ from Sec. IVA and $r_{vir}$ outputs from ROCKSTAR, finding the results agree within a factor of 2. We note that Kohri et al. [9] have previously shown that the WIMP-annihilation signal is extremely sensitive to the exponent $\alpha$ when they assumed a truncated power-law profile of the form $\rho \propto (1 + r/r_{cut})^\alpha$.

However, even for the largest values of $c$ observed (using a boost factor $A_b = 1000$), the WIMP-annihilation signal as far as we are able to resolve it remains about a factor of $10^4$ smaller than it would be for a UCMH with the same mass, which means they would have to be 100 times closer in order to be equally observable. Hence, we can only observationally rule out their existence in a volume

one-millionth as large as could be probed for halos with $r^{-9/4}$ profiles. The NFW halos are, however, much more common than the assumed abundance of UCMHs given the same initial conditions, meaning that the observational constraints on the power spectrum may not weaken as strongly as may be expected, but a detailed study of this issue goes beyond the scope of this paper.

In Fig. 18, we explore whether there is any connection between the combinations of three different parameters: the virial mass of a halo $M_{vir}$, its concentration parameter $c$, and a measure for the WIMP-annihilation signal. The WIMP-signal measure is higher for heavier halos and more compact ones. It is also significantly higher in simulations with a higher boost.

In Fig. 19, we show the WIMP-annihilation signal as a function of virial mass of all halos. The most interesting feature is that the WIMP-annihilation signal of the special seed in the peak-to-background ratio 15 simulation is an order of magnitude larger than any of the other halos in the same simulation with comparable mass. If we had better resolution, we would expect that the WIMP-annihilation signal from this halo would become even larger than what is plotted. The integral of the density squared is dominated by the very center of the halo where the density is largest, and this halo has a steeper profile toward the center than all of the others; see Fig. 5.

## VI. CONCLUSIONS

We have performed the first 3D N-body simulations of UCMH formation. Starting with an isolated spherical overdensity, we have shown that this would form a steep power-law profile with the density scaling close to
$\rho \propto r^{-9/4}$ [37], and we are able to resolve this profile by up to 3 orders of magnitude of length scales at low redshift. When we include random density fluctuations with an amplitude typically 15 times smaller than the special halo, we observe that the halo descended from the special seed becomes somewhat disrupted, with the density profile becoming shallower toward the center (compare Figs. 2 and 5). This flattening occurs on such small scales that it does not make a significant difference to the numerical values of the halo properties we extract, such as the power-law steepness $\alpha$ or the WIMP-annihilation signal. Since the WIMP-annihilation signal is dominated by the square of the density in the core of the halo, our estimates for the WIMP-annihilation signal present only the lower bound on this quantity, due to the limitations of our resolution.

However, even though the special halo flattens a little, it still remains exceptional relative to the other halos that form in the simulation box. Its steepness and WIMP-annihilation signal are both much larger than any other comparable mass halo in the simulation (see Figs. 4 and 17). It is therefore likely that the inner part of this halo, unresolvable by our simulations, remains close to a power law.

When we further increase the size of background fluctuations to be five times smaller than the special halo, we see the special halo become further disrupted. In fact, in this simulation, even for the special halo, we find an NFW profile is a better fit than a power law (see Fig. 10). Moreover, although it starts as equivalent to a $5-\sigma$ fluctuation of the background, it has ceased to be exceptional even by $z = 300$. These two facts strongly suggest that, even below the scales, we can resolve the special halo will not have a steep UCMH-like profile.

Note that in all of these situations the initial special seed remains the same size; it is only the background fluctuations that change. Therefore, we can conclude that in order to form this UCMH-like power-law profile an initial fluctuation must be substantially larger than five times its background. Ricotti and Gould estimated that the critical density threshold for collapse into a UCMH was $\delta_c = 10^{-3}$ [15], and this estimate was refined to include a scale dependence by Bringmann et al. [8]. We have shown here that in realistic cases $\delta_c$ is strongly dependent on the environment. The subsequent growth, profile steepness, and WIMP-annihilation signal are all affected by the neighboring perturbations.

We note that, although almost all constraints from UCMHs have been made assuming a steep $r^{-9/4}$ profile to hold down to extremely small radii, the expectation from N-body simulations is that the emergence of NFW profiles is generic due to radial instabilities, which will always exist and grow within realistic simulations [39–41]. However, we perform the first 3D simulations with boosted initial power spectra, probing small scales and starting deep in the radiation-dominated era. We confirm that the NFW profile remains generic. In Sec. II D, we have estimated an upper

bound on the typical distance between an extremely large overdensity and the nearest overdensity with a comparable size and shown that even a rare $5-\sigma$ fluctuation will not be isolated from the evolution of neighboring halos over a long period of time.

Despite providing evidence that UCMHs of the form typically considered in the literature are not very likely to be realized in nature, we stress that our simulations with the initial power spectrum boosted on the kpc scale do form a significant number of dense NFW-like halos with masses around $10^4–10^5 M_\odot/h$. Furthermore, these halos have extremely large values of the concentration parameter $c = r_{vir}/r_s$, which can grow to over 100 by redshift 10. Assuming the central density remains constant and the concentration grows roughly as $c \propto 1/(1+z)$ [46], this would imply concentrations as large as 1000 by redshift zero. Subsequently, it would be very interesting to extend our simulations to redshift zero to study the stability of the halos we simulate with large concentration parameters. However, this would also require the inclusion of baryonic effects on these extremely small scales.

Given that the presence of halos with profiles close to $r^{-9/4}$ is unrealistic in the Universe, the observational upper bounds on the primordial power spectrum on small scales, derived from the nonobservation of UCMHs, will have to be reevaluated. However, given that we have also shown that when the primordial power spectrum is boosted even typical halos form at much higher densities and are thus much more compact than usual there are likely to be other means to constrain these small scales, which remains to be explored.

## ACKNOWLEDGMENTS

We gratefully acknowledge Camilla Penzo for providing her version of RAMSES that includes the contribution of radiation density in the Hubble function. The authors thank Jens Chluba, Peter Coles, Neal Dalal, Oliver Hahn, Massimo Ricotti, Patrick Scott, David Seery, and Teruaki Suyama for interesting discussions. J. A., C. B., and M. G. would like to thank the Department of Physics at the University of Auckland for kind hospitality. S. H. thanks Sussex University for hosting him on a return visit. C. B. is supported by a Royal Society University Research Fellowship. This work was supported by a grant from the Swiss National Supercomputing Centre (CSCS) under Project No. s710.

# APPENDIX: NOTES ABOUT THE SIMULATIONS AND HALO ANALYSIS

## 1. Initial conditions and quasilinear evolution

To set up our initial conditions, we use a customized version of the initial condition generator that is implemented in gevolution [53]. In the public version of gevolution, the initial conditions are set up in terms of a linear
displacement field and velocity potential that are obtained by multiplying a Gaussian realization of the primordial curvature perturbation by the respective linear transfer functions. The latter can be computed using a Boltzmann code such as CAMB [51] or CLASS [54]. As usual, the particle initial conditions are obtained by displacing particles from a regular lattice according to the gradient of the displacement field and assigning velocities according to the gradient of the velocity potential. In our modification of this algorithm, we introduce the possibility of creating an initial configuration corresponding to a spherical Gaussian overdensity that can act as a seed for a UCMH. Using the superposition principle that applies in the linear regime, we allow in general for any linear combination of the random realization and this spherical overdensity.

The initial conditions are set deep inside the radiation-dominated era, at redshift $z = 10000$ for the simulations with the special central halo or $z = 5000000$ for the simulations with the boosted power spectrum. We ignore the fact that the baryonic component of matter is strongly coupled at that time, assuming that at the small scales relevant for our study the perturbations in the baryon-photon plasma are irrelevant and we only have to consider CDM. Using gevolution, we follow the initial logarithmic growth of the matter perturbations up to the point where the density contrast approaches unity for the most extreme initial conditions. This happens around $z \approx 3000$.

## 2. Nonlinear evolution

The particle-mesh scheme of gevolution, which works at fixed spatial resolution, eventually becomes inadequate for tracking the detailed evolution of the very compact structures we are interested in. One common approach to deal with this issue is to use adaptive mesh refinement (AMR), that is to successively fine grain the mesh in regions of high particle density. A cell-based AMR algorithm is implemented in the public code RAMSES [55]. We use a slightly modified version in which the Hubble function takes into account radiation density and hence is more accurate at high redshift. To pass from one code to the other, we have gevolution write a snapshot that we then use as initial data for RAMSES. By running the two codes just a bit further and comparing snapshots at a later time (still before AMR is triggered in RAMSES), we convince ourselves that the change between codes did not introduce any unexpected issues.

For simplicity, we do not use the hydrodynamics modules of RAMSES, which means that baryons are effectively treated as dark matter. Baryonic effects are expected to play an important role at low redshift $z \lesssim 10$, but we are not interested in these aspects here and stop our simulations before they become a serious concern.

## 3. Spherical shells and smoothing

Once halos have formed, we identify them using the friends-of-friends halo finder algorithm ROCKSTAR [42]. For each halo, the profile of $1 + \delta$ is obtained by counting the number of particles $n$ inside the spherically symmetric shells around the center of the halo and dividing by the volume of each shell,

$$1 + \delta = \frac{\rho}{\bar{\rho}} = \frac{nV_{\text{shell}}}{N^3V}, \qquad (A1)$$

where $N^3$ is the number of all particles in the simulation and $V_{\text{shell}}$ and $V$ are the volumes of a chosen shell and the entire box, respectively. To reduce the noise at higher $r$, we gradually increase the size of the bins such that every bin is larger than the previous by a constant factor: $r_{n+1} = qr_n$, which is set to $q = 1.01$. Therefore, in the $i$th bin, the density contrast is

$$1 + \delta = \frac{3nL^3}{4\pi N^3 r_0^3 (1 - q^{i+1})^3 - (1 - q^i)^3}, \qquad (A2)$$

where $L^3$ represents the entire volume of the simulation. For the size of the starting bin $r_0$, we chose four times the size of a cell for the finest AMR level. In our simulations, the maximum refinement level is 16, and therefore the starting bin is $r_0 = 32/2^{14}$ kpc/h = 1.95 pc/h. This approach works as long as the halos are approximately spherical. To avoid empty bins in the shells in which there happen to be no particles, we apply Gaussian smoothing with a width small enough to retain the shape of a halo

<figure><img src="image_27.png" /><figcaption>FIG. 20. The profile of the Gaussian-seed halo in the peak-to-background ratio 15 simulation at $z = 30$ for two different numbers of particles. We also show the effect of smoothing—the profiles before applying smoothing are shown in dashed lines, and the ones after are continuous lines. The smoothing we applied affects the shape of the profile even less than changing the resolution. The only noticeable difference is in the first bin.</figcaption></figure>
profile. We checked that this does not degrade our
resolution.

4. Convergence tests

To understand up to what minimal radius we can trust
our results, we ran some convergence tests. In Fig. 20, we
show two simulations performed with RAMSES with the
same initial conditions but different numbers of particles:

256³ and 512³. This corresponds to a mass resolution of 0.169 $M_\odot/h$ and 0.021 $M_\odot/h$, respectively. In both cases, the maximum AMR level was 16. From the plot, we conclude that we can trust the lower-resolution simulation down to approximately $r \sim 2 \times 10^{-3}$ kpc/h. Unless indicated otherwise, our numerical results presented in Secs. III–V were obtained with simulations that had 256³ particles.

[1] P. A. R. Ade et al. (Planck Collaboration), Planck 2015 results. XIII. Cosmological parameters, Astron. Astrophys. **594**, A13 (2016).

[2] B. J. Carr, K. Kohri, Y. Sendouda, and J. Yokoyama, New cosmological constraints on primordial black holes, Phys. Rev. D **81**, 104019 (2010).

[3] A. S. Josan and A. M. Green, Constraints from primordial black hole formation at the end of inflation, Phys. Rev. D **82**, 047303 (2010).

[4] M. Kawasaki, A. Kusenko, Y. Tada, and T. T. Yanagida, Primordial black holes as dark matter in supergravity inflation models, Phys. Rev. D **94**, 083523 (2016).

[5] C. Pattison, V. Vennin, H. Assadullahi, and D. Wands, Quantum diffusion during inflation and primordial black holes, J. Cosmol. Astropart. Phys. 10 (2017) 046.

[6] B. Carr, T. Tenkanen, and V. Vaskonen, Primordial black holes from inflaton and spectator field perturbations in a matter-dominated era, Phys. Rev. D **96**, 063507 (2017).

[7] P. S. Cole and C. T. Byrnes, Extreme scenarios: the tightest possible constraints on the power spectrum due to primordial black holes, arXiv:1706.10288.

[8] T. Bringmann, P. Scott, and Y. Akrami, Improved constraints on the primordial power spectrum at small scales from ultracompact minihalos, Phys. Rev. D **85**, 125027 (2012).

[9] K. Kohri, T. Nakama, and T. Suyama, Testing scenarios of primordial black holes being the seeds of supermassive black holes by ultracompact minihalos and CMB μ-distortions, Phys. Rev. D **90**, 083514 (2014).

[10] K. J. Mack, J. P. Ostriker, and M. Ricotti, Growth of structure seeded by primordial black holes, Astrophys. J. **665**, 1277 (2007).

[11] B. C. Lacki and J. F. Beacom, Primordial black holes as dark matter: Almost all or almost nothing, Astrophys. J. **720**, L67 (2010).

[12] R. Saito and S. Shirai, Primordial black hole as a source of the boost factor, Phys. Lett. B **697**, 95 (2011).

[13] H. A. Clark, G. F. Lewis, and P. Scott, Investigating dark matter substructure with pulsar timing? I. Constraints on ultracompact minihaloes, Mon. Not. R. Astron. Soc. **456**, 1394 (2016); Erratum, Mon. Not. R. Astron. Soc. **464**, 2468 (2017).

[14] H. A. Clark, G. F. Lewis, and P. Scott, Investigating dark matter substructure with pulsar timing? II. Improved limits on small-scale cosmology, Mon. Not. R. Astron. Soc. **456**,

1402 (2016); Erratum, Mon. Not. R. Astron. Soc. **464**, 955 (2017).

[15] M. Ricotti and A. Gould, A new probe of dark matter and high-energy universe using microlensing, Astrophys. J. **707**, 979 (2009).

[16] E. Zackrisson et al., Hunting for dark halo substructure using submilliarcsecond-scale observations of macrolensed radio jets, Mon. Not. R. Astron. Soc. **431**, 2172 (2013).

[17] F. Li, A. L. Erickcek, and N. M. Law, A new probe of the small-scale primordial power spectrum: Astrometric microlensing by ultracompact minihalos, Phys. Rev. D **86**, 043519 (2012).

[18] D. Zhang, Impact of primordial ultracompact minihaloes on the intergalactic medium and first structure formation, Mon. Not. R. Astron. Soc. **418**, 1850 (2011).

[19] L. G. van den Aarssen, T. Bringmann, and Y. C. Goedecke, Thermal decoupling and the smallest subhalo mass in dark matter models with Sommerfeld-enhanced annihilation rates, Phys. Rev. D **85**, 123512 (2012).

[20] J. Chluba, A. L. Erickcek, and I. Ben-Dayan, Probing the inflaton: Small-scale power spectrum constraints from measurements of the CMB energy spectrum, Astrophys. J. **758**, 76 (2012).

[21] Y. Yang, G. Yang, X. Huang, X. Chen, T. Lu, and H. Zong, Contribution of ultracompact dark matter minihalos to the isotropic radio background, Phys. Rev. D **87**, 083519 (2013).

[22] Y. Yang, G. Yang, and H. Zong, Dark matter decay and the abundance of ultracompact minihalos, Eur. Phys. Lett. **101**, 69001 (2013).

[23] S. Shandera, A. L. Erickcek, P. Scott, and J. Y. Galarza, Number counts and non-gaussianity, Phys. Rev. D **88**, 103506 (2013).

[24] Y. Yang, G. Yang, and H. Zong, Neutrino signals from ultracompact minihalos and constraints on the primordial curvature perturbation, Phys. Rev. D **87**, 103525 (2013).

[25] V. S. Berezinsky, V. I. Dokuchaev, and Yu. N. Eroshenko, Formation and internal structure of superdense dark matter clumps and ultracompact minihaloes, J. Cosmol. Astropart. Phys. 11 (2013) 059.

[26] I. Ben-Dayan and T. Kalaydzhyan, Constraining the primordial power spectrum from SNIa lensing dispersion, Phys. Rev. D **90**, 083509 (2014).
[27] V. S. Berezinsky, V. I. Dokuchaev, and Y. N. Eroshenko, Small-scale clumps of dark matter, *Usp. Fiz. Nauk* **184**, 3 (2014) [Phys. Usp. **57**, 1 (2014)].

[28] Y. Yang, Constraints on the primordial power spectrum of small scales using the neutrino signals from the dark matter decay, *Int. J. Mod. Phys. A* **29**, 1450194 (2014).

[29] A. Natarajan, N. Zhu, and N. Yoshida, Probing the Small Scale Matter Power Spectrum through Dark Matter Annihilation in the Early Universe, arXiv:1503.03480.

[30] M. Anthonisen, R. Brandenberger, and P. Scott, Constraints on cosmic strings from ultracompact minihalos, *Phys. Rev. D* **92**, 023521 (2015).

[31] G. Aslanyan, L. C. Price, J. Adams, T. Bringmann, H. A. Clark, R. Easther, G. F. Lewis, and P. Scott, Ultracompact minihalos as probes of inflationary cosmology, *Phys. Rev. Lett.* **117**, 141102 (2016).

[32] H. A. Clark, N. Iwanus, P. J. Elahi, G. F. Lewis, and P. Scott, Heating of galactic gas by dark matter annihilation in ultracompact minihalos, *J. Cosmol. Astropart. Phys.* **05** (2017) 048.

[33] G. Beck and S. Colafrancesco, Hunting Dark Matter in ultra-compact structures within the Milky-Way, arXiv: 1612.00169.

[34] R. Emami and G. Smoot, Observational Constraints on the Primordial Curvature Power Spectrum, arXiv:1705.09924.

[35] K.-Y. Choi and T. Takahashi, New bound on low reheating temperature for dark matter in models with early matter domination, *Phys. Rev. D* **96**, 041301 (2017).

[36] A. S. Josan and A. M. Green, Gamma-rays from ultra-compact minihalos: potential constraints on the primordial curvature perturbation, *Phys. Rev. D* **82**, 083527 (2010).

[37] E. Bertschinger, Self-similar secondary infall and accretion in an Einstein-de Sitter universe, *Astrophys. J. Suppl. Ser.* **58**, 39 (1985).

[38] Y. Lithwick and N. Dalal, Self-similar solutions of triaxial dark matter halos, *Astrophys. J.* **734**, 100 (2011).

[39] R. E. Angulo, O. Hahn, A. Ludlow, and S. Bonoli, Earth-mass haloes and the emergence of NFW density profiles, *Mon. Not. R. Astron. Soc.* **471**, 4687 (2017).

[40] G. Ogiya and O. Hahn, What sets the central structure of dark matter haloes?, arXiv:1707.07693.

[41] N. Dalal, Y. Lithwick, and M. Kuhlen, The Origin of Dark Matter Halo Profiles, arXiv:1010.2539.

[42] P. S. Behroozi, R. H. Wechsler, and H.-Y. Wu, The Rockstar phase-space temporal halo finder and the velocity offsets of cluster cores, *Astrophys. J.* **762**, 109 (2013).

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

[44] J. F. Navarro, C. S. Frenk, and S. D. M. White, The structure of cold dark matter halos, *Astrophys. J.* **462**, 563 (1996).

[45] J. F. Navarro, C. S. Frenk, and S. D. M. White, A Universal density profile from hierarchical clustering, *Astrophys. J.* **490**, 493 (1997).

[46] J. S. Bullock, T. S. Kolatt, Y. Sigad, R. S. Somerville, A. V. Kravtsov, A. A. Klypin, J. R. Primack, and A. Dekel, Profiles of dark haloes. Evolution, scatter, and environment, *Mon. Not. R. Astron. Soc.* **321**, 559 (2001).

[47] S. L. Lumsden, A. F. Heavens, and J. A. Peacock, The clustering of peaks in a random Gaussian field, *Mon. Not. R. Astron. Soc.* **238**, 293 (1989).

[48] J. M. Bardeen, J. R. Bond, N. Kaiser, and A. S. Szalay, The statistics of peaks of Gaussian random fields, *Astrophys. J.* **304**, 15 (1986).

[49] D. H. Lyth, The hybrid inflation waterfall and the primordial curvature perturbation, *J. Cosmol. Astropart. Phys.* **05** (2012) 022.

[50] M. Vogelsberger, S. D. M. White, R. Mohayaee, and V. Springel, Caustics in growing cold dark matter haloes, *Mon. Not. R. Astron. Soc.* **400**, 2174 (2009).

[51] A. Lewis, A. Challinor, and A. Lasenby, Efficient computation of CMB anisotropies in closed FRW models, *Astrophys. J.* **538**, 473 (2000).

[52] T. Bringmann and C. Weniger, Gamma ray signals from dark matter: Concepts, status and prospects, *Phys. Dark Universe* **1**, 194 (2012).

[53] J. Adamek, D. Daverio, R. Durrer, and M. Kunz, General relativity and cosmic structure formation, *Nat. Phys.* **12**, 346 (2016).

[54] D. Blas, J. Lesgourgues, and T. Tram, The cosmic linear anisotropy solving system (CLASS) II: Approximation schemes, *J. Cosmol. Astropart. Phys.* **07** (2011) 034.

[55] R. Teyssier, Cosmological hydrodynamics with adaptive mesh refinement: A new high resolution code called Ramses, *Astron. Astrophys.* **385**, 337 (2002).