Constraints on evaporating primordial black holes from
the AMS-02 positron data

Astrophysical and cosmological observations suggest that the dominant part of the matter in the present Universe is in the form of non-luminous dark matter (DM). In many particle physics models, DM consists of a new type of elementary particle beyond the Standard Model (SM) of particle physics, which may participate non-gravitational interactions with baryonic matter, such as weakly interacting massive particles (WIMPs), sterile neutrinos, QCD axions and ultralight bosonic particles. Despite great experimental efforts in recent decades, so far there is no confirmed signals of particle DM from direct, indirect and collider DM search experiments.

Primordial black hole (PBH) is an alternative DM candidate which does not require any new physics beyond the SM. PBHs are believed to have formed after the inflation, and subsequently evolved through accretion, mergers, and Hawking radiation. Depending on the time of formation, the masses of PBHs can vary in a large range. In general, the initial mass $M_{\text{PBH}}$ of a PBH should be close to the Hubble horizon mass at the production time $t$, $M_{\text{PBH}}\sim c^{3}t/G\simeq 10^{15}(t/10^{-23}\,\text{s})\,\text{g}$, where $c$ is the speed of light and $G$ is the Newton constant. For two typical formation time $t$ of the Planck time $t\sim 10^{-43}\,\mathrm{s}$ and the time just before the big-bang nucleosynthesis (BBN) $t\sim 1\,\mathrm{s}$, the initial PBH mass are around $10^{-5}\,\mathrm{g}$ and $10^{5}\,M_{\odot}$, respectively.

PBH can be whole or a fraction of DM. The energy fraction of PBHs relative to that of whole DM is defined as $f_{\text{PBH}}\equiv\Omega_{\text{PBH}}/\Omega_{\text{DM}}$, where $\Omega_{\text{DM}}$ and $\Omega_{\text{PBH}}$ are the energy density parameters of DM and PBHs relative to the critical density of the present Universe, respectively. For heavy PBHs with $M_{\text{PBH}}\gg 10^{17}\text{g}$, the value of $f_{\text{PBH}}$ can be constrained by the gravitational effects from PBHs such as microlensing MACHO:2000qbb ; Wyrzykowski:2010mh ; Wyrzykowski:2011tr ; Macho:2000nvd ; Griest:2013aaa ; CalchiNovati:2013jpj ; Griest:2013esa , dynamical constraint from globular clusters, galaxy disruption and other observables Koushiappas:2017chw ; Monroy-Rodriguez:2014ula ; Carr:2019bel ; MUSE:2020qbo (for recent reviews, see e.g.Carr:2020gox ; Carr:2020xqk ).

Light PBHs are expected to emit SM particles through Hawking radiations Hawking:1974rv . PBHs in the mass range from $10^{13}-10^{17}\,\text{g}$ are expected to emit SM particles with typical energies from a few GeV down to a few hundreds of keV. PBHs loss their masses through Hawking radiation at a rate $\text{d}M_{\text{PBH}}/\text{dt}\propto M_{\text{PBH}}^{-2}$ Baldes:2020nuv ; Cheek:2021odj , which suggests that heavier PBHs evaporate slower. It has been shown that PBHs with masses $M_{\text{PBH}}\gtrsim 5\times 10^{14}\text{ g}$ have lifetimes larger the age of the Universe, and can be considered as stable sources of photons and cosmic-ray (CR) particles in the Galaxy MacGibbon:1991tj . This type of evaporating PBHs can be searched by current space-borne detectors. For PBHs in this mass region, the value of $f_{\mathrm{PBH}}$ can be constrained by the data of extragalactic and galactic diffuse $\gamma$-rays Carr:2009jm ; Carr:2016hva ; Arbey:2019vqx , CMB Auffinger:2022khh , neutrinos Dasgupta:2019cae ; Wang:2020uvi , CR electrons Boudaud:2018hqb , CR positron annihilation into 511 keV gamma-ray lines Dasgupta:2019cae ; Laha:2019ssq , and CR antiprotons Maki:1995pa ; Barrau:2001ev , etc.

