The clustering dynamics of primordial black boles in $N$-body simulations

We explore the possibility that Dark Matter (DM) may be explained by a non-uniform background of approximately stellar-mass clusters of Primordial Black Holes (PBHs), by simulating the evolution them from recombination to the present with over 5000 realisations using a Newtonian $N$-body code. We compute the cluster rate of evaporation, and extract the binary and merged sub-populations along with their parent and merger tree histories, lifetimes and formation rates; the dynamical and orbital parameter profiles, the degree of mass segregation and dynamical friction, and power spectrum of close encounters. Overall, we find that PBHs can constitute a viable DM candidate, and that their clustering presents a rich phenomenology throughout the history of the Universe. We show that binary systems constitute about 9.5\% of all PBHs at present, with mass ratios of $\bar{q}_{\rm B} = 0.154$, and total masses of $\bar{m}_{\rm T,\,B} = 303\,M_\odot$. Merged PBHs are rare, about 0.0023\% of all PBHs at present, with mass ratios of $\bar{q}_{\rm B}= 0.965$ with total and chirp masses of $\bar{m}_{\rm T,\,B}= 1670\,M_\odot$ and $\bar{m}_{c,{\rm M}} = 642\,M_\odot$ respectively. We find that cluster puffing up and evaporation leads to bubbles of these PBHs of order 1 kpc containing at present times about 36\% of objects and mass, with hundred pc sized cores. We also find that these PBH sub-haloes are distributed in wider PBH haloes of order hundreds of kpc, containing about 63\% of objects and mass, coinciding with the sizes of galactic halos. We find at last high rates of close encounters of massive Black Holes ($M \sim 1000\,M_\odot$), with $\Gamma^{\mathrm{S}} = (1.2^{+5.9}_{-0.9}) \times 10^{7} \mathrm{yr^{-1} Gpc^{-3}}$ and mergers with $\Gamma^{\mathrm{M}} = 1337 \pm 41 \mathrm{yr^{-1} Gpc^{-3}}$.