# Superradiant Production of
Heavy Dark Matter from
Primordial Black Holes

Nicolás Bernal  Yuber F. Perez-Gonzalez  and Yong Xu

###### Abstract

Rotating black holes (BHs) can efficiently transfer energy to the surrounding
environment via superradiance. In particular, when the Compton length of a
particle is comparable to the gravitational radius of a BH, the particle’s
occupation number can be exponentially amplified. In this work, we investigate
the effect of the primordial-black-hole (PBH) superradiant instabilities on
the generation of heavy bosonic dark matter (DM) with mass above $\sim$ 1 TeV.
Additionally, we analyze its interplay with other purely gravitational and
therefore unavoidable DM production mechanisms such as Hawking emission and
the ultraviolet freeze-in. We find that superradiance can significantly
increase the DM density produced by PBHs with respect to the case that only
considers Hawking emission, and hence lower initial PBH densities are
required.

PI/UAN-2022-717FT

IPPP/22/33

##  1 Introduction

Superradiance is a radiation enhancement phenomena, ubiquitous in many
physical systems. In general relativity, a spinning black hole (BH) can
feature superradiant instabilities, allowing efficient transfer of energy and
angular momentum to the surrounding environment (for a review, see e.g. Ref.
[1]). In particular, the frequency of a wave that scatters with a Kerr BH can
be enhanced if it is smaller than the BH angular velocity [2], similar to the
Penrose process in the ergoregion of a Kerr BH [3]. Additionally, Press and
Teukolsky proposed the BH bomb Gedankenexperiment [4]: if one surrounds the BH
with a reflecting mirror, the frequency can be further amplified when the wave
bounces back and forth. In reality, the gravitational interaction with the
massive field that scatters can naturally act as a mirror [5].

For boson fields, the occupation number can be exponentially enhanced because
of superradiant instabilities, leading to the formation of gravitational
clouds around BHs, which could leave observable astrophysical signals. For
example, these clouds could source gravitational waves (GWs) via $i)$ direct
emission due to a dipole or quadrupolar configuration, $ii)$ through the
transition between different harmonic quantum levels, and $iii)$ from boson
condensate annihilation; for more details, see e.g. Refs. [1, 6]. This gives
rise to a potential new avenue to constrain boson fields captured by
superradiance.

There are two key parameters in superradiant instabilities: $i)$ the particle
mass $m_{\text{dm}}$ (corresponding to a Compton length
$\lambda_{C}=1/m_{\text{dm}}$), and $ii)$ the BH mass $M_{\text{BH}}$
(associated to a gravitational radius $r_{g}=G\,M_{\text{BH}}$, $G$ being the
gravitational constant). The detailed dynamics of the gravitational clouds
invokes a solution of the Klein-Gordon equation in a Kerr metric. Approximate
solutions are reported in the literature in the limits $\alpha\gg 1$ [7],
$\alpha\ll 1$ [8] and $\alpha\simeq 1$ [9], where $\alpha\equiv
r_{g}/\lambda_{C}$ is the dimensionless “fine-structure constant”.
Superradiance effects are maximized when $\alpha\simeq 1$ [9], corresponding
to a boson mass
$m_{\text{dm}}\simeq\mathcal{O}(10^{-11})\times\left(M_{\odot}/M_{\text{BH}}\right)$
eV [6]. For astrophysical solar-mass BHs, it is then possible to probe
ultralight bosons: see, e.g. the phenomenology of superradiance with QCD
axions [10, 11], axion-like particles (ALPs) [12, 13], spin-1 vector fields
[14, 15, 16, 12], spin-2 fields [17, 18], and self-interacting fields [19].
Moreover, the effect of superradiant instabilities with even lighter bosons
has recently been explored in the context of supermassive BHs [20, 21, 22].

Complementary to previous studies that consider astrophysical BHs, in this
work we instead investigate the phenomenology of superradiant instabilities in
the context of primordial BHs (PBH), which could have been generated in the
early Universe due to large density perturbations [23, 24]. PBHs with a large
angular momentum could also have been formed, for example, $i)$ from
collapsing domain walls [25], $ii)$ during a matter-dominated era [26, 27], or
$iii)$ due to the assembly of matter-like objects such as Q-balls or oscillons
[28]. Alternatively, PBHs originally produced with a small angular momentum
can gain spin111Note that here we mean the dimensionless spin parameter
$a_{\star}$. by emitting light scalar particles [29, 30, 31]. Here, we focus
on PBHs with mass $M_{\text{BH}}\lesssim\mathcal{O}(10^{9})~{}\text{g}$, which
have fully evaporated before the onset of Big Bang nucleosynthesis (BBN), so
that they would not jeopardize its successful predictions. During the
evaporation process, any particle including dark matter (DM) with mass small
than the horizon temperature can be sourced [32, 33]. If superradiant
instabilities exist for supermassive astrophysical BHs, so do those BHs with
smaller mass and primordial origin [7, 8]. Such instabilities with PBHs could
naturally give rise to interesting cosmological signatures. For example, it
has recently been shown that pion clouds can be formed around PBHs with mass
around $\sim 10^{12}$ kg due to the superradiant instability [34], offering a
potential new avenue to constrain the abundance of PBHs by considering
subsequent decay or annihilation of the pions.

In addition to studying the role of PBH superradiant instabilities in DM
production, we pay particular attention to its interplay with other
unavoidable gravitational channels such as UV freeze-in [35, 36, 37], and
Hawking evaporation [38, 39, 40, 41, 42, 43, 36, 44, 45, 46, 47, 48, 49]. For
illustration purposes, we first consider a spin-0 scalar field as the DM
candidate; the extensions to the spin-1 and spin-2 cases will be investigated
elsewhere. We show that PBH superradiance gives rise to an unexplored
nonthermal channel for generating heavy DM with mass above $\sim$ 1 TeV, and
could source the entire observed abundance of DM relics. We carefully track
the evolution of DM by solving a system of coupled equations for the evolution
of its number density, together with the PBH mass and spin. We find that
sizable reductions on the initial PBH abundance are required in the regimes
where superradiance is efficient.

The rest of the paper is organized as follows. We first briefly revisit the
formation of and constraints on PBHs in Sect. 2. Then in Sect. 3 we set up the
Hawking evaporation formalism and then discuss superradiant instabilities with
a particular focus on the evolution of the number of particles gravitationally
bounded to PBHs. In Sect. 4 we review the gravitational UV freeze-in. We
collect the complete numerical results with all the previously-discussed
gravitational contributions in Sect. 5. Finally, we summarize our findings in
Sect. 6. Along the paper we denote the reduced Planck mass as $M_{P}\equiv
1/(\sqrt{8\pi\,G})$, and we use natural units where $\hbar=c=k_{\rm B}=1$.

##  2 Primordial Black Holes

Different mechanisms could lead to the formation of PBHs in the early Universe
[32, 50]. For example, if above a threshold, density fluctuations can collapse
into a BH when reentering the Hubble horizon [51]. Other scenarios for PBH
formation include the collapse of cosmic string loops, the collapse of domain
walls, or the collision of bubbles [52]. Here, we remain agnostic about the
formation mechanism, and assume that an initial PBH density was formed in a
radiation-dominated era. For a plasma temperature $T=T_{\text{in}}$, it is
assumed that the initial mass of PBH $M_{\text{in}}$ is proportional to the
horizon mass [23, 24],

$M_{\text{in}}\equiv
M_{\text{BH}}(T_{\text{in}})=\frac{4\pi}{3}\,\kappa\,\frac{\rho_{R}(T_{\text{in}})}{H^{3}(T_{\text{in}})}\,,$
(2.1)

with numerical efficiency factor $\kappa\simeq 0.2$. Here $\rho_{R}$ and $H$
denote the radiation energy density and the Hubble rate, respectively.
Moreover, we adopt a monochromatic mass spectrum such that all PBHs were
generated with same mass.222Extended PBH mass functions arise naturally if
PBHs are created from inflationary fluctuations or cosmological phase
transitions (see, e.g., Refs. [53, 54, 55, 56, 57, 58, 59, 60, 61, 62]), and a
log-normal distribution represents a good approximation for models with a
symmetric peak in the power spectrum of curvature perturbations [63, 64, 65].
The extension to more realistic nonmonochromatic mass distributions will be
investigated elsewhere. We characterize the PBH initial energy density
$\rho_{\text{BH}}(T_{\text{in}})$ through the $\beta$ parameter, defined by

$\beta\equiv\frac{\rho_{\text{BH}}(T_{\text{in}})}{\rho_{R}(T_{\text{in}})}\,.$
(2.2)

Since PBHs behave as nonrelativitic matter, $\rho_{\text{BH}}$ redshifts
slower than radiation, and a nonstandard cosmological expansion era could
occur [66]. A phase dominated by PBH with $\rho_{\text{BH}}>\rho_{R}$ can be
triggered if the initial fraction is above a critical value $\beta_{c}$,

$\beta_{c}\equiv\frac{T_{\text{ev}}}{T_{\text{in}}}\,,$ (2.3)

where $T_{\text{ev}}$ is given by [67, 68]

$T_{\text{ev}}\simeq\left(\frac{g_{\star}(T_{\text{in}})}{640}\right)^{1/4}\left(\frac{M_{P}^{5}}{M_{\text{in}}^{3}}\right)^{1/2},$
(2.4)