Recently, low-energy CR all-electron ($e^{+}+e^{-}$) flux data from Voyager-1 have been used to set constraints on $f_{\mathrm{PBH}}$ Boudaud:2018hqb . The obtained limits turned out to be competitive with that derived from extragalactic $\gamma$-rays. As Voyager-1 is now outside the heliopuase, the electron flux measured by Voyager-1 can be considered as the true local interstellar (LIS) flux, and the derived constraints are expected to be relatively robust against the influence of the solar activity, the so called solar modulation effect. Note, however, that although the low-energy electron data from Voyager-1 can be considered as the true LIS flux, the theoretical prediction for CR flux from PBH evaporation involves a number of parameters for the CR propagation in the Galaxy, such as the diffusion coefficient, re-acceleration coefficient and convection velocity. These parameters are determined through fitting to the CR data (e.g the Boron to Carbon flux ratio, B/C) measured at the top of the atmosphere (TOA) deep inside the heliospher, which are strongly affected by the solar activity. Thus, the constraints on the PBH fraction from the Voyager data are inevitabally affected the solar modulation effect. Actually, the LIS fluxes measured by Voyager-1,2 are more useful in improving the modeling and calibrating of the solar modulation effect itself. In order to derive robust constraints on exotic contributions, it is necessary to consider both the LIS and TOA CR flux data simultaneously and consistently calculate the CR propagation in the Galaxy and heliosphere, as the parameters of the two processes are strongly correlated.

In this work we explore the possibility of using CR positron flux to constrain the abundance of PBHs. CR positrons are believed to be of secondary origin and relatively rare compared with the CR primaries such as protons and electrons. Thus CR positron flux should be sensitive to any exotic contributions. The currently measured CR positron flux already showed an unexpected rise starting at high energies above 20 GeV and a peak at around 300 GeV, which may hint at nearby astrophysical sources or halo DM contributions (see e.g. Bergstrom:2008gr ; Cirelli:2008pk ; Bergstrom:2009fa ; Lin:2014vja ; Jin:2013nta ; Jin:2014ica ). Alternatively, the low-energy CR positron flux below $\sim 10$ GeV can be used to set constraints on exotic contributions as they are roughly consistent with the expected backgrounds. In this work, we shall use the low-energy CR positron data with kinetic energy from 0.5–5.2 GeV from AMS-02 to constrain the PBH fraction $f_{\text{PBH}}$. Although the initial energy of the CR positrons from evaporating PBHs with mass above $5\times 10^{14}\text{ g}$ should be quite low around $\mathcal{O}(10)$ MeV, the positron energy can be further enhanced during the propagation process through interacting with the galactic magnetic fields, which is known as the re-acceleration process. The existence of re-acceleration is favored by a large number of independent analyses (see, e.g. Strong:1998pw ; Moskalenko:2001ya ; Strong:2001fu ; Moskalenko:2002yx ; Ptuskin:2005ax ; Jin:2014ica ; Johannesson:2016rlh ; Boschini:2019gow ; Boschini:2020jty ; Evoli:2016xgn ). We use the state-of-the-art models for CR propagation in the Galaxy (based on the numerical code ${\tt{Galprop}}$ Strong:1998pw ; Moskalenko:2001ya ; Strong:2001fu ; Moskalenko:2002yx ; Ptuskin:2005ax ) and in the heliospher (using the code ${\tt{Helmod}}$ Bobik:2011ig ; Bobik:2016 ; Boschini:2017gic ; Boschini:2019ubh ; Boschini:2022 which is based on the numerical solution of the Parker equation and the force-field approximation). In calculating the CR positron flux from PBH evaporation the numerical code ${\tt{BlackHawk}}$ Auffinger:2022sqj is used which includes both the primary and secondary contributions. The results show that for typical diffusive re-acceleration models with best-fit parameters determined from current proton, B/C and light CR nuclei data, the CR positron flux measured by AMS-02 can provide very stringent limits, which can be stronger than the previous constraints derived from the Voyager electron flux by around an order of magnitude.

The remaining part of this paper is organized as follows. In section II, we give a brief overview on the particle spectrum from PBH evaporation. In section III, we discuss the CR propagation in the Galaxy and the models for solar modulation. The constraints on PBH abundances from the AMS-02 positron data are discussed in section IV for the cases without background and with backgrounds. The results of this work is summarized in section V.

