Dark Matter Density Profile Around a Newborn First Star

TL;DR

This work addresses how the DM density profile around the first stars shapes DM-induced GW dephasing by performing ultra-high-resolution cosmological hydrodynamics simulations of Pop III star formation in a $\sim 3\times10^{5}\,M_\odot$ minihalo. Using hierarchical zoom-in initial conditions and particle splitting, the authors achieve $m_{\rm dm,min}=3.66\times10^{-4}\,M_\odot$ and $\varepsilon_{\min}=0.022$ pc, tracking collapse from $z=49$ to $z\simeq22$ and resolving down to $\sim 0.01$ pc. They find a three-layer DM structure: a central rotational core ($R_{\rm core} \approx 3$ pc) surrounded by a cusp ($R_{\rm core} < r < R_{\rm cusp}$) and an outer halo, with the inner slope evolving to $\rho_{ m dm} \propto r^{-0.6}$ inside $\sim 1$ pc due to adiabatic contraction; however, the inner slope varies across halos due to Lyman-Werner irradiation and streaming velocities, yielding a broad distribution that can alter GW dephasing by up to an order of magnitude compared to simplistic spike models. The results stress the importance of environment and resolution-tested DM profiles in GW predictions and provide a realistic input for modeling DM effects in the first black hole binaries.

Abstract

Ambient dark matter (DM) around binary black holes can imprint characteristic signatures on gravitational waves emitted from their merger. The exact signature depends sensitively on the DM density profile around the black holes. We run very high resolution cosmological hydrodynamics simulations of first star formation that follow the collapse of a $3\times10^{5}\,M_\odot$ mini-halo from $z=49$ to $z\simeq22$. Our flagship model achieves a DM particle mass of $3.7\times10^{-4}\,M_\odot$ and resolves the inner-most structure down to $0.02\,$pc. We show that the halo experiences a two-stage gravitational collapse, where a rotating, constant-density core with $r\lesssim3\,$pc is formed first, surrounded by an extended outskirts. Baryonic infall toward the center continues to raise the local Keplerian velocity and promotes adiabatic contraction of DM. The resulting density profile has an approximately power-law shape of $ρ_{\rm dm} \propto r^{-0.6}$ inside $\sim\!1\,$pc.We find that a piecewise power-law fit reproduces the simulation result to better than 10\%, and also find numerical convergence down to $\simeq\!0.01\,$pc. The DM density profile is typical for ordinary Pop~III halos, but our additional simulations reveal that inner slope varies significantly with halo-to-halo scatter, and the effect of Lyman-Werner irradiation and of supersonic baryon-DM streaming velocities, implying a wide distribution of slopes rather than a single universal curve. The large variation should be considered when calculating the predicted DM-induced dephasing of gravitational waves by up to an order of magnitude relative to the classical analytic model of the DM spike.

Figures (9)

• Figure 1: Averaged gas temperature of the collapsing cloud as a function of the gas particle number density for Model L3. The horizontal dotted line shows the CMB temperature floor $T_{\rm CMB}(z)=2.73(1+z)\,$K at $z=20.7$. Note the effect of HD cooling at $\rho_{\rm b} > 10^5$. The temperature does not increase gradually unlike in the usual case with H$_2$ cooling only.
• Figure 2: Density distributions of DM ($\rho_{\rm dm}$; top panels) and baryonic ($\rho_{\rm b}$; bottom panels) components projected on a cross-section through the density center. The left, center, and right panels show plots over $1000$, $10$, and $1$ pc per side, respectively.
• Figure 3: Radial density profiles of the DM and baryonic components for Model L3. Panel (a): the thick lines are simulation results when the maximum baryon density is $\rho_{\rm b,max}=10^9\,{\rm GeV\,cm^{-3}}$ and the thin lines are fitting functions for the two components (eqs. \ref{['eq:rhodm_fit']} and \ref{['eq:rhob_fit']}). The vertical dotted line is $R_{\rm III}=10\,R_\odot$, the radius of the typical Pop III star with mass $M_{\rm III} = 100\,M_\odot$HosokawaOmukai2009. Panel (b): enlarged view of the range indicated by the dotted rectangle in panel (a). The dashed lines are simulation results when baryon density exceeds DM density, $\rho_{\rm b,max}=10^3\,{\rm GeV\,cm^{-3}}$. Arrows indicate characteristic radii of the DM density profile: virial radius of the DM halo ($R_{\rm virial}$), scale radius where DM density profile steepens ($R_{\rm cusp}$), core radius at which the DM density flattens before adiabatic contraction due to gas condensation ($R_{\rm core}$), and gravitational softening length of the finest DM particle ($\varepsilon_{\rm min}$). Panel (c): relative deviation of our fitting functions from the simulation results. The vertical dotted lines are the radii at which the fitting function of the DM component switches, $r=1.0, 9.4, 66$ pc, respectively.
• Figure 4: Radial profiles of the enclosed masses of baryon and DM components (red and blue lines) for Model L3 when $\rho_{\rm b,max}=10^3$ and $10^9\,{\rm GeV\,cm^{-3}}$ (dashed and solid line styles). The long-dashed line shows a correlation as $M_{\rm enc} \propto r^{1.5}$.
• Figure 5: Radial profiles of the rotation velocity with respect to the Keplerian velocity for Model L3 when $\rho_{\rm b}=10^3$ and $10^9\,{\rm GeV\,cm^{-3}}$. The red lines show the ratio to the Keplerian velocity obtained from the total mass, $v_{\rm Kepler,dm+b} = \sqrt{G(M_{\rm enc,dm}(r)+M_{\rm enc,b}(r))/r}$, whereas the blue lines show the ratio to the Keplerian velocity obtained from the DM mass only, $v_{\rm Kepler,dm} = \sqrt{GM_{\rm enc,dm}(r)/r}$.