and corresponds to the SM temperature at which the PBHs have completely
evaporated. Here $g_{\star}(T)$ represents the number of relativistic degrees
of freedom contributing to the SM radiation energy density. Note that if
$\beta>\beta_{c}$, PBHs would start to dominate the total energy density at a
temperature $T=T_{\text{eq}}$, when the equality between SM radiation and the
PBH energy densities $\rho_{R}(T_{\text{eq}})=\rho_{\text{BH}}(T_{\text{eq}})$
occurred. Such a radiation-PBH equality temperature is given by

$T_{\text{eq}}=\beta\,T_{\text{in}}\left(\frac{g_{\star
s}(T_{\text{in}})}{g_{\star s}(T_{\text{eq}})}\right)^{1/3},$ (2.5)

where $g_{\star s}(T)$ is the number of relativistic degrees of freedom that
contribute to the SM entropy density.

In the evaporation process, all particles with a mass smaller than the horizon
temperature can be emitted. If such evaporation occurred during or after BBN,
the large particle injection could spoil the standard prediction of the
abundance of light elements. To avoid this, we consider PBHs that have fully
evaporated before the onset of BBN, that is,
$T_{\text{ev}}>T_{\text{BBN}}\simeq 4$ MeV [69, 70, 71, 72, 73, 74], which can
be translated into an upper bound on the initial mass of PBHs

$M_{\text{in}}\lesssim\mathcal{O}\left(10^{9}\right)\text{g}\,,$ (2.6)

as implied by Eq. (2.4). On the other hand, one can obtain a lower limit on
$M_{\text{in}}$ taking into account the upper bound on inflationary scale
$H_{I}$ from constraints on the inflationary tensor-to-scalar ratio. According
to the latest Planck collaboration results, one has $H_{I}\leq 2.5\times
10^{-5}~{}M_{P}$ [75], implying

$M_{\text{in}}\gtrsim 4\pi\,\kappa\,\frac{M_{P}^{2}}{H_{I}}\simeq
0.1~{}\text{g}\,.$ (2.7)

Finally, the GWs produced due to evaporation of PBH could also spoil BBN, as
they will contribute to the effective number of relativistic species $\Delta
N_{\rm eff}$ [76, 77]. This can be avoided if [77]

$\beta\lesssim 1.1\times
10^{-6}\left(\frac{\kappa}{0.2}\right)^{-\frac{1}{2}}\left(\frac{M_{\text{in}}}{10^{4}~{}\text{g}}\right)^{-\frac{17}{24}}.$
(2.8)

After describing the generic limits from evaporation, let us describe in
detail the properties of spinning PBHs. Since superradiant amplification takes
place only for such BHs, we will focus on Kerr PBHs in what follows.

##  3 Superradiance and Hawking Radiation

In this section, we first have a quick overview of Hawking radiation, then
focus on the particle number growth by BHs due superradiant instabilities.
After that we put things together, and describe the interplay between
superradiance and Hawking evaporation.

###  3.1 Hawking Radiation from Kerr Black Holes

PBHs could have been produced with a large initial angular momentum [26, 27,
25, 28, 29, 30, 31]. For simplicity, we consider a monochromatic spin
distribution, where all the PBHs had the same initial spin, without specifying
the initial mechanism that generated it. A Kerr BH is characterized by its
instantaneous mass $M_{\text{BH}}$ and dimensionless spin parameter
$a_{\star}\equiv 8\pi\,J\,M_{P}^{2}/M_{\text{BH}}^{2}\in[0,1]$, with $J$
denoting the total angular momentum of the BH.

The line element of the Kerr metric in Boyer-Lindquist coordinates is given by
[8]

$ds^{2}=-dt^{2}+\frac{2\,r_{g}\,r}{\Sigma}\left(dt-
a_{\star}\,r_{g}\,\sin^{2}\theta\,d\varphi\right)^{2}+\frac{\Sigma}{\Delta}\,dr^{2}+\Sigma\,d\theta^{2}+\left(r^{2}+a_{\star}^{2}\,r_{g}^{2}\right)\sin\theta^{2}\,d\varphi^{2},$
(3.1)

where $t$ denotes cosmic time and $(r,\,\theta,\,\varphi)$ are spherical
coordinates, $\Sigma\equiv r^{2}+a_{\star}^{2}\,r_{g}^{2}\,\cos^{2}\theta$,
and $\Delta\equiv r^{2}-2\,r_{g}\,r+a_{\star}^{2}\,r_{g}^{2}$, with
$r_{g}\equiv G\,M_{\text{BH}}$ denoting the gravitational radius. The event
horizon $r_{+}$ and the Cauchy horizon $r_{-}$ correspond to the coordinate
singularities $\Delta=0$, and are given by

$r_{\pm}=r_{g}\left(1\pm\sqrt{1-a_{\star}^{2}}\right).$ (3.2)

The region $r_{+}\leq r\leq r_{t}$, with

$r_{t}\equiv r_{g}\,\left(1+\sqrt{1-a_{\star}^{2}\cos^{2}\theta}\right),$
(3.3)

denotes the so-called ergosphere, in which particles can gain energy from
spinning BHs, and corresponds to the surface where the argument of $dt$
vanishes. Meanwhile, the horizon temperature for spinning BHs is given by

$T_{\text{BH}}=\frac{2M_{P}^{2}}{M_{\text{BH}}}\,\frac{\sqrt{1-a_{\star}^{2}}}{1+\sqrt{1-a_{\star}^{2}}}\,.$
(3.4)

As mentioned above, during the Hawking evaporation process, all particles with
masses smaller than $T_{\text{BH}}$ can be generated. Within a time interval
$[t,\,t+dt]$ and a momentum interval $[p,\,p+dp]$, the production rate for
particle species $i$ with spin $s_{i}$, mass $\mu_{i}$ and internal degrees of
freedom $g_{i}$ is [32, 47]

$\frac{d^{2}\mathcal{N}_{i}}{dp\,dt}=\frac{g_{i}}{2\pi^{2}}\sum_{l=s_{i}}\sum_{m=-l}^{l}\frac{\sigma_{s_{i}}^{lm}\left(M_{\text{BH}},p,a_{\star}\right)}{\exp\left[\left(E_{i}-m\,\Omega\right)/T_{\text{BH}}\right]-(-1)^{2s_{i}}}\,\frac{p^{3}}{E_{i}}\,,$
(3.5)

where $E_{i}(p)=\sqrt{p^{2}+\mu_{i}^{2}}$ denotes the energy, $\sigma_{s_{i}}$
describes the absorption cross-section, and $l$ and $m$ correspond to the
total and magnetic quantum numbers, respectively. Finally, in Eq. (3.5),
$\Omega$ is the BH horizon angular velocity, a crucial quantity for
superradiance, which is given by

$\Omega\equiv\frac{4\pi\,a_{\star}}{1+\sqrt{1-a_{\star}^{2}}}\,\frac{M_{P}^{2}}{M_{\text{BH}}}\,.$
(3.6)

Summing over all possible particle species and integrating over the phase
space the rate in Eq. (3.5), one has the following set of equations describing
the time evolution of the BH mass and spin, without including the effect of
superradiance [78]

$\displaystyle\frac{dM_{\text{BH}}}{dt}$
$\displaystyle=-\varepsilon\left(M_{\text{BH}},a_{\star}\right)\frac{M_{P}^{4}}{M_{\text{BH}}^{2}}\,,$
(3.7a) $\displaystyle\frac{da_{\star}}{dt}$
$\displaystyle=-a_{\star}\left[\gamma\left(M_{\text{BH}},a_{\star}\right)-2\varepsilon\left(M_{\text{BH}},a_{\star}\right)\right]\frac{M_{P}^{4}}{M_{\text{BH}}^{3}}\,,$
(3.7b)

where the evaporation functions
$\varepsilon\left(M_{\text{BH}},a_{\star}\right)\equiv\sum_{i}g_{i}\,\varepsilon_{i}(z_{i})$
and $\gamma\left(M_{\text{BH}}\right)\equiv\sum_{i}g_{i}\,\gamma_{i}(z_{i})$
are given by the dimensionless quantities

$\displaystyle\varepsilon_{i}\left(z_{i},a_{\star}\right)$
$\displaystyle=\frac{27}{2\pi^{2}}\int_{z_{i}}^{\infty}\sum_{l,m}\frac{\psi_{s_{i}}^{lm}\left(x,a_{\star}\right)\left(x^{2}-z_{i}^{2}\right)}{\exp\left(x^{\prime}/2f\left(a_{\star}\right)\right)-(-1)^{2s_{i}}}\,xdx\,,$
(3.8a) $\displaystyle\gamma_{i}\left(z_{i},a_{\star}\right)$
$\displaystyle=108\int_{z_{i}}^{\infty}\sum_{l,m}\frac{m\,\psi_{s_{i}}^{lm}\left(x,a_{\star}\right)\left(x^{2}-z_{i}^{2}\right)}{\exp\left(x^{\prime}/2f\left(a_{\star}\right)\right)-(-1)^{2s_{i}}}\,dx\,,$
(3.8b)

with $z_{i}=\mu_{i}/T_{\text{BH}}$, $\psi_{s_{i}}(E,\mu)\equiv 64\pi
M_{P}^{2}\sigma_{s_{i}}(E,\mu)/(27M_{\text{BH}}^{2})$,
$x^{\prime}=x-m\,\Omega^{\prime}$, $x=M_{\text{BH}}\,E_{i}/M_{P}^{2}$,
$\Omega^{\prime}=M_{\text{BH}}\,\Omega/M_{P}^{2}$, and
$f(a_{\star})=\sqrt{1-a_{\star}}/(1+\sqrt{1-a_{\star}})$. We refer the reader
to Refs. [78, 48, 47] for further details. We next describe the superradiant
instability for Kerr BHs.