In this work, unless otherwise stated, we adopt the natural system of units with $\hbar=k_{\text{B}}=c=1$, where $\hbar$ is the reduced Plank constant, $k_{\text{B}}$ is the Boltzmann constant, and $c$ is the speed of light. We consider a simple scenario where the spin of PBHs can be negligible, which can be obtained from a series of possible formation mechanisms DeLuca:2019buf ; 1704.06573 ; 1901.05963 (for the formation of PBHs with near-extremal spin, see e.g. Harada:2016mhb ; Harada:2017fjm ; Cotner:2017tir ). The emission rate of particle species $i$ per unit time $t$ and total energy $E$ from Hawking radiation is given by Hawking:1975iha

where $s_{i}$ and $g_{i}$ are the spin and the total degree of freedom of the particle $i$, respectively, and the temperature $T_{\text{PBH}}$ of a PBH with mass $M_{\text{PBH}}$ is Page:1976df

The graybody factor $\Gamma_{i}$ in Eq. (1) is determined by the equation of motion of the particle in curved space time near the horizon. It describes the probability that the particle $i$, created at the PBH horizon, finally escapes to spatial infinity. In the geometric optics limit (i.e., the high-energy limit), the graybody factor for electron can be approximated as $\Gamma_{e}\simeq 27\,G^{2}\,M_{\text{PBH}}^{2}\,E^{2}$. Note that Eq. (1) only describes the primary particles directly emitted from PBHs. The decays of the unstable primary particles can produce secondary stable particles. In the calculation of the graybody factor and the energy spectra of the emitted particles, we use the numerical code of ${\tt{BlackHawk}}$ Arbey:2019mbc ; Auffinger:2022sqj in which both the primary and secondary production processes are calculated. For the low-energy particle production and decay, we use the results from the Hazma code 1907.11846 which is now included in the updated version BlackHawk-v2.1 Auffinger:2022sqj .

The energy spectrum $\text{d}^{2}n_{e^{\pm}}/\text{d}t\,\text{d}E_{\text{kin}}$ from all the evaporating PBHs with different masses is given by

where $\text{d}n/\text{d}M_{\text{PBH}}$ is the mass distribution function of PBH. Depending on formation mechanisms, the mass function can be a peak theory distribution Tashiro:2008sf ; Germani:2018jgr , log-normal distribution Kannike:2017bxn , monochromatic distribution Carr:2017jsz or power-law distribution Carr:2017edp . In this work, we consider two widely used models, namely, monochromatic and log-normal mass functions. A nearly monochromatic mass function is naturally expected if all PBHs are formed at a same epoch, and a log-normal mass function can arise from inflationary fluctuations Dolgov:1992pu ; Clesse:2015wea . The two mass functions are given by

where $M_{\text{c}}$ is the characteristic mass, $\sigma$ is the width of the log-normal mass distribution, $A_{1}$ and $A_{2}$ are normalized factors with different unit, which are determined by $\rho_{\text{PBH}}$ the energy density of PBH as follows

Due to Hawking radiation, PBHs continually lose their masses. Analytical and numerical calculations both confirmed that the lighter PBHs, the faster evaporation MacGibbon:1990zk ; Chao:2021orr , suggesting PBHs with very small masses $(\lesssim 5\times 10^{14}\,\mathrm{g})$ are absent by now. In this work we only take into account the contribution from the existing PBHs, the lower bound $M_{\mathrm{min}}$ in Eq. (5) is fixed at $M_{\text{min}}=5\times 10^{14}\,\mathrm{g}$. From Eq. (5), the relation between the normalized factor and $\rho_{\mathrm{PBH}}$ can be obtained. It is evident that $A_{1}=\rho_{\mathrm{PBH}}/M_{\mathrm{c}}$ for $M_{\mathrm{c}}>M_{\mathrm{min}}$. The value of $A_{2}$ can be expressed as $A_{2}=\rho_{\mathrm{PBH}}/k$, where $0.5\lesssim k\lesssim 1$ from numerical calculations. The case of $k=1$ corresponds to the limit of $M_{\mathrm{min}}=0$.