###  3.2 Superradiance in a Nutshell

When the Compton length $\lambda_{C}=1/m_{\text{dm}}$ of a boson with mass
$m_{\text{dm}}$ is comparable to the BH gravitational radius, a superradiant
instability can occur, leading to an enhancement of the particle occupation
number and a formation of bound states. Such bound states are generally
referred to as _gravitational atoms_ , characterized by a “fine-structure
constant” [11]

$\alpha\equiv\frac{r_{g}}{\lambda_{C}}=G\,M_{\text{BH}}\,m_{\text{dm}}\simeq
0.38\left(\frac{M_{\text{BH}}}{10^{7}~{}{\rm
g}}\right)\left(\frac{m_{\text{dm}}}{10^{7}~{}{\rm GeV}}\right).$ (3.9)

Similarly to the hydrogen atom, the energy spectrum is controlled by three
quantum numbers, namely the principal, orbital, and magnetic quantum numbers
$(n,\,l,\,m)$, and at leading order is given by [11, 34]

$\omega\simeq m_{\text{dm}}\left(1-\frac{\alpha^{2}}{2n^{2}}\right).$ (3.10)

Generically, in order to allow superradiance to occur, the energy or frequency
$\omega$ has to be smaller than the BH angular velocity. In other words, if
the latter is smaller than the former, PBHs would be unable to further amplify
$\omega$. Let us now state the superradiance condition in a more mathematical
way. At leading order $\omega\sim m_{\text{dm}}$ (cf. Eq. (3.10)), and
$\Omega\sim(2\,r_{g})^{-1}$ for $a_{\star}\simeq 1$ (cf. Eq. (3.6)). To allow
superradiance to occur, $\omega<m\,\Omega$, which translates into

$\alpha<\frac{1}{2}\,,$ (3.11)

for $m=l=1$. As mentioned above, due to superradiance, the frequencies can be
enhanced, which could also be understood as a growth of the particle
occupation number. The growth rate is defined to be the imaginary part of the
bound state frequency $\mathfrak{Im}(\omega)$.333Further details can be found
in Appendix A. In terms of the horizon radius $r_{+}$, for the dominant
unstable mode ($n=2$ and $l=m=1$) the growth rate $\Gamma_{\text{sr}}$ can be
approximated as [8, 1, 15, 34]

$\Gamma_{\text{sr}}=\frac{m_{\text{dm}}}{24}\left(\frac{M_{\text{BH}}\,m_{\text{dm}}}{8\pi\,M_{P}^{2}}\right)^{8}\left(a_{\star}-2\,m_{\text{dm}}\,r_{+}\right),$
(3.12)

where $r_{+}$ represents the event horizon as defined in Eq. (3.2). Note that
in the limit $a_{\star}\to 0$, $r_{+}\to 0$, and therefore the growth rate
$\Gamma_{\text{sr}}$ vanishes, as expected, since superradiant instabilities
do not occur for Schwarzschild BHs. This approximation models well the
behavior of $\Gamma_{\text{sr}}$ for small values of $\alpha$, however, fails
in the large $\alpha$ regime. A better determination of $\Gamma_{\text{sr}}$
requires a numerical solution.

For completeness, we next describe the numerical method that we have used to
compute the growth rates, following Ref. [9]. To determine the frequencies of
the bound states, it is required to solve the Klein-Gordon equation on a Kerr
background [9, 12]. After a standard variable separation, two power series
solutions are proposed to solve the radial and angular equations, which lead
to three-term recurrence relations between the different expansion
coefficients [79, 9]. Implicit conditions for continued fractions, which are
only satisfied by the bound-state frequencies, are obtained after performing
algebraic manipulations. Thus, one needs to employ numerical method to solve
such implicit conditions and obtain the growth rate. We use a scipy standard
root finder that applies the HYBR method [80] to simultaneously solve the
radial and angular implicit conditions. We have checked that our code gives
the same rates as those presented in Ref. [9], which we show in Fig. 1 for
$l=m=1$ (black) and $l=m=2$ (blue), for $a_{\star}=0.99$, 0.9, 0.8 and 0.6. In
what follows, we only consider the instability on the state $l=m=1$ since it
grows exponentially faster than the others. Contributions to the DM abundance
coming from other $l$ and $m$ states will be discussed elsewhere. Detailed
numerical simulations show that for $a_{\star}\simeq 1$, the superradiant
instability is largest when $\alpha\sim 0.42$ [9, 81], as shown in Fig. 1.

Figure 1: Normalized growth rate $\Gamma_{\text{sr}}$ as function of the fine-
structure constant $\alpha$, for $l=m=1$ (black) or $l=m=2$ (blue), and
different values of the initial BH spin.

Superradiant instabilities exponentially amplify the number of DM particles.
Defining $N_{\text{dm}}^{\text{sr}}$ to be the total number of DM particles
gravitationally bounded to a PBH and considering for the moment only
superradiance, the dynamics of the system is governed by the system of coupled
equations [11, 1, 82]

$\displaystyle\frac{dN_{\text{dm}}^{\text{sr}}}{dt}$
$\displaystyle=\Gamma_{\text{sr}}\,N_{\text{dm}}^{\text{sr}}\,,$ (3.13a)
$\displaystyle\frac{dM_{\text{BH}}}{dt}$
$\displaystyle=-m_{\text{dm}}\,\frac{dN_{\text{dm}}^{\text{sr}}}{dt}\,,$
(3.13b) $\displaystyle\frac{dJ}{dt}$
$\displaystyle=-\sqrt{2}\,\frac{dN_{\text{dm}}^{\text{sr}}}{dt}\,.$ (3.13c)

Apart from having $\alpha\sim 1$, superradiant exponential growth requires a
nonvanishing initial population of gravitationally bounded particles. However,
it is interesting to realize that the exact number of initial particles has a
small impact on its final number. Noticing that

$\frac{dN_{\text{dm}}^{\text{sr}}}{da_{\star}}=\frac{dN_{\text{dm}}^{\text{sr}}}{dt}\times\frac{dt}{da_{\star}}=\left[2\,\frac{m_{\text{dm}}}{M_{\text{BH}}}\,a_{\star}-8\pi\sqrt{2}\left(\frac{M_{P}}{M_{\text{BH}}}\right)^{2}\right]^{-1}\simeq-\frac{1}{8\pi\sqrt{2}}\left(\frac{M_{\text{BH}}}{M_{P}}\right)^{2},$
(3.14)

it can be seen that the final amount of DM particles
$N_{\text{dm}}^{\text{sr}}(\text{end})$ after the superradiant instability is

$N_{\text{dm}}^{\text{sr}}(\text{end})\simeq
N_{\text{dm}}^{\text{sr}}(\text{in})+1.5\times
10^{19}\left(\frac{M_{\text{BH}}}{10^{5}~{}\text{g}}\right)^{2}\left(a^{\text{ini}}_{\star}-a^{\text{end}}_{\star}\right)\simeq\mathcal{O}\left(10^{18}\right)\left(\frac{M_{\text{BH}}}{10^{5}~{}\text{g}}\right)^{2},$
(3.15)

featuring a small additive dependence on the initial number of DM particles
$N_{\text{dm}}^{\text{sr}}(\text{in})$.

###  3.3 Interplay of Superradiance and Hawking Radiation

As mentioned in Sect. 3.1, on top of the contribution from superradiance, DM
particles are also unavoidably sourced by Hawking evaporation. In this section
we investigate the interplay between the two mechanisms. The evolution of the
system is governed by the set of coupled equations:

$\displaystyle\frac{dN_{\text{dm}}^{\text{sr}}}{dt}$
$\displaystyle=\Gamma_{\text{sr}}\,N_{\text{dm}}^{\text{sr}}\,,$ (3.16a)
$\displaystyle\frac{dM_{\text{BH}}}{dt}$
$\displaystyle=-\varepsilon\,\frac{M_{P}^{4}}{M_{\text{BH}}^{2}}-m_{\text{dm}}\,\Gamma_{\text{sr}}\,N_{\text{dm}}^{\text{sr}}\,,$
(3.16b) $\displaystyle\frac{da_{\star}}{dt}$
$\displaystyle=-a_{\star}\left[\gamma-2\varepsilon\right]\frac{M_{P}^{4}}{M_{\text{BH}}^{3}}-8\pi\left[\sqrt{2}-2\alpha\,a_{\star}\right]\Gamma_{\text{sr}}\,N_{\text{dm}}^{\text{sr}}\,\frac{M_{P}^{2}}{M_{\text{BH}}^{2}}\,,$
(3.16c) $\displaystyle\frac{dn_{\text{dm}}}{dt}+3H\,n_{\text{dm}}$
$\displaystyle=n_{\text{BH}}\left[\Gamma_{\text{BH$\to$DM}}+\Gamma_{\text{sr}}\,N_{\text{dm}}^{\text{sr}}\right],$
(3.16d)

with $n_{\text{dm}}$ and $n_{\text{BH}}$ denoting the DM and PBH number
density, respectively. The first line describes the usual evolution of total
bounded DM numbers due to superradiant instability. The second and third lines
are a combination of Eqs. (3.7) and (3.13), and describe the evolution of PBH
mass and spin. Finally, the last line denotes the modified Boltzmann equation
for the DM number density with contributions from both the superradiant
instability and the Hawking evaporation. In the first term of the right-hand
side (RHS), $\Gamma_{\text{BH$\to$DM}}$ denotes the DM production rate from BH
through evaporation, which can be obtained by integrating Eq. (3.5) over the
phase space for a scalar DM field [47]. In the second term,
$\Gamma_{\text{sr}}\,N_{\text{dm}}^{\text{sr}}$ corresponds to the DM
production rate per BH. Note that we are assuming that after the evaporation
of the PBHs, DM behaves as free particles. This allows us to explore the
maximum reach of superradiant DM production. The details the fate of the DM
clouds will be explored in a future work. Finally, Eqs. (3.16) reduce to Eqs.
(3.7) or (3.13) in the limit where superradiance or Hawking emission are
artificially switched off.

To allow superradiance to take place, the PBH should not evaporate before the
instability becomes significant. Thus, we can compare the typical timescale
for superradiance to the lifetime of PBHs. The PBH lifetime is roughly [44]

$\tau=\frac{160}{\pi g_{\star}}\frac{M_{\text{in}}^{3}}{M_{P}^{4}}\,,$ (3.17)

whereas the typical timescale for superradiance $t_{\text{sr}}$ corresponds to
the inverse of the superradiance rate in Eq. (3.12),

$t_{\text{sr}}\sim\frac{1}{\Gamma_{\text{sr}}}\simeq\frac{24}{m_{\text{dm}}}\left(\frac{8\pi\,M_{P}^{2}}{M_{\text{in}}\,m_{\text{dm}}}\right)^{8}\frac{1}{a_{\star}}\,.$
(3.18)

To develop superradiance, $t_{\text{sr}}<\tau$ is required, so we can estimate
a lower bound on the initial PBH mass, depending on the DM mass and initial
PBH spin

$M_{\text{in}}\gtrsim
70\left(\frac{1}{a_{\star}}\right)^{1/11}\left(\frac{10^{11}~{}\text{GeV}}{m_{\text{dm}}}\right)^{9/11}\text{g}\,.$
(3.19)

However, we stress that Kerr BHs lose angular momentum much faster than their
mass [78], so they effectively become Schwarzschild BHs before completely
evaporating. As the BH should have enough angular momentum for the
superradiant instability to rise, the estimate in Eq. (3.19) is conservative,
since it does not consider the time when the BH would have shed their spin.
Still, we have verified that this characteristic of Kerr BHs only affects the
estimated lower bound for light PBHs $M_{\text{in}}\lesssim 10$ g.

Figure 2: Evolution of the PBH mass (left), its spin (center), and the DM
comoving number density (right) for $M_{\text{in}}=3.2\times 10^{4}$ g,
$a_{\star}=0.999$, $\beta=4.4\times 10^{-22}$, and $m_{\text{dm}}=10^{9}$ GeV.
The dotted lines only take into account Hawking emission, whereas the solid
lines also account for superradiance. The latter case reproduces the observed
DM abundance. Note that SR represents superradiance.

The coupled system (3.16), together with the Friedmann equations for the PBH
and SM radiation energy densities are numerically solved using the Python
package ULYSSES [83].444https://github.com/earlyuniverse/ulysses An example of
the evolution of the PBH mass and spin, and the DM number density is shown in
Fig. 2, for $M_{\text{in}}=3.2\times 10^{4}$ g, $a_{\star}=0.999$,
$\beta=4.4\times 10^{-22}$ and $m_{\text{dm}}=10^{9}$ GeV, a benchmark that
reproduces the observed relic abundance once superradiance is included. The
dotted lines only take into account Hawking emission, whereas the solid lines
also contain the contribution from superradiance. The effect of the
superradiant instability manifests itself significantly when $t\sim
1/\Gamma_{\text{sr}}$, which corresponds to
$a/a_{\text{ini}}\sim\mathcal{O}\left(10^{5}\right)$ as depicted in Fig. 2,
with $a_{\text{ini}}$ being the scale factor at $T=T_{\text{in}}$. For
simplicity, we define $a_{\text{ini}}=1$ in the following. At the end of the
superradiant era, the BH mass is reduced only by a factor of $\sim 5\%$, while
the spin decreases to $\sim 40\%$ of its initial value. The right panel of
Fig. 2 tracks the comoving DM number density, namely
$N_{\text{dm}}=a^{3}n_{\text{dm}}$. The superradiant instability exponentially
amplifies the DM number density by $\sim\mathcal{O}(11)$ orders of magnitude
in a short period of time. Finally, once superradiance ends, the PBH evolution
is governed by the usual Hawking evaporation. Interestingly, compared to the
case where superradiance is absent, the final DM number density gets boosted
by factor of $\sim 21$. Hawking radiation produces an additional contribution
to the total DM number density, but its impact is negligible compared to the
abundance produced by superradiance.

Figure 3: The solid and dotted black lines show the regions of the parameter
space where ratio of the DM relic abundance with and without the inclusion of
superradiance is 15 and 5, respectively. The left and right panels correspond
to $a_{\star}=0.999$ and $a_{\star}=0.6$. In the gray area the universe could
be overclosed due to the gravitational freeze-in, unless there is a large
entropy dilution. The red regimes are excluded by different observables
described in the text.

The relative increase of the DM relic abundance due to superradiance can also
be seen in Fig. 3. The solid and dotted black lines show the regions of the
parameter space where ratio of the DM relic abundances with and without the
inclusion of superradiance is 15 and 5, respectively. The left and right
panels correspond to $a_{\star}=0.999$ and $a_{\star}=0.6$. As expected, the
maximum boost occurs for high values of the PBH spin and values of
$\alpha\sim\mathcal{O}(1)$. Additionally, the red bands are excluded by the
CMB or BBN data. For $a_{\star}=0.999$, superradiance increases the relic
abundance by at most a factor of $\sim 20$, for a large region of the
parameter space, where the conditions are met. For initial PBH masses
$M_{\text{in}}\lesssim 10~{}{\rm g}$ and DM masses $m_{\text{dm}}\gtrsim
10^{12}~{}{\rm GeV}$, the region where there is an enhancement due to
superradiance shrinks, since the depletion of angular momentum via Hawking
radiation is fast enough so that the instability does not have time to
develop. In the case of heavier PBHs and lighter DM, superradiance can
increase the relic abundance even if the gravitational coupling is small. For
example, for $M_{\text{in}}=10^{7}$ g and $m_{\text{dm}}=1.5\times 10^{5}$ GeV
(corresponding to $\alpha\simeq 6\times 10^{-3}$), the relic abundance is
enhanced by a factor $\sim 7$. Furthermore, we observe that the region where
the increase in the relic abundance is $\gtrsim 5$ grows for larger initial
PBH masses. In this case, PBHs have longer lifetimes that allow for the
development of the instability, although at a lower rate since $\alpha\lesssim
0.1$. As expected, for PBHs rotating with a smaller angular momentum, the
increment on the final relic abundance is reduced. For $a_{\star}=0.6$ (cf.
the right panel of Fig. 3), we find an increase by a factor of at most $\sim
7$. Moreover, we observe that for initial PBH masses smaller than $\sim 70$ g,
the superradiant instability does not develop (as implied by Eq. (3.19)), and
the final relic abundance comes only from Hawking evaporation.

##  4 Gravitational UV Freeze-in

Independently from the PBH evaporation, there is an irreducible DM production
channel which is particularly efficient in the gray region in Fig. 3 (i.e.
large values for the DM mass and the reheating temperature), and corresponds
to the gravitational UV freeze-in [84, 85, 86, 35]. DM can be generated by
2-to-2 annihilations of SM particles, mediated by the exchange of massless
gravitons in the $s$ channel. We emphasize that this gravitational production
mechanism is unavoidable due to universal couplings between the metric and the
energy-momentum tensor $\sim g_{\mu\nu}\,h^{\mu\nu}$.

The evolution of the DM number density $n$ from gravitational freeze-in is
given by the Boltzmann equation

$\frac{dn}{dt}+3\,H\,n=\gamma\,,$ (4.1)

where $\gamma$ corresponds to the DM production rate density, given by [84,
85, 35]

$\gamma(T)=\delta\,\frac{T^{8}}{M_{P}^{4}}\,,$ (4.2)

where $\delta\simeq 1.1\times 10^{-3}$ for scalar DM. Light DM particles
(lighter than the reheating temperature) are mainly produced after the end of
reheating. The DM yield from gravitational freeze-in $Y\equiv n(T)/s(T)$, with
entropy density $s(T)=\frac{2\pi^{2}}{45}g_{\star s}T^{3}$, is obtained by
integrating Eq. (4.1),

$Y_{0}=\frac{45\,\delta}{2\pi^{3}g_{\star
s}}\sqrt{\frac{10}{g_{\star}}}\left(\frac{T_{\text{rh}}}{M_{P}}\right)^{3},\qquad\text{
for }m_{\text{dm}}\ll T_{\text{rh}}\,.$ (4.3)

It is interesting to note that away from the instantaneous decay
approximation, the DM yield is typically only boosted by a small factor of
order $\mathcal{O}(1)$ [87, 88, 89].

Alternatively, in the case where DM is heavier than the reheating temperature
(but still lighter than the highest temperature reached by the SM thermal bath
[90]), it could not be produced after but during the reheating era. Equation
(4.1) therefore yields