In the left panel of Fig. 1, we show the calculated emission spectrum $\text{d}^{2}N_{e^{\pm}}/\text{d}t\,\text{d}E_{\text{kin}}$ for CR $e^{\pm}$ as a function of kinetic energy $E_{\text{kin}}=E-m_{e}$ in the monochromatically distributed PBHs with typical masses $M_{\text{PBH}}=5\times 10^{14}\ \text{g}$, $1\times 10^{16}\ \text{g}$ and $5\times 10^{16}\ \text{g}$, respectively. For light PBHs with $M_{\text{PBH}}=5\times 10^{14}\ \text{g}$, the contributions from secondary particles are important, the secondaries can be dominant in the low-energy region below 10 MeV and can change the spectral shape significantly. For heavier PBHs with $M_{\text{PBH}}\gtrsim 1\times 10^{16}\ \text{g}$, the secondary contribution is negligible, so only the total contribution is shown. In the right panel of Fig. 1, the energy spectra for the case of log-normal PBH mass distribution is shown for fixed $M_{\text{c}}=1\times 10^{16}\ \text{g}$ with three different widths of $\sigma=0.5$, 1.0 and 2.0. Compared with the left panel, for larger width such as $\sigma=2.0$, the spectra can extend to much higher energies, as the possibility for having lighter PBHs is significantly increase.

In this section, we derive the constraints on $f_{\text{PBH}}$ from the current CR electron and positron data. In the first step, we derive the constraints under the assumption of null astrophysical backgrounds, which should be more conservative and robust compared with that with astrophysical backgrounds included. In calculating the constraints, we adopt a simple approach by requiring that the predicted CR fluxe from PBH evaporation should not exceed the measured value by more than $2\leavevmode\nobreak\ \sigma$ uncertainty in any energy bin, which was adopted in Ref. Boudaud:2018hqb , and allows for a direct comparison with their results. The constraints can be obtained by simply rescaling the value of $f_{\text{PBH}}$ and compare the predictions with the experimental data.

In Fig. 3, we show the obtained individual constraints from Voyager-1,-2 all-electron data, AMS-02 all-electron data and AMS-02 positron data for monochromatically distributed PBHs in the mass range $5\times 10^{14}-10^{17}\leavevmode\nobreak\ \text{g}$. The results for four different propagation models GH, MIN, MED and MAX listed in Tab. 1 are shown in the four panels of Fig. 3. It can be seen from the figure that for all the four propagation models, the constraints from the low-energy AMS-02 positron data are the most stringent. In the GH propagation model, we obtain conservative upper limit of

which is stronger than that from the Voyager all-electron data by an order of magnitude. The AMS-02 positron constraints are also more stringent than that from the AMS-02 all-electron data. The reason is that the initial spectra of electrons and positrons from PBH evaporation are the same as the evaporation process is charge symmetric. Thus the prediction for all-electron flux from PBH evaporation should be approximately twice the positron flux. However, the all-electron flux measured by AMS-02 is around an order of magnitude higher than of positron flux. Thus the corresponding constraints are significantly weaker in the null background case. The constraints derived in Fig. 3 are in a good agreement with the values previously estimated from the left panel of Fig. 2. In the left panel of Fig. 4 we show the constraints for log-normally distributed PBHs with a typical width $\sigma=1$ in the GH propagation model. Similar to the case with monochromatic mass function, the AMS-02 positron data provide the most stringent limits in the case with log-normal mass function. In the right panel of Fig. 4, the constraints for different widths $\sigma=0.5,1.0$ and 2.0 for the log-normal mass function are compared. For the log-normal mass function constraints becomes more stringent as the value of $\sigma$ increases, which is consistent with the flux predictions in the right panel of Fig. 2.

In the next step, we consider the CR-positron constraints including the secondary positron backgrounds. Including the astrophysical background in general results in more stringent constraints. Note, however, that so far the theoretical predictions for low-energy (sub-GeV) positron flux still suffer from significant uncertainties. In addition to the uncertainties in the propagation model such as the diffusion halo half-height $z_{h}$, large uncertainties arise from the poorly known low-energy $pp$ and $pA$ inelastic scattering cross sections for positron production. At present, the cross sections can be obtained either from analytical parameterizations to the experimental data  Orusa:2022pvp ; Kafexhiu:2014cua ; Kamae:2006bf or from QCD-based Monte-Carlo event generators  DelaTorreLuque:2023zyd ; Koldobskiy:2021nld ; Bierlich:2022pfr . The difference can easily reach a factor of two. Other uncertainties specifically related to CR electron and positron energy loss involve that in the Galactic magnetic fields and gas distribution (for recently analyses, see e.g. DiMauro:2023oqx ; DelaTorreLuque:2023zyd ). Thus the constraints including the astrophysical backgrounds are less robust compared with that without backgrounds. To account for the uncertainties in the secondary positron flux , we add two additional parameters $N_{e^{+}}$ and $\delta_{e^{+}}$ to the calculated CR positron flux $\Phi$ from a given propagation model, namely, $\Phi(E_{\text{kin}})\to N_{e^{+}}\,E_{\text{kin}}^{\delta_{e^{+}}}\,\Phi(E_{\text{kin}})$. The two parameters are determined by fitting to the experimental data. In the GH propagation model, a $\chi^{2}$ fitting to the AMS-02 positron flux data in the low-energy region $0.5-5.2$ GeV (16 data points in total) gives $N_{e^{+}}=0.255\pm 0.003$, and $\delta_{e^{+}}=0.995\pm 0.010$ with $\chi^{2}/\text{d.o.f}=32.51/14$. For other propagation models MIN, MED and MAX similar results are found.