$Y_{0}=\frac{45\,\delta}{2\pi^{3}g_{\star
s}}\sqrt{\frac{10}{g_{\star}}}\frac{T_{\text{rh}}^{7}}{M_{P}^{3}\,m_{\text{dm}}^{4}}\,,\qquad\text{
for }m_{\text{dm}}\gg T_{\text{rh}}\,.$ (4.4)

We would like to stress that once produced by the UV freeze-in gravitational
mechanism, the DM could suffer from a dilution due to the PBH evaporation, in
the case where they manage to dominate the total energy density of the
universe (i.e., if $\beta>\beta_{c}$). The total entropy dilution is given by
[67]

$\frac{S(T_{\text{in}})}{S(T_{\text{ev}})}\simeq\begin{dcases}1&\text{ for
}\beta<\beta_{c}\,,\\\
\frac{T_{\text{ev}}}{T_{\text{eq}}}\simeq\left(\frac{g_{\star}}{640}\right)^{\frac{1}{4}}\left(\frac{M_{P}}{M_{\text{in}}}\right)^{\frac{3}{2}}\frac{M_{P}}{\beta\,T_{\text{in}}}&\text{
for }\beta>\beta_{c}\,,\end{dcases}$ (4.5)

with the entropy defined as $S(T)=a^{3}s(T)$. To reproduce the observed DM
relic density, one requires

$m_{\text{dm}}\,Y_{0}\,\frac{S(T_{\text{in}})}{S(T_{\text{ev}})}=\Omega_{\text{dm}}h^{2}\,\frac{1}{s_{0}}\,\frac{\rho_{c}}{h^{2}}\simeq
4.3\times 10^{-10}~{}\text{GeV},$ (4.6)

with $\rho_{c}\simeq 1.05\times 10^{-5}\,h^{2}$ GeV/cm3 being the critical
energy density, $s_{0}\simeq 2.69\times 10^{3}$ cm-3 the entropy density at
present [91], and $\Omega_{\text{dm}}h^{2}\simeq 0.12$ the observed DM relic
abundance [92].

The regions where the gravitational UV freeze-in is effective and can
reproduce and even overshoot the whole observed DM abundance are overlaid in
gray in Fig. 3. Here, two assumptions were made: $i)$ as PBHs are assumed to
be formed in a radiation-dominated universe, $T_{\text{rh}}\geq
T_{\text{in}}$. Here we take $T_{\text{rh}}\geq T_{\text{in}}$. And $ii)$ no
entropy dilution was considered, which can be achieved if $\beta<\beta_{c}$.

##  5 Results

In this section, we summarize the comprehensive parameter space by taking into
account the interplay of the aforementioned three mechanisms. Figure 4 shows
with thick colored lines the parameter space $[M_{\text{in}},\,\beta]$
reproducing the entire observed DM density for different DM masses:
$m_{\text{dm}}=10^{6}$ GeV (black), $10^{9}$ GeV (blue), $10^{12}$ GeV (red)
and $10^{15}$ GeV (green), and four initial values for the spin of PBH
$a_{\star}=0$ (top left panel), $a_{\star}=0.3$ (top right panel),
$a_{\star}=0.6$ (lower left panel) and $a_{\star}=0.999$ (lower right panel).
Different line styles correspond to Hawking radiation only (dotted), Hawking
radiation and superradiance (dashed), and Hawking radiation, gravitational
freeze-in plus superradiance (thick solid). The shaded regions represent areas
constrained by different observables: $M_{\text{in}}\lesssim 10^{-1}$ g, cf.
Eq. (2.7), is disfavored by CMB, while large values for $\beta$ lead to
modifications of BBN predictions due to the energy density stored in GWs, cf.
Eq. (2.8). For the case of the BBN, we present the bound obtained by studying
in detail the evolution of the proton-to-neutron ratio and primordial
abundances of light elements, including the evaporation of PBH, taken from
Ref. [23]. Finally, the dotted line ($\beta=\beta_{c}$) shows the transition
between radiation (lower part) and PBH-dominated eras (upper part).

Figure 4: Parameter space that reproduces the whole observed DM abundance, for
Hawking radiation only (dotted), Hawking radiation and superradiance (dashed),
and Hawking radiation, superradiance and gravitational freeze-in (solid). The
colors correspond to different DM masses $m_{\text{dm}}=10^{6}$ GeV (black),
$10^{9}$ GeV (blue), $10^{12}$ GeV (red) and $10^{15}$ GeV (green). The
different panel correspond to $a_{\star}=0.999$, 0.6, 0.3 and 0. The shaded
areas are excluded by different observables described in the text. We assumed
that $T_{\text{rh}}=T_{\text{in}}$.

Figure 4 can be understood as follows. First, the dotted lines correspond to
the case where DM is only produced by Hawking radiation, and therefore this
requirement shows three different slopes: $i)$ If PBHs dominated the energy
density of the universe (for $\beta\gg\beta_{c}$), the DM yield is independent
of $\beta$, and therefore the lines are vertical. However, a dependence of
$\beta$ appears if, during the entire lifetime of BH, the universe was always
radiation dominated (for $\beta\ll\beta_{c}$). $ii)$ If DM is lighter than the
initial BH temperature $\beta\propto T_{\text{in}}$, while in the opposite
case $iii)$ $\beta\propto T_{\text{in}}^{-3}$ [43, 36, 37]. In the present
case where DM self-interactions are not efficient, DM has to be heavier than
$\mathcal{O}(1)$ MeV in order not to be hot erasing seeds for structure
formation [40, 43, 41], and can be as heavy as $M_{P}$ [93, 94].555Note that
we are not considering a BH evaporation process stopping at $T_{\text{BH}}\sim
M_{P}$, with the associated production of Planck mass relics [95, 96, 54, 97,
98, 39].

However, let us remember that such light PBHs are assumed to be created after
the end of the reheating process, when the universe was SM-radiation
dominated, at a temperature $T_{\text{in}}\leq T_{\text{rh}}$. The
gravitational freeze-in production has a strong dependence on the reheating
temperature, so that it can be minimized by taking
$T_{\text{in}}=T_{\text{rh}}$, which is the assumption used hereafter. In that
sense, this analysis is conservative, as we assume minimal possible
gravitational freeze-in DM production. In Fig. 4 the dashed lines show the
impact on the DM relic abundance once Hawking radiation and gravitational UV
freeze-in are simultaneously taken into account. The total abundance of DM can
be produced by UV freeze-in for $M_{\text{in}}\lesssim 10^{5}$ g and large
values of DM masses, $m_{\text{dm}}\gtrsim 10^{4}$ GeV, cf. Fig. 3. Therefore,
a strong reduction in the $\beta$ parameter is needed to avoid DM
overproduction.

Finally, the solid lines in Fig. 4 take into account the three DM production
mechanisms presented previously: Hawking radiation, gravitational freeze-in,
and superradiance. Superradiance is effective when the gravitational coupling
$\alpha$ is close to unity. As mentioned earlier and also shown in the right
panel of Fig. 2, there is a boost effect in the final DM number density once
superradiance gets triggered. This must be balanced with a reduction of
$\beta$ in the regime when superradiant instability is efficient. Since the
superradiance rate $\Gamma_{\text{sr}}\propto a_{\star}$, cf. Eq. (3.12), it
is expected that the larger $a_{\star}$, the greater the reduction in $\beta$,
which would be zero in the case $a_{\star}=0$ (upper left panel of Fig. 4).

##  6 Conclusions

When the Compton wave length of a particle is comparable to the size of a
spinning black hole (BH), energy can be efficiently transferred to the
particle, leading to an amplification of occupation number in the ergoregion.
For astrophysical BHs with mass $M_{\text{BH}}\geq M_{\odot}$, it is possible
to form clouds around BHs, which could give rise to some observational
signatures; this has been utilized as a potential new avenue to constrain
ultralight particles beyond the standard model of particle physics [1].

In this work, complementary to the astrophysical BHs as widely discussed in
literature, we pay particular attention to the phenomenology of the
superradiant instability with primordial black holes (PBHs) lighter than $\sim
10^{9}$ g. In particular, we investigate the impact of superradiance on the
generation of heavy ($m_{\text{dm}}\gtrsim 1$ TeV) dark matter (DM) particles
in the early Universe. To this end we minimally set up a framework with the
inevitable scenario where DM is generated pure gravitationally via either UV
freeze-in or Hawking emission. Due to the exponential enhancement effect, we
find that in the regime where superradiance manifests, the required BH
abundance $\beta$ to yield the correct DM relics is reduced.

Figure 5 summarizes our main conclusions, where we show the parameter space
that reproduces the whole observed DM abundance, taking into account all the
gravitational production channels described previously (Hawking radiation,
superradiance and gravitational freeze-in). The thickness of the bands
brackets all possible values of the PBH spin: $0\leq a_{\star}<1$. The color
code is the same used before, i.e. $m_{\text{dm}}=10^{3}$ GeV (orange),
$10^{6}$ GeV (black), $10^{9}$ GeV (blue), $10^{12}$ GeV (red) and $10^{15}$
GeV (green). We find that for a close-to-maximally rotating PBHs
($a_{\star}=0.999$) the relic abundance is increased by a factor of $\sim 20$,
when the gravitational coupling between the DM and the PBH is $\sim{\cal
O}(1)$. For lower values of $a_{\star}$, the superradiant effect decreases and
the observed relic abundance can be explained by purely Hawking evaporation.

Figure 5: Parameter space that reproduces the whole observed DM abundance,
taking the three contributions: Hawking radiation, superradiance and
gravitational freeze-in. The colors correspond to different DM masses
$m_{\text{dm}}=10^{3}$ GeV (orange), $10^{6}$ GeV (black), $10^{9}$ GeV
(blue), $10^{12}$ GeV (red) and $10^{15}$ GeV (green). The thickness of the
bands corresponds to the different PBH spins: $0\leq a_{\star}\leq 0.999$. The
shaded areas are excluded by different observables described in the text.

## Acknowledgments

The authors would like to thank Manuel Drees for discussions and comments on
the draft. We also would like to thank Jamie Mc Donald for his comments on the
first version of this manuscript. NB received funding from the Patrimonio
Autónomo - Fondo Nacional de Financiamiento para la Ciencia, la Tecnología y
la Innovación Francisco José de Caldas (MinCiencias - Colombia) grants
80740-465-2020 and 80740-492-2021. NB is also funded by the Spanish
FEDER/MCIU-AEI under grant FPA2017-84543-P. This project has received
funding/support from the European Union’s Horizon 2020 research and innovation
program under the Marie Skłodowska-Curie grant agreement No 860881-HIDDeN.
YFPG would like to thank the Instituto de Fisica Teorica (IFT UAM-CSIC) in
Madrid for support via the Centro de Excelencia Severo Ochoa Program under
Grant CEX2020-001007-S, during the Extended Workshop “Neutrino Theories”,
where part of this work was completed. This work has made use of the Hamilton
HPC Service of Durham University.

##  Appendix A Massive Scalar Fields around Kerr Black Holes

We revisit the solution of a massive scalar $\phi$ in the vicinity of a Kerr
BH and further present the approximation of the superradiance rate
$\Gamma_{\text{sr}}$. We closely follow Ref. [34] where $G=\hbar=c=1$ is used,
however we will restore the units at the end. The Klein-Gordon equation for a
massive scalar $\phi$ with mass $\mu$ in a Kerr metric, defined in Eq. (3.1),
is given by

$\left(g_{\mu\nu}\nabla^{\mu}\nabla^{\nu}-\mu^{2}\right)\phi=0\,.$ (A.1)

Separating variables in a standard manner,

$\phi(t,r,\theta,\varphi)=\sum_{\omega,n,l.m}e^{-i\omega
t+im\varphi}S_{lm}(\theta)R_{nlm}(r)\,,$ (A.2)

then one has

$\displaystyle\frac{1}{\sin\theta}\partial_{\theta}\left(\sin\theta\,\partial_{\theta}S\right)+\left[a^{2}_{\star}r_{g}^{2}\left(\omega^{2}-\mu^{2}\right)\cos^{2}\theta-\frac{m^{2}}{\sin^{2}\theta}+\lambda\right]S=0\,,$
(A.3)
$\displaystyle\Delta\partial_{r}\left(\partial_{r}R\right)-\Delta\left[\mu^{2}r^{2}+a^{2}_{\star}r_{g}^{2}\omega^{2}-2\omega\,m\,a_{\star}r_{g}\,r+\left(\omega\left(r^{2}+a^{2}_{\star}r_{g}^{2}\right)-m\,a_{\star}r_{g}\right)^{2}+\lambda\right]R=0\,,$
(A.4)

where $\lambda$ denotes a separation constant, which is the eigenvalue of Eq.
(A.3). Here our main interest is the solution of the radial component $R$, Eq.
(A.4), with a particular focus on the behavior of the frequency $\omega$
therein. Note that if the imaginary part $\mathfrak{Im}(\omega)<0$, the wave
would decay exponentially. However, the amplitude of the field grows
exponentially if $\mathfrak{Im}(\omega)>0$.

Generically, there are no analytical solutions for Eq. (A.4). However, we can
consider the limit $\mu M\ll 1$, where $\lambda\simeq l(l+1)$ [8]. With a
dimensionless coordinate $x=(r-r_{+})/r_{+}$, it turns out that there are
analytical solutions in two limit regimes with near-horizon region $x\ll
l/\bar{\omega}$ and far region $x\gg 1$, where we have defined
$\bar{\omega}=r_{+}\,\omega$. One can then match the two solutions to obtain
the field behavior in the overlapping region, $1\ll x\ll l/\bar{\omega}$,
which includes the ergosphere. There, an exponential growing mode appears if
$\mathfrak{Im}(\omega)$ is positive. This matching procedure leads to a
constraint among the frequency $\omega$, BH mass $M$, as well as the quantum
numbers, from which one can solve for $\mathfrak{Im}(\omega)$ [34]

$\mathfrak{Im}(\omega
M)=-\frac{1}{2}A_{nl}\left(\frac{\omega^{\prime}M}{\tau}\right)(\mu
M)^{4l+5}\left(\frac{r_{+}-r_{-}}{r_{+}+r_{-}}\right)^{2l+1},$ (A.5)

where $\omega^{\prime}=(2-\tau)(\omega-m\Omega)$, $\tau=(r_{+}-r_{-})/r_{+}$
and

$A_{nl}=\left(\frac{l!}{(2l+1)!(2l)!}\right)^{2}\frac{(l+n)!}{(n-l-1)!}\frac{4^{2l+2}}{n^{2l+4}}\prod_{k=1}^{l}\left[k^{2}+16\left(\frac{\omega^{\prime}M}{\tau}\right)^{2}\right].$
(A.6)

For the mode of interest with $n=2$ and $m=l=1$, one has

$A_{21}=\frac{1}{6}\left[1+16\left(\frac{M(a_{\star}-2\omega\,r_{+})}{2r_{+}\sqrt{1-a_{\star}^{2}}}\right)^{2}\right]\,,$
(A.7)

and

$\displaystyle\mathfrak{Im}(\omega M)$
$\displaystyle=-\frac{1}{2}A_{21}\left(\frac{\omega^{\prime}M}{\tau}\right)(\mu
M)^{9}\left(\frac{r_{+}-r_{-}}{r_{+}+r_{-}}\right)^{3}$
$\displaystyle=\frac{1}{24}(\mu
M)^{9}(a_{\star}-2\omega\,r_{+})\frac{(1-a_{\star}^{2})}{\left(1+\sqrt{1-a_{\star}^{2}}\right)}+\mathcal{O}\left(\alpha-\frac{1}{2}\right)^{3}.$
(A.8)

Note that in the Schwarzschild limit with $a_{\star}\to 0$, one has $r_{+}\to
0$, and then $\mathfrak{Im}(\omega M)$ vanishes, as expected. At leading order
$\omega\sim\mu$, one has

$a_{\star}-2\omega r_{+}\simeq a_{\star}-2\mu
r_{+}\simeq-2\left(1+\sqrt{1-a_{\star}^{2}}\right)\left(\alpha-\frac{1}{2}\right),$
(A.9)

where we have considered $\Omega\sim 1/(2\,r_{g})$, cf. Eq. (3.6). Due to the
superradiance condition, $\alpha<1/2$, higher order terms of
$\left(\alpha-1/2\right)$ in Eq. (A) can be omitted. For $n=2$ and $m=l=1$,
one has

$\mathfrak{Im}(\omega M)\simeq\frac{1}{24}(\mu
M)^{9}(a_{\star}-2\omega\,r_{+})\,,$ (A.10)

where
$(1-a_{\star}^{2})/\left(1+\sqrt{1-a_{\star}^{2}}\right)\sim\mathcal{O}(1)$
has been assumed.666We note that this approximation fails for large values of
$a_{\star}$. Defining the superradiance growth rate as
$\Gamma_{\text{sr}}\equiv\mathfrak{Im}(\omega)$, one has

$\Gamma_{\text{sr}}=\frac{\mu}{24}(\mu M)^{8}(a_{\star}-2\mu\,r_{+})\,,$
(A.11)

which corresponds to Eq. (3.12) once restoring $G^{-1}=8\pi M_{P}^{2}$.