In the presence of backgrounds, we derive the corresponding constraint on $f_{\text{PBH}}$ using the standard minimal chi-squred method (or, the log likelihood ratio method). For a given mass function parameters of $M_{\text{c}}$ or a pair of $(M_{\text{c}}\,,\sigma)$, we first find the best-fit parameters $(f_{\text{PBH}},N_{e}^{+},\delta_{e}^{+})$ and the corresponding minimal value of $\chi^{2}_{\text{min}}$. Then, we increase the value of $f_{\text{PBH}}$ and redo the fitting with respect to the other nuisance parameters $N_{e}^{+},\delta_{e}^{+}$ and obtain another value of $\chi^{2}$. This process is repeated until a particular $f_{\text{PBH}}$ satisfying $\Delta\chi^{2}\equiv\chi^{2}(f_{\text{PBH}},\hat{N}_{e^{+}}\,,\hat{\delta}_{e^{+}})-\chi^{2}_{\text{min}}=3.84$ is found. Here, the parameters with a single hat minimized the $\chi^{2}$ function for the given $f_{\text{PBH}}$. The obtained value of $f_{\text{PBH}}$ corresponds to the value excluded at the 95% confidence level. In the left panel of Fig. 5, we show the obtained constraints on $f_{\text{PBH}}$ as a function of $M_{\text{c}}$ in the case of monochromatic mass function in four propagation models. The difference between different propagation models are typically within two orders of magnitude. Among these models, GH model give the most stringent limits. In the right panel of Fig. 5, the constraints for the log-normal mass function with three width $\sigma=0.5\,,1.0\,,2.0$ in the GH propagation model are shown. Similar to the right panel of Fig. 4, increasing $\sigma$ leads to more stringent constraints. In all the case, we have verified that the constraints from the low-energy AMS-02 positron data are always the most stringent compared with that from the AMS-02 all-electron or the Voyager-1,2 all-electron data.

In Fig. 6, we compare the constraints from the AMS-02 positron data obtained in the GH model with a selection of constraints from other observables such as that from the extragalactic $\gamma$-ray Carr:2009jm , CMB Auffinger:2022khh and $511\,\mathrm{keV}$ $\gamma$-ray line Dasgupta:2019cae . We also compare our results with that previously obtained from the Voyager-1 all-electron data in a diffusive-reacceleration propagation model (model A) in Ref. Boudaud:2018hqb . The figure shows that the constraints from AMS-02 positron data are competitive with all these previously obtained constraints, and more stringent than that from the Voyager all-electron data in both the cases with and without including the astrophysical backgrounds.

In summary, in this work we have explored the possibility of using the low-energy CR positron data to constrain the abundance of primordial black holes. The advantage of using the CR positron data is based on its secondary origin and thus sensitive to any exotic contributions. We have shown that in some typical diffusive re-acceleration CR propagation models, the current AMS-02 data can place stringent constraints on $f_{\text{PBH}}$. As an example, we have shown that in the Galprop+Helmod model for CR propagation in the Galaxy and heliosphere, a conservative upper limit of $f_{\text{PBH}}\lesssim 3\times 10^{-4}$ at $M_{\mathrm{PBH}}\simeq 2\times 10^{16}$ g can be obtained, which improves the previous constraints from the Voyager-1 data of all-electrons by around an order of magnitude. Compared with other observables such as diffusive $\gamma$-rays, using the CR flux data to constrain PBH abundance in general suffers from uncertainties in the propagation models, which are expected to be improved in the near future with more precise data from the experiments such as AMS-02, DAMPE and HERD.