## References

  * [1] R. Brito, V. Cardoso and P. Pani, _Superradiance: New Frontiers in Black Hole Physics_ , _Lect. Notes Phys._ 906 (2015) pp.1 [1501.06570].
  * [2] S.A. Teukolsky and W.H. Press, _Perturbations of a rotating black hole. III - Interaction of the hole with gravitational and electromagnet ic radiation_ , _Astrophys. J._ 193 (1974) 443.
  * [3] R. Penrose, _Gravitational collapse: The role of general relativity_ , _Riv. Nuovo Cim._ 1 (1969) 252.
  * [4] W.H. Press and S.A. Teukolsky, _Floating Orbits, Superradiant Scattering and the Black-hole Bomb_ , _Nature_ 238 (1972) 211.
  * [5] V. Cardoso, O.J.C. Dias, J.P.S. Lemos and S. Yoshida, _The Black hole bomb and superradiant instabilities_ , _Phys. Rev. D_ 70 (2004) 044039 [hep-th/0404096].
  * [6] R. Brito and P. Pani, _Black-Hole Superradiance: Searching for Ultralight Bosons with Gravitational Waves_ , in _Handbook of Gravitational Wave Astronomy_ (2021), DOI.
  * [7] T.J.M. Zouros and D.M. Eardley, _Instabilities of Massive Scalar Perturbations of a Rotating Black Hole_ , _Annals Phys._ 118 (1979) 139.
  * [8] S.L. Detweiler, _Klein-Gordon Equation and Rotating Black Holes_ , _Phys. Rev. D_ 22 (1980) 2323.
  * [9] S.R. Dolan, _Instability of the massive Klein-Gordon field on the Kerr spacetime_ , _Phys. Rev. D_ 76 (2007) 084001 [0705.2880].
  * [10] A. Arvanitaki, S. Dimopoulos, S. Dubovsky, N. Kaloper and J. March-Russell, _String Axiverse_ , _Phys. Rev. D_ 81 (2010) 123530 [0905.4720].
  * [11] A. Arvanitaki, M. Baryakhtar and X. Huang, _Discovering the QCD Axion with Black Holes and Gravitational Waves_ , _Phys. Rev. D_ 91 (2015) 084011 [1411.2263].
  * [12] V. Cardoso, Ó.J.C. Dias, G.S. Hartnett, M. Middleton, P. Pani and J.E. Santos, _Constraining the mass of dark photons and axion-like particles through black-hole superradiance_ , _JCAP_ 03 (2018) 043 [1801.01420].
  * [13] M.J. Stott and D.J.E. Marsh, _Black hole spin constraints on the mass spectrum and number of axionlike fields_ , _Phys. Rev. D_ 98 (2018) 083006 [1805.02016].
  * [14] S. Endlich and R. Penco, _A Modern Approach to Superradiance_ , _JHEP_ 05 (2017) 052 [1609.06723].
  * [15] M. Baryakhtar, R. Lasenby and M. Teo, _Black Hole Superradiance Signatures of Ultralight Vectors_ , _Phys. Rev. D_ 96 (2017) 035019 [1704.05081].
  * [16] V.P. Frolov, P. Krtouš, D. Kubizňák and J.E. Santos, _Massive Vector Fields in Rotating Black-Hole Spacetimes: Separability and Quasinormal Modes_ , _Phys. Rev. Lett._ 120 (2018) 231103 [1804.00030].
  * [17] R. Brito, V. Cardoso and P. Pani, _Massive spin-2 fields on black hole spacetimes: Instability of the Schwarzschild and Kerr solutions and bounds on the graviton mass_ , _Phys. Rev. D_ 88 (2013) 023514 [1304.6725].
  * [18] R. Brito, S. Grillo and P. Pani, _Black Hole Superradiant Instability from Ultralight Spin-2 Fields_ , _Phys. Rev. Lett._ 124 (2020) 211101 [2002.04055].
  * [19] M. Baryakhtar, M. Galanis, R. Lasenby and O. Simon, _Black hole superradiance of self-interacting scalar fields_ , _Phys. Rev. D_ 103 (2021) 095019 [2011.11646].
  * [20] H. Davoudiasl and P.B. Denton, _Ultralight Boson Dark Matter and Event Horizon Telescope Observations of M87*_ , _Phys. Rev. Lett._ 123 (2019) 021102 [1904.09242].
  * [21] H. Davoudiasl, P.B. Denton and J. Gehrlein, _Supermassive Black Holes, Ultralight Dark Matter, and Gravitational Waves from a First Order Phase Transition_ , _Phys. Rev. Lett._ 128 (2022) 081101 [2109.01678].
  * [22] Y. Chen, R. Roy, S. Vagnozzi and L. Visinelli, _Superradiant evolution of the shadow and photon ring of Sgr A ∗_, 2205.06238.
  * [23] B.J. Carr, K. Kohri, Y. Sendouda and J. Yokoyama, _New cosmological constraints on primordial black holes_ , _Phys. Rev. D_ 81 (2010) 104019 [0912.5297].
  * [24] B. Carr, K. Kohri, Y. Sendouda and J. Yokoyama, _Constraints on primordial black holes_ , _Rept. Prog. Phys._ 84 (2021) 116902 [2002.12778].
  * [25] Y.N. Eroshenko, _Spin of primordial black holes in the model with collapsing domain walls_ , _JCAP_ 12 (2021) 041 [2111.03403].
  * [26] T. Harada, C.-M. Yoo, K. Kohri and K.-I. Nakao, _Spins of primordial black holes formed in the matter-dominated phase of the Universe_ , _Phys. Rev. D_ 96 (2017) 083517 [1707.03595].
  * [27] F. Kuhnel, _Enhanced Detectability of Spinning Primordial Black Holes_ , _Eur. Phys. J. C_ 80 (2020) 243 [1909.04742].
  * [28] M.M. Flores and A. Kusenko, _Spins of primordial black holes formed in different cosmological scenarios_ , _Phys. Rev. D_ 104 (2021) 063008 [2106.03237].
  * [29] C.M. Chambers, W.A. Hiscock and B. Taylor, _Spinning down a black hole with scalar fields_ , _Phys. Rev. Lett._ 78 (1997) 3249 [gr-qc/9703018].
  * [30] B.E. Taylor, C.M. Chambers and W.A. Hiscock, _Evaporation of a Kerr black hole by emission of scalar and higher spin particles_ , _Phys. Rev. D_ 58 (1998) 044012 [gr-qc/9801044].
  * [31] M. Calzà, J. March-Russell and J.a.G. Rosa, _Evaporating primordial black holes, the string axiverse, and hot dark radiation_ , 2110.13602.
  * [32] S.W. Hawking, _Particle Creation by Black Holes_ , _Commun. Math. Phys._ 43 (1975) 199.
  * [33] B.J. Carr, _Some cosmological consequences of primordial black-hole evaporations_ , _Astrophys. J._ 206 (1976) 8.
  * [34] P.B. Ferraz, T.W. Kephart and J.a.G. Rosa, _Superradiant pion clouds around primordial black holes_ , 2004.11303.
  * [35] N. Bernal, M. Dutra, Y. Mambrini, K. Olive, M. Peloso and M. Pierre, _Spin-2 Portal Dark Matter_ , _Phys. Rev. D_ 97 (2018) 115020 [1803.01866].
  * [36] N. Bernal and Ó. Zapata, _Gravitational dark matter production: primordial black holes and UV freeze-in_ , _Phys. Lett. B_ 815 (2021) 136129 [2011.02510].
  * [37] N. Bernal, _Gravitational Dark Matter and Primordial Black Holes_ , in _Beyond Standard Model: From Theory to Experiment_ , 5, 2021 [2105.04372].
  * [38] O. Lennon, J. March-Russell, R. Petrossian-Byrne and H. Tillim, _Black Hole Genesis of Dark Matter_ , _JCAP_ 04 (2018) 009 [1712.07664].
  * [39] D. Hooper, G. Krnjaic and S.D. McDermott, _Dark Radiation and Superheavy Dark Matter from Black Hole Domination_ , _JHEP_ 08 (2019) 001 [1905.01301].
  * [40] I. Masina, _Dark matter and dark radiation from evaporating primordial black holes_ , _Eur. Phys. J. Plus_ 135 (2020) 552 [2004.04740].
  * [41] I. Baldes, Q. Decant, D.C. Hooper and L. Lopez-Honorez, _Non-Cold Dark Matter from Primordial Black Hole Evaporation_ , _JCAP_ 08 (2020) 045 [2004.14773].
  * [42] P. Gondolo, P. Sandick and B. Shams Es Haghi, _Effects of primordial black holes on dark matter models_ , _Phys. Rev. D_ 102 (2020) 095018 [2009.02424].
  * [43] N. Bernal and Ó. Zapata, _Self-interacting Dark Matter from Primordial Black Holes_ , _JCAP_ 03 (2021) 007 [2010.09725].
  * [44] N. Bernal and Ó. Zapata, _Dark Matter in the Time of Primordial Black Holes_ , _JCAP_ 03 (2021) 015 [2011.12306].
  * [45] I. Masina, _Dark Matter and Dark Radiation from Evaporating Kerr Primordial Black Holes_ , _Grav. Cosmol._ 27 (2021) 315 [2103.13825].
  * [46] S. Jyoti Das, D. Mahanta and D. Borah, _Low scale leptogenesis and dark matter in the presence of primordial black holes_ , _JCAP_ 11 (2021) 019 [2104.14496].
  * [47] A. Cheek, L. Heurtier, Y.F. Perez-Gonzalez and J. Turner, _Primordial black hole evaporation and dark matter production. I. Solely Hawking radiation_ , _Phys. Rev. D_ 105 (2022) 015022 [2107.00013].
  * [48] A. Cheek, L. Heurtier, Y.F. Perez-Gonzalez and J. Turner, _Primordial black hole evaporation and dark matter production. II. Interplay with the freeze-in or freeze-out mechanism_ , _Phys. Rev. D_ 105 (2022) 015023 [2107.00016].
  * [49] P. Sandick, B.S. Es Haghi and K. Sinha, _Asymmetric reheating by primordial black holes_ , _Phys. Rev. D_ 104 (2021) 083523 [2108.08329].
  * [50] B.J. Carr and S.W. Hawking, _Black holes in the early Universe_ , _Mon. Not. Roy. Astron. Soc._ 168 (1974) 399.
  * [51] P. Villanueva-Domingo, O. Mena and S. Palomares-Ruiz, _A brief review on primordial black holes as dark matter_ , _Front. Astron. Space Sci._ 8 (2021) 87 [2103.12087].
  * [52] B. Carr and F. Kuhnel, _Primordial Black Holes as Dark Matter: Recent Developments_ , _Ann. Rev. Nucl. Part. Sci._ 70 (2020) 355 [2006.02838].
  * [53] S.W. Hawking, I.G. Moss and J.M. Stewart, _Bubble Collisions in the Very Early Universe_ , _Phys. Rev. D_ 26 (1982) 2681.
  * [54] B.J. Carr, J.H. Gilbert and J.E. Lidsey, _Black hole relics and inflation: Limits on blue perturbation spectra_ , _Phys. Rev. D_ 50 (1994) 4853 [astro-ph/9405027].
  * [55] J. Garcia-Bellido, A.D. Linde and D. Wands, _Density perturbations and black hole formation in hybrid inflation_ , _Phys. Rev. D_ 54 (1996) 6040 [astro-ph/9605094].
  * [56] J. Yokoyama, _Cosmological constraints on primordial black holes produced in the near critical gravitational collapse_ , _Phys. Rev. D_ 58 (1998) 107502 [gr-qc/9804041].
  * [57] J.C. Niemeyer and K. Jedamzik, _Dynamics of primordial black hole formation_ , _Phys. Rev. D_ 59 (1999) 124013 [astro-ph/9901292].
  * [58] M. Drees and E. Erfani, _Running-Mass Inflation Model and Primordial Black Holes_ , _JCAP_ 04 (2011) 005 [1102.2340].
  * [59] K. Kohri, C.-M. Lin and T. Matsuda, _Primordial black holes from the inflating curvaton_ , _Phys. Rev. D_ 87 (2013) 103527 [1211.2371].
  * [60] F. Kühnel, C. Rampf and M. Sandstad, _Effects of Critical Collapse on Primordial Black-Hole Mass Spectra_ , _Eur. Phys. J. C_ 76 (2016) 93 [1512.00488].
  * [61] H. Deng and A. Vilenkin, _Primordial black hole formation by vacuum bubbles_ , _JCAP_ 12 (2017) 044 [1710.02865].
  * [62] N. Bernal, V. Muñoz Albornoz, S. Palomares-Ruiz and P. Villanueva-Domingo, _Current and future neutrino limits on the abundance of primordial black holes_ , 2203.14979.
  * [63] A. Dolgov and J. Silk, _Baryon isocurvature fluctuations at small scales and baryonic dark matter_ , _Phys. Rev. D_ 47 (1993) 4244.
  * [64] A.M. Green, _Microlensing and dynamical constraints on primordial black hole dark matter with an extended mass function_ , _Phys. Rev. D_ 94 (2016) 063530 [1609.01143].
  * [65] K. Kannike, L. Marzola, M. Raidal and H. Veermäe, _Single Field Double Inflation and Primordial Black Holes_ , _JCAP_ 09 (2017) 020 [1705.06225].
  * [66] R. Allahverdi et al., _The First Three Seconds: a Review of Possible Expansion Histories of the Early Universe_ , 2006.16182.
  * [67] N. Bernal, F. Hajkarim and Y. Xu, _Axion Dark Matter in the Time of Primordial Black Holes_ , _Phys. Rev. D_ 104 (2021) 075007 [2107.13575].
  * [68] N. Bernal, Y.F. Perez-Gonzalez, Y. Xu and Ó. Zapata, _ALP dark matter in a primordial black hole dominated universe_ , _Phys. Rev. D_ 104 (2021) 123536 [2110.04312].
  * [69] S. Sarkar, _Big bang nucleosynthesis and physics beyond the standard model_ , _Rept. Prog. Phys._ 59 (1996) 1493 [hep-ph/9602260].
  * [70] M. Kawasaki, K. Kohri and N. Sugiyama, _MeV scale reheating temperature and thermalization of neutrino background_ , _Phys. Rev. D_ 62 (2000) 023506 [astro-ph/0002127].
  * [71] S. Hannestad, _What is the lowest possible reheating temperature?_ , _Phys. Rev. D_ 70 (2004) 043506 [astro-ph/0403291].
  * [72] F. De Bernardis, L. Pagano and A. Melchiorri, _New constraints on the reheating temperature of the universe after WMAP-5_ , _Astropart. Phys._ 30 (2008) 192.
  * [73] P.F. de Salas, M. Lattanzi, G. Mangano, G. Miele, S. Pastor and O. Pisanti, _Bounds on very low reheating scenarios after Planck_ , _Phys. Rev. D_ 92 (2015) 123534 [1511.00672].
  * [74] T. Hasegawa, N. Hiroshima, K. Kohri, R.S.L. Hansen, T. Tram and S. Hannestad, _MeV-scale reheating temperature and thermalization of oscillating neutrinos by radiative and hadronic decays of massive particles_ , _JCAP_ 12 (2019) 012 [1908.10189].
  * [75] Planck collaboration, _Planck 2018 results. X. Constraints on inflation_ , _Astron. Astrophys._ 641 (2020) A10 [1807.06211].
  * [76] K. Inomata, M. Kawasaki, K. Mukaida, T. Terada and T.T. Yanagida, _Gravitational Wave Production right after a Primordial Black Hole Evaporation_ , _Phys. Rev. D_ 101 (2020) 123533 [2003.10455].
  * [77] G. Domènech, C. Lin and M. Sasaki, _Gravitational wave constraints on the primordial black hole dominated early universe_ , _JCAP_ 04 (2021) 062 [2012.08151].
  * [78] D.N. Page, _Particle Emission Rates from a Black Hole: Massless Particles from an Uncharged, Nonrotating Hole_ , _Phys. Rev. D_ 13 (1976) 198.
  * [79] E.W. Leaver, _An Analytic representation for the quasi normal modes of Kerr black holes_ , _Proc. Roy. Soc. Lond. A_ 402 (1985) 285.
  * [80] P. Virtanen, R. Gommers, T.E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau et al., _SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python_ , _Nature Methods_ 17 (2020) 261.
  * [81] H. Witek, V. Cardoso, A. Ishibashi and U. Sperhake, _Superradiant instabilities in astrophysical systems_ , _Phys. Rev. D_ 87 (2013) 043513 [1212.0551].
  * [82] R. Roy, S. Vagnozzi and L. Visinelli, _Superradiance evolution of black hole shadows revisited_ , _Phys. Rev. D_ 105 (2022) 083002 [2112.06932].
  * [83] A. Granelli, K. Moffat, Y.F. Perez-Gonzalez, H. Schulz and J. Turner, _ULYSSES: Universal LeptogeneSiS Equation Solver_ , _Comput. Phys. Commun._ 262 (2021) 107813 [2007.09150].
  * [84] M. Garny, M. Sandora and M.S. Sloth, _Planckian Interacting Massive Particles as Dark Matter_ , _Phys. Rev. Lett._ 116 (2016) 101302 [1511.03278].
  * [85] Y. Tang and Y.-L. Wu, _On Thermal Gravitational Contribution to Particle Production and Dark Matter_ , _Phys. Lett. B_ 774 (2017) 676 [1708.05138].
  * [86] M. Garny, A. Palessandro, M. Sandora and M.S. Sloth, _Theory and Phenomenology of Planckian Interacting Massive Particles as Dark Matter_ , _JCAP_ 02 (2018) 027 [1709.09688].
  * [87] M.A.G. Garcia, Y. Mambrini, K.A. Olive and M. Peloso, _Enhancement of the Dark Matter Abundance Before Reheating: Applications to Gravitino Dark Matter_ , _Phys. Rev. D_ 96 (2017) 103510 [1709.01549].
  * [88] N. Bernal, F. Elahi, C. Maldonado and J. Unwin, _Ultraviolet Freeze-in and Non-Standard Cosmologies_ , _JCAP_ 11 (2019) 026 [1909.07992].
  * [89] B. Barman, N. Bernal, Y. Xu and Ó. Zapata, _Ultraviolet Freeze-in with a Time-dependent Inflaton Decay_ , 2202.12906.
  * [90] G.F. Giudice, E.W. Kolb and A. Riotto, _Largest temperature of the radiation era and its cosmological implications_ , _Phys. Rev. D_ 64 (2001) 023508 [hep-ph/0005123].
  * [91] Particle Data Group collaboration, _Review of Particle Physics_ , _PTEP_ 2020 (2020) 083C01.
  * [92] Planck collaboration, _Planck 2018 results. VI. Cosmological parameters_ , _Astron. Astrophys._ 641 (2020) A6 [1807.06209].
  * [93] D.J.H. Chung, E.W. Kolb, A. Riotto and I.I. Tkachev, _Probing Planckian physics: Resonant production of particles during inflation and features in the primordial power spectrum_ , _Phys. Rev. D_ 62 (2000) 043508 [hep-ph/9910437].
  * [94] G.F. Giudice, M. Peloso, A. Riotto and I. Tkachev, _Production of massive fermions at preheating and leptogenesis_ , _JHEP_ 08 (1999) 014 [hep-ph/9905242].
  * [95] J.H. MacGibbon, _Can Planck-mass relics of evaporating black holes close the universe?_ , _Nature_ 329 (1987) 308.
  * [96] J.D. Barrow, E.J. Copeland and A.R. Liddle, _The Cosmology of black hole relics_ , _Phys. Rev. D_ 46 (1992) 645.
  * [97] A.D. Dolgov, P.D. Naselsky and I.D. Novikov, _Gravitational waves, baryogenesis, and dark matter from primordial black holes_ , astro-ph/0009407.
  * [98] D. Baumann, P.J. Steinhardt and N. Turok, _Primordial Black Hole Baryogenesis_ , hep-th/0703250.