# Dark Matter Capture in Celestial Objects:
Treatment Across Kinematic and Interaction Regimes

Rebecca K. Leane  Juri Smirnov

###### Abstract

Signatures of dark matter in celestial objects have become of increasing
interest due to their powerful detection prospects. To test any of these
signatures, the fundamental quantity needed is the rate in which dark matter
is captured by celestial objects. Depending on whether dark matter is light,
heavy, or comparable in mass to the celestial-body scattering targets, there
are different considerations when calculating the capture rate. Furthermore,
if dark matter has strong or weak interactions, the physical behaviour
important for capture varies. Using both analytic approximations and
simulations, we demonstrate how to treat dark matter capture in a range of
celestial objects for arbitrary dark matter mass and interaction strength. We
release our calculation framework as a public package available in both Python
and Mathematica versions, called Asteria .

##  1 Introduction

Capture of dark matter (DM) in celestial objects gives rise to wide-ranging
and exciting opportunities for DM discovery. Assuming interactions with the
Standard Model (SM) particles, DM may scatter with the celestial-body matter,
lose energy, and become gravitationally captured. This can lead to the build
up of a DM population inside the celestial object, with readily detectable
signals. For annihilating DM, these include heating signals, neutrino signals,
gamma-ray signals, or electron signals. For non-annihilating DM, very large
amounts of DM can accumulate over time due to a lack of depletion through
annihilation, leading to interesting observables such as eventual destruction
via black holes from over-accumulation, or changes to stellar evolution. The
range of objects considered in the past includes the Earth and the Sun [1, 2,
3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
Jupiter [1, 24, 25, 26, 27], Brown Dwarfs [28, 29], Uranus [30], Exoplanets
[28], White Dwarfs and Neutron Stars [31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59,
60, 61, 62, 63, 64, 65, 66, 67, 68, 29, 69], and other stars [70, 71, 72, 73,
74, 75].

To test any of these signatures, the fundamental quantity needed is the DM
capture rate. Formalisms to calculate the capture rate of DM go back a very
long time [76, 77], and have recently been improved predominantly in the
context of heavy DM in high-escape velocity objects, such as neutron stars or
white dwarfs [78, 62, 66, 79, 80]. However, in recent years, interest has
grown in signatures of objects with low escape velocities, such as the Earth
or Jupiter, due to highly detectable signals [24, 25, 26, 27, 81, 82, 83, 84,
85, 86, 87, 88, 89, 90, 91]. Interest has also recently grown in the prospects
for light DM discovery, due to strengthening direct detection bounds for GeV-
scale DM [92, 93]. This is especially pertinent as it was recently shown that
even sub-MeV DM masses can be retained in a wide range of celestial objects,
depending on the DM model [94]. Capture of light DM with strong interactions
in the Earth was analytically approximated in Ref. [81]. However, any
framework considering comparable DM and SM target masses, or light DM with
cross sections near the strong-weak interaction boundary, have not been
considered. There are also a number of other regimes which have not yet been
explored in full generality, as well as some inconsistencies in the existing
literature.

The goal of this work is to establish capture rate calculations which are
valid for a wide range of celestial objects, in arbitrary DM mass regimes, and
arbitrary DM-SM scattering cross sections, covering all the regions in our
Fig. 1, and focusing on DM-nucleon interactions. To do this, we will present
both analytic and simulated results, and present some conditions on how to
treat DM capture. We find our improved treatment can lead to results that are
very different to other frameworks. Our framework is released in a public
package available in both Python and Mathematica versions, called Asteria
[95], named after the Greek goddess of stars and planets.

Figure 1: Schematic of the DM capture regimes of interest, as a function of
DM-SM scattering cross section $\sigma_{{\chi}{\rm SM}}$ and DM mass
$m_{\chi}$. Each regime can have unique features requiring different capture
treatments. The transition cross section $\sigma_{\rm tr}$ is the cross
section corresponding to a mean free path of the size of the celestial object,
and therefore the transition between the single and multi-scatter interaction
regime. The DM mass relative to the celestial body SM target mass also affects
the capture kinematics and features.

Our paper is organized as follows. In Section 2, we give an overview as well
as some improvements in the foundations of the main existing multiscatter
formalism, and detail some treatments for the capture summation in Section 3.
In Section 4 we briefly discuss the weak interaction regime. In Section 5, we
discuss the diffusive DM regime, corresponding to DM with light or comparable
to SM target masses, with large interaction cross sections. This includes new
simulation results. In Section 6, we discuss the intermediate “transition”
interaction regime, where the DM-SM interaction rate switches between the
single and multi-scatter treatments. In Section 7, we discuss treatment of the
heavy ballistic DM regime, and in Section 8 we apply our capture formalism to
the Earth, the Sun, Jupiter, and a Brown Dwarf, and describe the physical
features observed throughout each capture regime. We summarize our findings in
Section 9.

##  2 Overview of Capture Formalism

###  2.1 Previous Multiscatter Formalism Assumptions

DM capture can occur via single or multiple scatters in the celestial body,
depending on the kinematic regime. A treatment covering mulitscatter of DM was
first studied in Ref. [77], and a currently widely used framework is from Ref.
[96], though also see Ref. [97]. Further improvements were also shown for
ultra-dense objects in Refs. [57, 62, 66, 80]. However, as the framework of
Ref. [96] was created for heavy DM in neutron stars, it is not generically
applicable to all celestial objects in all DM mass regimes. In particular, the
assumptions of Ref. [96] are that:

  * •

The celestial-body is a neutron star, such that its escape velocity $v_{\rm
esc}$ is much greater than the mean incoming halo DM velocity $v_{\chi}$, i.e.
$v_{\rm esc}\gg v_{\chi}$,

  * •

The DM is heavy compared to the SM targets in the celestial body. This is
because the assumption is that DM travels in a straight line through the
object, corresponding to ballistic forwards motion, which does not apply to
light DM.

For strongly interacting light DM in the Earth, an analytic modification was
proposed in Ref. [81]. This can be used for DM much lighter than the SM target
matter, but does not apply in the weak interaction transition nor weak
interaction regime. No framework exists for comparable DM and SM masses
outside the weak-interaction regime. We will address these cases shortly, but
first we now discuss the foundation for the capture calculation, which will
serve as a backbone for a few of our improved treatments.

###  2.2 Foundation for the Capture Calculation

As we aim to detail capture for a range of objects, DM masses, and DM-SM
interaction regimes, we start by giving an overview of the minimal ingredients
for capture, under the assumption of forward linear DM motion; we will address
the modifications required for specific regimes in the upcoming sections.

The probability of capture after $N$ scatters, assuming isotropic scattering
(i.e. that the differential cross section is independent of the scattering
angle), is given by

$\displaystyle g_{N}(u)$
$\displaystyle=\int_{0}^{1}dz_{1}\,...\int_{0}^{1}dz_{N}\,\Theta\left(1-\sqrt{1+w^{2}}\,\prod_{i=1}^{i=N}\sqrt{1-z_{i}\beta}\right)\,,$
(2.1)

with $z_{i}$ being the angular scattering variables, $w=u/v_{\rm esc}$, where
$u$ is the velocity of the incoming DM particle, and

$\beta=\frac{4m_{\rm SM}\,m_{\chi}}{(m_{\rm SM}+m_{\chi})^{2}},$ (2.2)

where $m_{\rm SM}$ is the SM target mass, and $m_{\chi}$ is the DM mass. We
define $A=\prod_{i=1}^{i=N-1}\sqrt{1-z_{i}\beta}$, and use it rewrite the
$N$th integral. This integral can be evaluated analytically, however, to do
this the required condition is that

$\displaystyle A\,\sqrt{u^{2}+v_{\rm esc}^{2}}(1-\beta)^{1/2}<v_{\rm esc}\,,$
(2.3) $\displaystyle\text{and }A\,\sqrt{u^{2}+v_{\rm esc}^{2}}>v_{\rm esc}.$

The condition in the first line of Eq. (2.3) implies that the maximal energy
transfer possible in the $N$th scatter with $z=1$ brings the DM velocity below
the escape velocity such that DM is captured. The condition in the second line
implies that the DM particle is still above the escape velocity after the
previous $(N-1)$th scatter. Thus, as pointed out in Ref. [89], the evaluated
integral gives the probability for getting captured after exactly $N$
collisions, and reads

$\displaystyle\int_{0}^{1}dz_{N}\,\theta\left(1-\sqrt{1+w^{2}}\,A\,\sqrt{1-z_{N}\beta}\right)=1-\frac{1}{\beta}+\frac{1}{\beta\,(1+w^{2})A^{2}}\,.$
(2.4)

In Ref. [59], this expression has been evaluated, and after further $N-1$
integrations is

$\displaystyle\left.g_{N}(u)\right|_{\rm Exactly\ N}$
$\displaystyle=\int_{0}^{1}dz_{1}\,...\int_{0}^{1}dz_{N-1}\left(1-\frac{1}{\beta}+\frac{1}{\left(1+w^{2}\right)\prod_{i=1}^{i=N-1}\left(1-z_{i}\beta\right)\beta}\right)$
(2.5)
$\displaystyle=1-\frac{1}{\beta}+\frac{1}{\beta^{N}\left(1+w^{2}\right)}\log{\left[\frac{1}{1-\beta}\right]}^{N-1}\,.$
(2.6)

The capture rate after $N$ or less scatterings is then obtained by summation
over $N$ and integration over the DM velocity distribution [89]

$\displaystyle C_{N}=\pi
R^{2}p_{N}(\tau)\,\sum_{i=1}^{N}\int_{0}^{\infty}f(u)u\left(1+w^{-2}\right)\,\left.g_{i}(u)\right|_{(\rm
Exactly\ i)}\,du\,.$ (2.7)

To get the total capture rate one needs to sum the coefficients in Eq. (2.7)
up to a maximal number of scatterings $N_{\rm max}$,

$\displaystyle C_{\rm total}=\sum_{i=1}^{N_{\rm max}}C_{i}\,.$ (2.8)

$N_{\rm max}$ is a truncation condition used to allow numerical evaluation,
which we discuss in the next section. Note that while Eq. (2.6) was evaluated
in Ref. [59], that reference misses the summation in Eq. (2.7) and therefore
produces unphysical behaviour; see Appendix B of Ref. [89] for details.

The expression in Eq. (2.7) is greatly simplified in the large interaction
regime, where the DM undergoes many scatterings. The nested integral in Eq.
(2.1) can be evaluated using the assumption that the scattering variable is on
average $\langle z_{i}\rangle\approx 1/2$, which yields an approximate result
for $g_{N}(u)\approx g_{N}(u)^{\rm avg}$, as discussed in Ref. [96]. The
expression $g_{N}(u)^{\rm avg}$ corresponds to the probability for the DM
particle to drop below escape velocity at any of the $N$ scatterings, as no
additional kinematic conditions are needed to evaluate the integral.
Therefore, no additional summation is required and integrating $g_{N}(u)^{\rm
avg}$ over the DM velocity distribution

$\displaystyle C_{N}\approx\pi
R^{2}p_{N}(\tau)\int_{0}^{\infty}f(u)u\left(1+w^{-2}\right)g_{N}(u)^{\rm
avg}\,du\,,$ (2.9)

we obtain a compact analytic expression for the capture coefficient

$C_{N}=C_{\rm
geo}\,p_{N}(\tau)\,\left(1-\frac{\exp\left[-\frac{3}{2}\frac{v_{\rm
esc}^{2}}{v_{\chi}^{2}}\left(\alpha_{\mu}^{-N}-1\right)\right]\left(1+\frac{3}{2}\frac{v_{\rm
esc}^{2}}{v_{\chi}^{2}}\alpha_{\mu}^{-N}\right)}{1+\frac{3}{2}\frac{v_{\rm
esc}^{2}}{v_{\chi}^{2}}}\right)\,,$ (2.10)

where the geometric capture rate $C_{\rm geo}$ is the total rate of DM
particles passing through the celestial body

$C_{\rm geo}=\frac{\pi
R^{2}\rho_{\chi}v_{\chi}}{m_{\chi}}\sqrt{\frac{8}{3\pi}}\left(1+\frac{3}{2}\frac{v_{\rm
esc}^{2}}{v_{\chi}^{2}}\right)\,,$ (2.11)

where $\rho_{\chi}$ is the DM density in the surrounding environment, as found
in Ref [96]. The energy fraction remaining after each scatter is given by

$\alpha_{\mu}=1-\frac{4\langle z\rangle\mu}{(1+\mu)^{2}},$ (2.12)

with $\mu=m_{\chi}/m_{\rm SM}$.

We use the full expression for arbitrary scattering angle in the weak
interaction limit, as the averaged scattering angle is only a good
approximation for multi-scatter interactions. In the case of large
interactions and multi-scatter, the $\langle z_{i}\rangle\approx 1/2$
approximation speeds up and the calculation without a relevant loss of
accuracy (see Fig. 5 in App. B of Ref. [89]), and so we use Eq. (2.10) for the
capture coefficient, to get an expression for the capture rate after $N$ or
less scatters, which is then summed over per Eq. (2.8).

In the above form of Eq. (2.10) it is easy to see that the capture rate for
all $C_{N}$ drops to zero if $\mu\rightarrow 0$ or $\mu\rightarrow\infty$, as
$\alpha_{\mu}\rightarrow 1$, i.e. the remaining energy fraction after each
scattering approaches one. This shows that the energy loss efficiency is poor
both for increasingly light or heavy DM mass compared to the target mass, such
that more scatters are increasingly required in these limits. Note that we
neglect the relative motion of the celestial object to the DM halo, which is a
small correction. We also neglect thermal motion of the targets, which is
justified for interactions with nuclei in most celestial objects given
sufficient DM mass.

The probability of a single DM particle undergoing $N$ scatters is [96]

$p_{N}(\tau)=2\int_{0}^{1}dy\frac{ye^{-y\tau}(y\tau)^{N}}{N!},$ (2.13)

where $y$ is an angular impact variable and $\tau$ is the optical depth,

$\displaystyle\tau=\sum_{i}\,\frac{3}{2}\frac{\sigma_{\chi_{A_{i}}}}{\sigma_{\text{tr}_{A_{i}}}}=\sum_{i}\tau_{A_{i}},$
(2.14)

where $\sigma_{\text{tr}_{A}}=\pi R^{2}/N_{A}$ is the transition cross section
for DM capture, with $N_{A}=f_{A}M_{\rm obj}/m_{A}$ being the number of atoms
with atomic mass number $A$ in the object, with $f_{A}$ being the mass
fraction of the element. The transition cross section is when the mean free
path of the DM is the diameter of the celestial object, and therefore is at
the transition of the single and multiscatter regimes. Note that $\sigma_{\rm
tr}$ has been called the “saturation cross section” $\sigma_{\rm sat}$ in
previous works, and we advocate for renaming it as calling it the saturation
cross section is not physically accurate for all objects; we discuss this in
more detail in the upcoming Sec. 6.1.

For spin-independent scattering rates in celestial objects, we convert the
nucleon to nucleus cross sections using the Born approximation,

$\displaystyle\sigma_{\chi\rm
A}=A^{2}\left(\frac{\mu_{R}(A)}{\mu_{R}(N)}\right)^{2}{\sigma_{\chi N}},$
(2.15)

where $A$ is the atomic mass number of the relevant nucleus, and
$\mu_{R}(A\text{ or }N)$ are the reduced masses of the DM with the nucleus or
nucleon. Note that the Born approximation breaks down for large cross
sections, and in this regime particle DM models should be used directly [98,
99]. For practical purposes an approximation procedure allows preservation of
unitarity of the interaction by assuming $\sigma_{\chi A}^{\rm
tot}=\min\left(\sigma_{\chi A},4\pi r_{N}^{2}\right)$, where
$r_{N}\approx(1.2\,\text{fm})\,A^{1/3}$ [98]. However, if DM is an extended
object with $r_{\chi}>r_{N}$, larger cross sections are possible. To remain
model agnostic we do not limit the cross section size, but note that depending
on the model, the treatment of large cross sections ($\sigma_{\chi N}\gtrsim
10^{-31}$ [98]) can substantially deviate from our setup here. In our package
Asteria [95] there is the option to include the $A$ scaling, however we
emphasize that this will not be valid for large cross sections. We therefore
also provide the possibility for the input of the total $\sigma_{\chi A}$
cross section, such that no scaling assumption with atomic mass number is
enforced internally.

###  2.3 Multi-Element Treatment

Treating capture for multi-component celestial objects can quickly become
cumbersome, due to the large number of scattering combinatorics. We make a few
simplifying assumptions, without sacrificing much accuracy. Firstly, we assume
a homogeneous distribution of the chemical elements throughout the celestial
body. We discuss our additional assumptions below depending on the regime.

####  2.3.1 Multi-Scatter Case

As discussed in Ref. [100] in the context of a two component model in the
optically thick regime, the scattering combinatorics can start to play a role.
We find that a sufficient treatment is to switch at cross sections above the
transition cross section i.e. $\tau>3/2$, to the effective treatment using an
average mass, dominated by the element that is most likely to be encountered
by DM

$m_{\rm SM}=\frac{1}{\tau}\sum_{i}\tau_{A_{i}}\,m_{i}\,.$ (2.16)

This treatment greatly accelerates the computation compared to evaluating all
the combinatorics, and produces results in an agreement within a few percent
with the full approach in Ref. [100].

####  2.3.2 Single-Scatter Case

In the case that $\tau<3/2$ the scattering is dominated by the first, or first
few, scattering events. In this regime it is sufficient to calculate the
capture rates per element $C_{A_{i}}$, which are dependent on the
$\tau_{A_{i}}$ for each element independently, and sum them to obtain a total
capture rate

$\displaystyle C_{\rm all\,elements}=\sum_{i}C_{A_{i}}\,.$ (2.17)

This treatment agrees within percent level accuracy with the full calculation
that takes into account combinatorics, see Ref. [100].

##  3 How to Treat the Capture Summation

The multi-scatter sum has formally an infinite number of elements, thus in
order to evaluate it numerically we use a truncation condition, and several
simplifications to accelerate the calculation in certain regimes. We will
first formulate and justify the truncation condition, defining the maximum
number of scatters, and discuss regime specific simplifications in the
following dedicated sections.

In order to formulate the truncation condition it is crucial to understand the
behaviour of the probability of DM to scatter $N$ times in the object. While
the expression for the capture rate at large cross section in Eq. (2.8) is
correct, it can be non-trivial to numerically compute, particularly due to the
$p_{N}(\tau)$ as written in Eq. (2.13), which can require evaluating large
$N!$. We now discuss the behavior of the scattering probability at large
optical depth, that will allow us to formulate a truncation condition for the
capture rate sum.

To improve numerical stability in the strong interaction regime, it is
advantageous to write Eq. (2.13) as [97]

$p_{N}(\tau)=\frac{2(N+1)}{\tau^{2}}\left(1-\frac{\Gamma[N+2,\tau]}{\Gamma[N+2]}\right),$
(3.1)

where $\Gamma[a]$ is Gamma function, and $\Gamma[a,b]$ is the upper incomplete
Gamma function; gamma functions are used as they are known as the generalized
factorial, as $\Gamma[N+1]=N!$ for natural numbers $N$.

Figure 2: Probability $p_{N}$ of exactly $N$ scatters for different choices of
the optical depth $\tau$.

Figure 2 shows the scatter probability as a function of the number of scatters
for different choices of optical depth $\tau$, 3000 and 5000. We see that the
probability to scatter $N$ times peaks around $N\sim\tau$, which is true for
all $\tau\gg 3/2$. This means that the probability to scatter many more times
than the optical depth is negligible, in the strongly interacting regime. We
make this point more explicit by considering an expansion of Eq. (3.1) in two
regimes

$\displaystyle
p_{N}(\tau)\approx\frac{2(N+1)}{\tau^{2}}\times\begin{cases}\left(1-\dfrac{e^{-\tau}\tau^{1+N}}{\Gamma[N+2]}\right)&(N<\tau),\\\
\\\
\left(\dfrac{e^{N-\tau}N^{-5/2-N}\tau^{N+2}}{\sqrt{2\pi}}\right)&(N>\tau)\,.\end{cases}$
(3.2)

In the limit that $\tau\gg 3/2$, the right bracket for the $N<\tau$ case
approaches one, while the bracket for $N>\tau$ approaches zero, i.e.

$p_{N}(\tau)\approx\frac{2(N+1)}{\tau^{2}}\Theta(\tau-N),$ (3.3)

in agreement with Ref. [97]. This shows that for large cross sections only
$N<\tau$ meaningfully contributes to the total capture rate. If instead we are
in the weakly interacting regime, i.e. if $\tau$ is smaller than about $3/2$,
the probability function $p_{N}(\tau)$ peaks at $N=1$ (meaning one scatter)
and rapidly decays with growing $N$.

In the package that we provide along with this paper, the numerical accuracy
goal for $p_{N}$ is at the percent level, and so we find that it is
computationally sufficient to truncate the summation in Eq. (2.8) at

$N_{\rm max}=\max\\{10,\lfloor e\tau\rfloor\\}\,.$ (3.4)

This is an arbitrary choice, which is justified because at $N\sim e\tau$, Eq.
(3.2) for $\tau<N$ leads to a scaling for $p_{N}(\tau)\sim
e^{-\tau}/\tau^{5/2}$, which is negligible for large $\tau$ values. For
smaller optical depth, the capture is controlled by the first few scatters,
since the at small $\tau<3/2$ Eq. (3.1) leads to a rapidly falling scattering
probability $p_{N}(\tau)\sim\tau^{N}/N!$. We make thus a conservative choice
if we sum to $N_{\rm max}=10$ in this case.

There is one more caveat concerning the optical depth, which is that there are
a finite number of SM targets to scatter with. This is generally only relevant
in the ultra-heavy DM limit; we discuss it further in Sec. 7.

By implementing the above conditions we find that the numerical evaluation of
Eq. (2.8) is greatly accelerated, and is within a sub-percent level accuracy
of the fully converged general result.

We now discuss regime-specific calculations for DM capture rates.

##  4 Weak Interaction Regime

Figure 3: The weak regime discussed in this section.

The optically thin regime, also known as the weak interaction regime, has been
the most commonly studied. The original results of Gould [77] (see also
Spergel and Press [76]) hold in this regime, but are also approximately
reproduced in the limiting case of weak cross sections in the multiscatter
calculation in Eq. (2.8). Here, there is no need to deviate from the simple
non-linear forward scattering picture described in Sec. 2, as the DM may
change direction, but it only does so roughly once (as far as capture is
concerned), making the direction it turns largely irrelevant.

###  4.1 Numerical Simplifications

When taking the single scatter limit of Ref. [96], one point that potentially
requires some care is the evaluation of the scattering probability. It can be
numerically advantageous to use the leading order expansion of Eq. (3.2) for
$\tau\ll 1$,

$\displaystyle p_{N}(\tau)\approx\frac{2\,\tau^{N}(N+1)}{(N+2)\Gamma[N+2]}\,.$
(4.1)

In addition, when the mass mismatch between the DM and the SM target is very
large i.e. $\mu\ll 1$ or $\mu\gg 1$ the capture rate coefficients are
$C_{N}/G_{\rm geo}\ll 1$ and can become numerically unstable in the weak
regime. In this case an approximate expression can be used to alleviate the
problem,

$\displaystyle C_{N}\approx 9\,C_{\rm geo}\,p_{N}(\tau)N\,\frac{v_{\rm
esc}^{4}}{v_{\chi}^{4}}\frac{\mu}{(1+\mu)^{2}}\,\left(2+3\frac{v_{\rm
esc}^{2}}{v_{\chi}^{2}}\right)^{-1}\,.$ (4.2)

This expression is valid for $\mu\ll 1$ and $\mu\gg 1$, which is the case for
inefficient energy transfer, given a large mass mismatch. When in the small
optical depth regime, we replace $C_{N}$ by the approximate expression in Eq.
(4.2) when $\mu/(1+\mu^{2})\times v_{\rm esc}^{2}/v_{\chi}^{2}<10^{-5}$ which
improves the numerical stability and allows us to evaluate even highly
suppressed capture rates.

##  5 Diffusive Dark Matter Regime

Figure 4: The strong light and strong comparable regime discussed in this
section.

As discussed above, the capture ingredients in Sec. 2 assume that all DM
trajectories are a straight forward line through the object. This holds for DM
much heavier than the SM target, or weakly interacting DM that is much lighter
than the SM target.

On the other hand, when DM has strong interactions and is light or comparable
in mass to the SM target, it has diffusive motion while being captured, and
can be more easily have its trajectory altered backwards or sideways [101,
102, 81]. In the multi-scatter regime (i.e., the strong interaction regime),
these non-forwards trajectories need to be tracked. To include this effect in
Eq. (2.10), we will modify the expression with an additional factor $f_{\rm
cap}$, which takes into account reflection of light DM out of the celestial
object. This is required as Eq. (2.10) only applies for heavy DM in celestial
objects with high escape velocities, and we want to extend the setup to
include light DM masses in objects with low escape velocities. In addition, as
we will show, this diffusive motion factor can also actually be important for
high-escape velocity objects such as the Sun when DM becomes very light, as
the DM is sped up to the escape velocity and can still end up reflected back
out of the object.

###  5.1 Required Number of Scatters for Maximal Dark Matter Capture

To calculate the amount of reflection, we will need an estimate of the
expected number of scatters before capture. The number of scatters at which
the majority of DM particles are moving at an average velocity was derived in
Ref. [103], however, the only the halo velocity was considered for the
incoming DM velocity. We extend the formula to take into account the effect of
the gravitational infall i.e. $v_{\rm in}=\sqrt{v_{\chi}^{2}+v_{\rm
esc}^{2}}$, which is important for heavy objects. We therefore find the number
of scatterings that is required for the bulk of the DM particles to be
captured,

$\displaystyle N_{\rm req}(\mu,v_{\chi},v_{\rm esc})=\frac{\log\left(v_{\rm
esc}^{2}/(v_{\chi}^{2}+v_{\rm
esc}^{2})\right)}{\log\left(\alpha_{\mu}\right)}\,,$ (5.1)

which depends on the escape velocity and average DM velocity, and the DM mass
and the average target mass, and where again the average of a scattering
variable $\langle z\rangle\sim 1/2$ has been taken assuming isotropic
scattering. From Eq. (5.1), it is clear that if the DM mass is significantly
larger or smaller than the SM target mass, multiple scatterings have to be
taken into account in order to correctly compute the DM capture rate.

###  5.2 Reflection Coefficient in the Low Dark Matter Mass Limit

The simplest case to consider including reflection is the limit that the DM is
much lighter than the SM target. In this case, the reflection factor is given
by

$\displaystyle f_{\rm cap}^{\rm light}\approx\frac{2}{\sqrt{\pi\,N_{\rm
req}}}=\left[\frac{4}{\pi}\frac{\log\left(\alpha_{\mu}\right)}{\log\left(v_{\rm
esc}^{2}/(v_{\chi}^{2}+v_{\rm esc}^{2})\right)}\right]^{1/2},$ (5.2)

where $N_{\rm req}$ is the total number of scatters required for capture
(distinct from $N$, which is instead the total number of scatters the DM
actually undertakes). A similar expression was pointed out in Ref. [81] for
the Earth, though there the DM speed up from the escape velocity was not
included, which makes a non-negligible difference for higher escape velocity
objects. In any case, as the DM mass approaches the SM target mass, this
approximation becomes less reliable, and so far there was no analytic
expression to be used. We now investigate more general behavior of the
reflection coefficient $f_{\rm cap}$.

###  5.3 New Simulations and Analytic Fits

####  5.3.1 Simulation Setup

We perform simple simulations to investigate the generic behaviour of the
reflection coefficient. Our simulation setup is similar to that performed in
Ref. [104] for the DM distribution within the Sun, although we will apply our
simulation to a different problem which is the fraction of reflected DM
particles, rather than the internal DM distribution.

We perform a 3D Monte Carlo Simulation assuming the DM particles enter the
celestial object with a velocity of $v_{\rm in}=\sqrt{v_{\rm
esc}^{2}+v_{\chi}^{2}}$, and at an impact vector randomly drawn from a
distribution uniform in angles. In order to determine the next interaction
point, including the first point after entry, we draw the realized free path
from an exponential distribution proportional to $e^{-\tau}$, where $\tau$ is
the optical depth. We perform the DM particle scattering in the center of mass
frame, which conserves the magnitude of the velocity vector. We assume
isotropic scattering, thus drawing the new direction angles from a flat
distribution, and boost to the lab frame of the object. The simulation is run
until the DM particle either leaves the object, or its velocity drops below
$v_{\rm esc}$, in which case it is considered captured. We sample over
$10^{4}$ particles per mass point, as this number of particles exhibits well
converged behavior. This gives us our data for our analytic fit to the
reflection coefficient. The simulations are valid in the optically thick
regime, as the mean free path is limited to near the surface of the object. We
therefore simplify the setup by only taking the $v_{\rm esc}$ value at the
surface.

Figure 5: Simulations for DM in Jupiter in the multi-scatter regime;
qualitatively similar results are obtained for the Earth. The color of the
lines indicates the DM velocity at the end of the run, as a fraction of the
incoming DM velocity. At the end of the simulation DM particles have either
dropped below the escape velocity to become captured, or have escaped from the
object. Left: The diffusive regime where the DM mass is half that of the SM
targets. Right: Beginning the departure from the diffuse regime, with DM five
times heavier than the SM target mass.

####  5.3.2 Simulation Results

Figure 5 shows two examples of our simulated outcomes, for DM in Jupiter. In
the left panel, the DM is lighter than the SM target, and we clearly see
diffusive behavior. That is, upon entering the celestial object, the DM
wanders around and does not exhibit dominantly forward motion. In the right
plot, we simulated DM that is slightly heavier than the DM target. In this
case, the DM more often travels forwards in a straight line, though as the
masses are not too different, it is still sometimes reflected, albeit not as
often as the lighter DM scenario. As the escape velocity of Jupiter
corresponds to about 22 percent of the average incoming DM halo velocity, any
velocities below the 0.22 remaining DM velocity fraction (dark blue/violet
trajectories) are captured in both panels.

To make our simulations easily applicable to arbitrary celestial objects, we
fit analytic functions to the simulation results. The functional form of the
fit functions is chosen in order to describe the simulation data with a
minimal number of parameters, and the coefficients are determined in a common
best-fit minimization across six different celestial objects, ranging in
escape velocities from that of the Earth to the Sun. We find therefore that
the reflection coefficient can be generalized to

$\displaystyle f_{\rm cap}(\mu)=\begin{cases}f_{\rm cap}^{\rm
light}(\mu)&(\mu<\mu_{T})\\\ \\\ \dfrac{\mu_{M}f_{\rm cap}^{\rm
light}(\mu_{T})-\mu_{T}f_{M}+\mu\left(f_{M}-f_{\rm cap}^{\rm
light}(\mu_{T})\right)}{\mu_{M}-\mu_{T}}&(\mu_{T}<\mu<\mu_{M})\\\ \\\
\dfrac{\mu}{(\mu-\mu_{M})+\mu_{M}/f_{M}}&(\mu>\mu_{M})\end{cases}$ (5.3)

for the DM mass regimes as indicated. Here the DM-SM mass ratio $\mu_{T}$ is
determined by $N_{\rm req}(\mu_{T})=N_{T}$, where the critical number of
scatterings needed for capture $N_{T}$ is determined by a fit to the
simulations and is given by

$N_{T}\approx 12+1.8\log\left(\sqrt{1+\frac{v_{\chi}^{2}}{v_{\rm
esc}^{2}}}\right).$ (5.4)

This implies a value for $\mu_{T}$, given by

$\mu_{T}=\frac{\sqrt{2\left(v_{\rm esc}^{2}/(v_{\chi}^{2}+v_{\rm
esc}^{2})\right)^{\frac{1}{N_{T}}}-1}-\left(v_{\rm
esc}^{2}/(v_{\chi}^{2}+v_{\rm esc}^{2})\right)^{\frac{1}{N_{T}}}}{\left(v_{\rm
esc}^{2}/(v_{\chi}^{2}+v_{\rm esc}^{2})\right)^{\frac{1}{N_{T}}}-1}$ (5.5)

The mass ratio at the second regime change is found to be

$\mu_{M}=1.56\,\left[1-\frac{1}{1+0.52\log\left(\sqrt{1+v_{\chi}^{2}/v_{\rm
esc}^{2}}\right)}\right],$ (5.6)

and the capture fraction at the second regime change is given by

$f_{M}=0.22\,\left[1+\frac{3.58}{1+0.23\log\left(\sqrt{1+v_{\chi}^{2}/v_{\rm
esc}^{2}}\right)}\right].$ (5.7)

Our treatment is applicable across a wide range of objects with different
escape velocities, up to $v_{\rm esc}\sim 1000\text{ km/s}$. Objects with
larger escape velocities, such as white dwarfs and neutron stars, can be
approximated with the formalism we present here, however their high densities
can require a more detailed treatment of nuclear effects to obtain more
accurate results, and so are beyond the scope of our treatment.

Figure 6: Captured DM fraction taking into account reflection for the Earth, a
Super Earth, Jupiter, a Super Jupiter, a Brown Dwarf, and the Sun, with escape
velocities as labeled. Our simulations are shown as dots labeled “Simulation”,
and our analytic model fits to our simulations are shown as solid lines
labeled “Best-fit”. The previous approximation for the low mass limit for the
Earth is shown as a gray dashed line, labeled “Low Mass Approximation”, which
we match well in the limit of its validity.

Note that the reflection coefficients we approximate here are only relevant
for $\tau\gtrsim N_{\rm req}$. At smaller optical depth, the capture is
dominated by the first few scatters, and reflection is not relevant. The
reflection coefficient is also only relevant in the diffusive regime, which
corresponds to the light DM or comparable to target mass regime. However, in
the heavy mass regime where $m_{\chi}\gg m_{\rm SM}$, $f_{\rm cap}\sim 1$ only
provided that the DM mass does not become so large that it cannot be
kinematically stopped. We discuss the separate modification for strongly-
interacting and ultra-heavy DM in Sec. 7.

Figure 6 shows a comparison of our simulated data to our functional fits in
Eq. (5.3) for the reflection coefficient. These are shown for six benchmark
celestial objects. The only assumption about the celestial object needed for
these plots is the escape velocity of the object as labeled. In the top-left
panel, we also compare our results to the analytic expression in Eq. (5.2) as
pointed out for the Earth in Ref. [81] in the light DM limit, shown as “Low
Mass Approximation”. Comparing to our simulations, we find that the analytic
approximation in Eq. (5.2) is accurate in the limit that the DM is much
lighter than the SM target as expected, and corresponds to a critical number
of scatterings of $\mathcal{O}(10)$ being required. For DM mass closer to the
SM target mass, our analytic fits to our simulations are better approximations
to quantify reflection. We also show that high-escape velocity objects such as
the Sun have DM capture efficiency loss due to reflection, which was not
previously considered.

###  5.4 Connecting Reflection Coefficients to Light Dark Matter Capture Rates

The maximum capture rate in any scenario is given by

$C_{\rm max}=C_{\rm geo}\,f_{\rm cap},$ (5.8)

where $f_{\rm cap}$ is the reflection factor. However, even at very large
scattering cross sections, $C_{\rm max}$ will not necessarily be reached due
to scattering kinematics. Therefore, the total light DM capture rate is given
as a minimum of

$C_{\rm tot}^{\rm light}=\textrm{Min}\left(C^{\rm light}\,,\,C_{\rm
max}\right).$ (5.9)

To evaluate $C^{\rm light}$, we use different numerical tricks depending the
cross section size. For $\tau>3/2$ (the minimum condition for multi-scatter),
and $\tau<100$, we simply use Eq. (2.9) for $C^{\rm light}$. To accelerate the
computation in the case of very large cross sections above $\tau\approx 100$,
we use the fact that the number of scatterings is sufficiently large in order
for the sum in Eq (2.8) to be replaced by an integral over $N$. Given the
simple form of the expression for $p_{N}(\tau)$ in Eq. (3.3), which is
justified at such large optical depth, the integral evaluation is very
efficient. Finally, at even larger optical depth above $\tau\approx 10^{6}$,
the evaluation of the integrand leads to a numerical instability. This can be
avoided, however, by performing a Taylor expansion of the integrand. In the
case of light DM an expansion in $\mu\ll 1$ to order $\mathcal{O}(\mu^{4})$ is
used,

$\displaystyle\frac{C_{N}^{\rm light}}{C_{\rm
geo}p_{N}(\tau)}\approx\min\Biggl{\\{}1,\Biggl{|}\,3\mu
Nv_{\text{esc}}^{4}\Biggl{(}24+\mu^{3}N^{3}\Bigl{(}64-3\Bigl{(}9\frac{v_{\text{esc}}^{4}}{v_{\chi}^{4}}\Bigl{(}\frac{v_{\text{esc}}^{2}}{v_{\chi}^{2}}-6\Bigr{)}+76\frac{v_{\text{esc}}^{2}}{v_{\chi}^{2}}\Bigr{)}\Bigr{)}$
$\displaystyle-4\mu^{2}\Bigl{(}3\mu-1\Bigr{)}N^{2}\Bigl{(}\frac{9v_{\text{esc}}^{4}}{v_{\chi}^{4}}-\frac{30v_{\text{esc}}^{2}}{v_{\chi}^{2}}+16\Bigr{)}-4\mu\Bigl{(}\mu\Bigl{(}5\mu-6\Bigr{)}+3\Bigr{)}N\Bigl{(}\frac{3v_{\text{esc}}^{2}}{v_{\chi}^{2}}-4\Bigr{)}$
$\displaystyle+8\mu\Bigl{(}\mu-3\Bigr{)}\Biggr{)}\times\Biggl{(}8\Bigl{(}3v_{\text{esc}}^{2}v_{\chi}^{2}+2v_{\chi}^{4}\Bigr{)}\Biggr{)}^{-1}\Biggr{|}\Biggr{\\}}.$
(5.10)

This procedure results in a sub-percent accuracy at a substantial increase in
calculation speed.

##  6 Discussion of the Intermediate Interaction Regime

Figure 7: The transition between weak and strong interaction regimes, as
discussed in this section.

We now consider the intermediate interaction regime, which we mark by the
transition cross section $\sigma_{\rm tr}$, which has been previously labeled
$\sigma_{\rm sat}$ in some works in the literature. We first discuss the
meaning of this cross section and advocate for the term of a “transition cross
section” rather than the previous term “saturation cross section”, then
discuss the capture behavior when the cross sections switch from weak to
strong interactions.

###  6.1 The Meaning of the “Saturation” Cross Section

We emphasize that the physical meaning of $\sigma_{\rm tr}$ (previously
sometimes called $\sigma_{\rm sat}$ in the literature), where the optical
depth is $\tau=3/2$, is simply that the mean free path becomes of order the
size of the object, and therefore larger cross sections than $\sigma_{\rm tr}$
enter the multiscatter regime, and smaller cross sections than $\sigma_{\rm
tr}$ correspond to the single scatter regime. $\sigma_{\rm tr}$ should not be
interpreted as the cross section generically resulting in the maximum capture
rate, as we now discuss.

There is some ambiguity in the literature about the meaning of the saturation
cross section, $\sigma_{\rm sat}$ (which we are now calling $\sigma_{\rm
tr}$). This quantity was first referred to as a saturation cross section in
Ref. [78], which studied DM capture in neutron stars. Because neutron stars
have extreme escape velocities, they are very efficient at capturing DM.
Provided the DM mass is not too different compared to the SM target mass this
means that the maximum capture rate is usually achieved with a single scatter,
such that the maximal capture and the “saturation” of the cross section is
reached right before the multi-scatter regime kicks in. This is because
neutron stars’ capture efficiency implies adding any scatters more than
roughly once will not help them capture DM any further. Therefore, calling
$\sigma_{\rm tr}$ a saturation cross section is valid when applied to its
original use in Ref. [78] for neutron stars, although it is not the physical
meaning of the cross section. In fact, for any object where the escape
velocity is higher than the DM halo velocity, $\sigma_{\rm tr}$ can correspond
to the maximum capture rate provided the DM mass is not too far from the SM
target mass.

However, we emphasize that this is not true in objects with lower escape
velocities compared to the DM velocity, such as Jupiter or the Earth. For
objects with escape velocities below the DM halo velocity, because they are
not very efficient at capturing DM, usually many scatters are required to slow
down and capture the DM, and so larger cross sections than $\sigma_{\rm tr}$
are often required to capture the bulk of the DM.

Figure 8 makes this point explicit, and shows the fraction of DM particles
passing through Jupiter that are captured as functions of the DM-SM scattering
cross section. Here Jupiter is chosen as an example of a low-escape velocity
object to illustrate the meaning of $\sigma_{\rm tr}$; qualitatively similar
results are obtained for other low escape velocity objects such as the Earth.
Fig. 8 shows a range of DM masses: the solid lines are for DM heavier than
target, and the dashed lines are the inverse of the solid lines, making them
lighter than the SM target and symmetrically far away from the target mass. As
the DM scattering is most kinematically efficient when the DM is comparable to
the SM mass, taking DM masses that are symmetrically lighter or heavier than
the target requires multiscattering for capture. The symmetry in the behaviour
of DM mass both the heavier and lighter than the target show the identical
scaling of the capture rate, apart from the fact that at very high cross
sections light DM is reflected, which explains the departure of the dashed
curves for the light DM at high cross sections.

In Fig. 8, we see that compared to the transition cross section $\sigma_{\rm
tr}\approx 10^{-34}~{}$cm2 for Jupiter, the maximum capture rate is actually
achieved at cross sections higher than the transition cross section, simply
because Jupiter is not very efficient at capturing DM due to its escape
velocity being lower than the local DM halo velocity, and therefore requires
the multiscatter regime in order to capture the bulk of the DM. This clearly
illustrates that the “saturation cross section” is not the physical meaning of
$\sigma_{\rm tr}={\pi R^{2}}/{N_{n}}$, and taking the “saturation cross
section” as the cross section providing the maximum DM capture rate, would
yield capture rates orders of magnitude too small.

Figure 8: Fraction of DM particles passing through Jupiter that are captured
as functions of the DM-SM scattering cross section, for several DM masses as
labeled. The transition cross section $\sigma_{\rm tr}$ for Jupiter is shown
as a vertical line, and the gray dashed line at $C/C_{\rm geo}=1$ corresponds
to all the DM passing through Jupiter being captured.

Given the broadened interested in DM in a variety of celestial objects, and
the true physical meaning of this cross section, we advocate for calling this
cross section instead the transition cross section, which in this work we are
calling $\sigma_{\rm tr}$. As we described above, $\sigma_{\rm tr}={\pi
R^{2}}/{N_{n}}$ marks the transition between the single and multiscatter
regimes for DM capture in arbitrary celestial objects, it is not generically
the cross section where the maximum capture rate will be found. Indeed, note
that already in some works, e.g. Ref. [80], the term “threshold cross section”
and denoted also by $\sigma_{\rm tr}$ was also already used, instead of
“saturation cross section”.

##  7 Ballistic Dark Matter Regime

Figure 9: The heavy DM mass regime, with strong interactions discussed in this
section.

As discussed above, when DM is heavy, the straight forward-line motion is
accurate for the DM trajectory inside a celestial object. However, if DM
becomes substantially heavier than the SM target mass, it can have another
problem. This is simply that it can have so much momentum that the celestial-
body has insufficient stopping power, and the ultra-heavy DM blasts out the
other end of the celestial body without being captured, even in the strong
interaction regime. This can occur if the optical depth of DM in a given
object does not allow the required number of scatters for capture. In Sec. 3,
we included an improved condition, similar to the one introduced in Ref. [105]
for the Earth, which covers this scenario.

In addition, when the DM mass is much heavier than the target mass, and the
scattering is fully forward dominated, there is a finite number of targets
that the particle can encounter on its way through the object. This number is
given by $N_{c}\sim N_{n}^{1/3}$, with $N_{n}$ being the total target number
in the object. Thus, we consider only physically relevant values of the
optical depth $\tau_{\rm physical}=\min\left[\tau,N_{c}\right]$ in this
regime.

Finally, to calculate DM capture in the case of heavy DM and very large cross
sections, we use similar approximations to calculating $C_{N}^{\rm light}$ in
the previous section. For optical depths above $\tau\approx 100$, we replace
the sum in Eq. (2.8) by the integral, and use the $p_{N}(\tau)$ from Eq.
(3.3). Again at optical depth above $\tau\approx 10^{6}$, the evaluation of
the integrand leads to a numerical instability, and is avoided by performing a
Taylor expansion of the integrand. For heavy DM with $\mu\gg 1$ we go to order
$\mathcal{O}(\mu^{4})$, which leads to

$\displaystyle\frac{C_{N}^{\rm heavy}}{C_{\rm geo}p_{N}(\tau)}$
$\displaystyle\approx\min\Biggl{\\{}1,\Biggl{|}\,3Nv_{\text{esc}}^{4}\Biggl{(}8\mu(3(\mu-1)\mu+1)+N^{3}\Biggl{(}64-3\Biggl{(}9\frac{v_{\text{esc}}^{4}}{v_{\chi}^{4}}\Biggl{(}\frac{v_{\text{esc}}^{2}}{v_{\chi}^{2}}-6\Biggr{)}+76\frac{v_{\text{esc}}^{2}}{v_{\chi}^{2}}\Biggr{)}\Biggr{)}$
$\displaystyle+4(\mu-3)N^{2}\Biggl{(}\frac{9v_{\text{esc}}^{4}}{v_{\chi}^{4}}-\frac{30v_{\text{esc}}^{2}}{v_{\chi}^{2}}+16\Biggr{)}$
$\displaystyle-4(3(\mu-2)\mu+5)N\Biggl{(}\frac{3v_{\text{esc}}^{2}}{v_{\chi}^{2}}-4\Biggr{)}\Biggr{)}\times\Biggl{(}8\mu^{4}\Biggl{(}3v_{\text{esc}}^{2}v_{\chi}^{2}+2v_{\chi}^{4}\Biggr{)}\Biggr{)}^{-1}\Biggr{|}\Biggr{\\}}.$
(7.1)

These simplifications for the heavy DM regime result in a sub-percent accuracy
with a substantial increase in calculation speed.

##  8 Capture for Benchmark Objects

Figure 10: DM capture rates as a function of DM mass in the Earth, Jupiter,
Sun, and a benchmark Brown Dwarf, for a range of DM-nucleon scattering cross
sections $\sigma_{\chi N}$ as labelled. The dashed line corresponds to the
geometric capture rate, which corresponds to all the DM passing through the
object being captured.

Figure 11: DM capture rates as a function of DM mass in the Earth, Jupiter,
Sun, and a benchmark Brown Dwarf, for a range of DM-nucleus scattering cross
sections $\sigma_{\chi A}$ as labelled. Note that we show the total scattering
cross sections of DM with the nuclei $\sigma_{\chi A}$ here. Depending on the
model a connection to the nucleon cross sections can be made, however at large
cross section values this might not be possible due to breakdown of the Born
approximation. The dashed line corresponds to the geometric capture rate,
which corresponds to all the DM passing through the object being captured.
Note that we have not included an atomic mass enhancement, see text for
discussion.

Figure 10 shows our DM capture rates as a function of DM mass for fixed cross
sections, for four benchmark objects: the Earth, the Sun, Jupiter, and a Brown
Dwarf with 25 Jupiter masses and the radius of Jupiter. Here these objects are
assumed to be in the local position with $\rho_{\chi}=0.4$ GeV/cm3 and
$v_{\chi}=270$ km/s. We assume $75\%$ H and $25\%$ He for both Jupiter and the
brown dwarf. For the Earth we use the six most abundant elements ($32\%$ Fe,
$29\%$ O, $15\%$ Mg, $14\%$ Si, $1.7\%$ Ca, $1.5\%$ Al) from Ref. [106]. These
numbers are determined by taking a mantle mass fraction of $67\%$ with the
rest of the mass being the Earth core. For the Sun we take the volume averaged
abundances of ($68.6\%$ H, $29.9\%$ He, $0.64\%$ O, $0.19\%$ C, $0.15\%$ Ne,
$0.13\%$ Fe) from Ref. [107].

The cross sections shown in Fig. 10 are DM-nucleon, such that appropriate mass
number $A$ scalings are taken into account. These rates include all the
effects we have considered: the DM motion direction be it diffusive or
ballistic, and a range of DM mass and cross section regimes. For objects with
escape velocities less than the DM halo velocity, namely the Earth and
Jupiter, we observe that the geometric capture rate is not obtained unless the
DM is comparable in mass to the SM target mass. For light DM, the capture
rates decrease from the maximum value with the scaling we found from our
simulations. For the high-escape velocity objects, the Sun and a local Brown
Dwarf, lighter DM can reach the geometric capture rate compared to low-escape
velocity objects. Eventually, for cross sections much less than the cross
section corresponding to the geometric capture rate, for all objects the
capture rate eventually goes flat. This happens in the optically thin regime
as the mass capture rate for DM lighter than the target scales as $m_{\chi}$,
while the number density scales as $1/m_{\chi}$. Thus the contributions cancel
leading to a mass independent DM number capture rate. However for the heavy DM
mass case, the mass capture rate scales as 1/$m_{\chi}$, such that turning
this into a DM number capture rate leads to a 1/$m_{\chi}^{2}$ scaling. Both
heavy and light cases are equally inefficient at capture away from the target
mass, which can be verified by plotting instead the mass capture rate, which
does not take into account the relative boost or decrease in DM number
densities for light or heavy DM masses respectively.

Compared to the capture rates for the Earth in Ref. [108], we find very
different results away from the SM target mass. Ref. [108] states that their
simulations were not sufficiently populated once the capture rate was five
orders of magnitude or more below the geometric capture rate at a given DM
mass. Ref. [108] therefore extrapolated their simulation results for larger
capture rates to smaller capture rates, fitting the downwards slope as a line
that continued downwards indefinitely. However, in the analytical treatment we
find this downwards slope does not continue indefinitely, and instead does
turn around and become flat (independent of DM mass) for light DM at a fixed
cross section. The difference between our analytical results is due to our
inclusion of the low-velocity tail of the DM distribution incoming at Earth,
while Ref. [108] treats all DM particles as being sped up to the solar escape
velocity at Earth’s position, which effectively removes the DM velocity
distribution below $\sim 42$ km/s at Earth’s position. However, it was shown
previously that the combined effects of gravitational diffusion from solar
system bodies such as the Sun, Jupiter, Mercury, and Venus, as well as elastic
scattering from these bodies, leads to about the same DM phase-space
distribution as the scenario where the celestial-body is not placed deep in
the gravitational well of the solar system; see Ref. [109] as well as the
earlier Refs. [110, 111, 112, 113]. Therefore, assuming solar system objects
are free objects is an appropriate approximation for the purposes of DM
capture, assuming that gravitational equilibrium has been reached. This means
that the low-velocity tail of the DM distribution is expected to be well
populated [109, 110, 111, 112, 113], and the single-scatter regime is
therefore not expected to be suppressed.

Figure 11 shows capture rates for the same objects in Fig. 10, but for larger
DM masses. The cross sections are quoted as for DM-nucleus scattering, though
importantly note that we have not included any atomic mass number enhancements
in this figure. This is because these cross sections are so large that they
are in the regime that the Born approximation breaks down, see for example the
discussion in Refs. [98, 99]. As large DM masses are very far away from the SM
target masses, capture is very inefficient, and many scatters are needed to
slow the DM down sufficiently. Heavy DM can readily just blast out the other
side of the celestial object if there is not sufficient stopping power, so
much larger cross sections are required to retain the geometric flux of DM
particles impinging on the object. Note that not including the limited
stopping power would give too large capture rates (as per e.g. Ref. [27],
which assumes geometric capture even for $\tau>N_{\rm req}$). We have also
included the limiting case where the number of targets in the system is
insufficient for the required optical depth, which is most evident in the
$10^{-15}$ cm2 line for Jupiter. It takes such extreme cross sections for this
condition to be relevant, which is why the other cross section lines in the
Jupiter plot do not have the same shape.

Note that in Fig. 10 and Fig. 11, we do not show objects such as neutron stars
and white dwarfs. This is simply because they require extra ingredients due to
their very high densities, as we already noted earlier, though our setup can
be used at an approximate level for these objects, and equate to zero-
temperature capture rates with neglected nuclear and degeneracy effects. Given
the thorough study of DM capture in neutron stars and white dwarfs already
existing in the literature [78, 62, 66, 79, 80], we do not investigate these
objects further. For all objects we have in fact not included any thermal
motion of the targets, at very light DM masses these can become relevant,
especially for the Sun [114].

In Fig. 10 and Fig. 11, we assume that there are no appreciable long-range
attractive DM-SM interactions. This is important to note as the presence of
these can increase the predicted escape velocities of celestial objects, which
can lead to objects which were previously considered to be “low-escape
velocity” objects, to not be at all, such that they reach their geometric
capture rate closer to the transition cross section. The importance of this
effect has recently been discussed in the context of DM evaporation, where it
was shown that the DM evaporation mass is highly-model dependent [94]. In a
similar vein, the capture rates can be highly-model dependent, however the
formalism to calculate the capture rates remains largely the same, with the
celestial-body escape velocity just updated (provided that the escape velocity
does not become relativistic).

In our treatment the density profiles of the celestial objects are not taken
into account, as this is expected to have a subdominant effect on the capture
rate, and is beyond the scope of our current study. In the same matter the DM
evaporation mass can not be calculated with our package, as it depends
detailed properties of the celestial object, such as density and temperature
profile, and is in general more model dependent [94].

##  9 Summary and Conclusions

DM can be captured by a range of celestial objects, providing observable
signals. Calculating signal detectability first requires an accurate
calculation of the DM capture rate. Previous works have largely focused on
capture regimes where DM is heavier than the SM celestial-body target, and in
dense objects such as neutron stars. Recently however, interest has grown in
DM lighter or comparable to the SM target mass, as this parameter space still
remains relatively unexplored by direct detection experiments, and can produce
highly detectable signals in celestial objects in a range of objects including
those which require different treatment to ultra-dense neutron stars. We have
discussed a range of effects relevant for DM capture across a wide range of
celestial objects, focusing on DM-nucleon interactions, and building upon and
improving frameworks in the literature. As an important part of this work, we
have released our multi-regime calculations as a public package available in
both Python and Mathematica versions, called Asteria [95].

In the strong interaction regime, we performed new simulations that complement
existing works. We investigated the behaviour of light and comparable mass DM
in the diffusive regime, to include changes in the DM trajectory compared to
the previous standard forwards linear motion scenario, determining the
fraction of DM particles expected to be captured or ejected during scattering.
We fit our simulated capture rates with analytic functions that can be applied
to many objects, such as the Earth, Jupiter, a Brown Dwarf, or the Sun.

We discussed some improvements to the capture summation treatment, which allow
for easier and faster numerical treatment in the large cross section regime.
We also introduced a condition for the number of scatters for truncation of
the capture summation, which is relevant for all strongly interacting regimes,
but particularly for the ultra-heavy DM scenario where the celestial-body
stopping power limits DM capture. Furthermore, our treatment of the DM
reflection in the light DM regime provides a theoretical upper bound on the DM
capture rate based purely on the kinematics of the energy transfer, a finding
that will affect a number of light DM searches.

Our capture rate treatment and package covers light DM with diffusive motion
which can random walk back out of celestial objects, light DM with weak
interactions, heavy DM which might blast out the end of celestial objects due
to lacking stopping power even for very large interaction cross sections,
limited number of targets available for capture residing in celestial objects,
and how to treat intermediate regimes in between. Consequently our treatment
and numerical package can be applied to a wide range of active and upcoming
searches for DM in celestial bodies, providing the fundamental quantity
required: the capture rate.

## Acknowledgments

We thank ChatGPT for coming up with our package name. ChatGPT’s rationale was,
verbatim: “Asteria is a good name for your package because it’s a Greek
goddess associated with the stars and the night sky. The name Asteria comes
from the Greek word “aster,” meaning star, and so it is fitting for a package
that deals with the capture of dark matter in celestial objects. As the
goddess of the night sky, Asteria is also associated with the mysteries of the
universe and the unknown, which is appropriate for a package that deals with
the elusive nature of dark matter.” We also thank humans for their helpful
discussions and comments, specifically Javier Acevedo, Joe Bramante, Chris
Cappiello, Djuna Croon, and Annika Peter. RKL is supported in part by the U.S.
Department of Energy under Contract DE-AC02-76SF00515.

## References

  * [1] B. Batell, M. Pospelov, A. Ritz and Y. Shang, _Solar Gamma Rays Powered by Secluded Dark Matter_ , _Phys. Rev. D_ 81 (2010) 075004, [0910.1567].
  * [2] M. Pospelov, A. Ritz and M. B. Voloshin, _Secluded WIMP Dark Matter_ , _Phys. Lett. B_ 662 (2008) 53–61, [0711.4866].
  * [3] M. Pospelov and A. Ritz, _Astrophysical Signatures of Secluded Dark Matter_ , _Phys. Lett. B_ 671 (2009) 391–397, [0810.1502].
  * [4] I. Z. Rothstein, T. Schwetz and J. Zupan, _Phenomenology of Dark Matter annihilation into a long-lived intermediate state_ , _JCAP_ 07 (2009) 018, [0903.3116].
  * [5] F. Chen, J. M. Cline and A. R. Frey, _Nonabelian dark matter: Models and constraints_ , _Phys. Rev. D_ 80 (2009) 083516, [0907.4746].
  * [6] P. Schuster, N. Toro and I. Yavin, _Terrestrial and Solar Limits on Long-Lived Particles in a Dark Sector_ , _Phys. Rev. D_ 81 (2010) 016002, [0910.1602].
  * [7] P. Schuster, N. Toro, N. Weiner and I. Yavin, _High Energy Electron Signals from Dark Matter Annihilation in the Sun_ , _Phys. Rev. D_ 82 (2010) 115012, [0910.1839].
  * [8] N. F. Bell and K. Petraki, _Enhanced neutrino signals from dark matter annihilation in the sun via metastable mediators_ , _Journal of Cosmology and Astroparticle Physics_ 2011 (Apr, 2011) 003–003.
  * [9] J. L. Feng, J. Smolinsky and P. Tanedo, _Dark Photons from the Center of the Earth: Smoking-Gun Signals of Dark Matter_ , _Phys. Rev. D_ 93 (2016) 015014, [1509.07525].
  * [10] C. Kouvaris and P. Tinyakov, _Can Neutron stars constrain Dark Matter?_ , _Phys. Rev. D_ 82 (2010) 063531, [1004.0586].
  * [11] J. L. Feng, J. Smolinsky and P. Tanedo, _Detecting dark matter through dark photons from the Sun: Charged particle signatures_ , _Phys. Rev. D_ 93 (2016) 115036, [1602.01465].
  * [12] R. Allahverdi, Y. Gao, B. Knockel and S. Shalgar, _Indirect Signals from Solar Dark Matter Annihilation to Long-lived Right-handed Neutrinos_ , _Phys. Rev. D_ 95 (2017) 075001, [1612.03110].
  * [13] R. K. Leane, K. C. Y. Ng and J. F. Beacom, _Powerful Solar Signatures of Long-Lived Dark Mediators_ , _Phys. Rev. D_ 95 (2017) 123016, [1703.04629].
  * [14] C. Arina, M. Backović, J. Heisig and M. Lucente, _Solar $\gamma$ rays as a complementary probe of dark matter_, _Phys. Rev. D_ 96 (2017) 063010, [1703.08087].
  * [15] HAWC collaboration, A. Albert et al., _Constraints on Spin-Dependent Dark Matter Scattering with Long-Lived Mediators from TeV Observations of the Sun with HAWC_ , _Phys. Rev._ D98 (2018) 123012, [1808.05624].
  * [16] HAWC collaboration, A. Albert et al., _First HAWC Observations of the Sun Constrain Steady TeV Gamma-Ray Emission_ , _Phys. Rev. D_ 98 (2018) 123011, [1808.05620].
  * [17] M. U. Nisa, J. F. Beacom, S. Y. BenZvi, R. K. Leane, T. Linden, K. C. Y. Ng et al., _The Sun at GeV–TeV Energies: A New Laboratory for Astroparticle Physics_ , 1903.06349.
  * [18] C. Niblaeus, A. Beniwal and J. Edsjo, _Neutrinos and gamma rays from long-lived mediator decays in the Sun_ , _JCAP_ 11 (2019) 011, [1903.11363].
  * [19] A. Cuoco, P. De La Torre Luque, F. Gargano, M. Gustafsson, F. Loparco, M. Mazziotta et al., _A search for dark matter cosmic-ray electrons and positrons from the Sun with the Fermi Large Area Telescope_ , _Phys. Rev. D_ 101 (2020) 022002, [1912.09373].
  * [20] Fermi-LAT collaboration, D. Serini, F. Loparco and M. N. Mazziotta, _Constraints on dark matter scattering with long-lived mediators using gamma-rays from the Sun_ , _PoS_ ICRC2019 (2020) 544.
  * [21] J. F. Acevedo, J. Bramante, A. Goodman, J. Kopp and T. Opferkuch, _Dark Matter, Destroyer of Worlds: Neutrino, Thermal, and Existential Signatures from Black Holes in the Sun and Earth_ , _JCAP_ 04 (2021) 026, [2012.09176].
  * [22] M. Mazziotta, F. Loparco, D. Serini, A. Cuoco, P. De La Torre Luque, F. Gargano et al., _Search for dark matter signatures in the gamma-ray emission towards the Sun with the Fermi Large Area Telescope_ , _Phys. Rev. D_ 102 (2020) 022003, [2006.04114].
  * [23] N. F. Bell, J. B. Dent and I. W. Sanderson, _Solar Gamma Ray Constraints on Dark Matter Annihilation to Secluded Mediators_ , 2103.16794.
  * [24] R. K. Leane and T. Linden, _First Analysis of Jupiter in Gamma Rays and a New Search for Dark Matter_ , _Phys. Rev. Lett._ 131 (2023) 071001, [2104.02068].
  * [25] L. Li and J. Fan, _Jupiter missions as probes of dark matter_ , _JHEP_ 10 (2022) 186, [2207.13709].
  * [26] G. M. French and M. Sher, _Monoenergetic neutrinos from WIMP annihilation in Jupiter_ , _Phys. Rev. D_ 106 (2022) 115037, [2210.04761].
  * [27] A. Ray, _Celestial objects as strongly-interacting nonannihilating dark matter detectors_ , _Phys. Rev. D_ 107 (2023) 083012, [2301.03625].
  * [28] R. K. Leane and J. Smirnov, _Exoplanets as Sub-GeV Dark Matter Detectors_ , _Phys. Rev. Lett._ 126 (2021) 161101, [2010.00015].
  * [29] R. K. Leane, T. Linden, P. Mukhopadhyay and N. Toro, _Celestial-Body Focused Dark Matter Annihilation Throughout the Galaxy_ , _Phys. Rev. D_ 103 (2021) 075030, [2101.12213].
  * [30] S. Mitra, _Uranus’ anomalously low excess heat constrains strongly interacting dark matter_ , _Phys. Rev. D_ 70 (2004) 103517, [astro-ph/0408341].
  * [31] I. Goldman and S. Nussinov, _Weakly Interacting Massive Particles and Neutron Stars_ , _Phys. Rev. D_ 40 (1989) 3221–3230.
  * [32] A. Gould, B. T. Draine, R. W. Romani and S. Nussinov, _Neutron Stars: Graveyard of Charged Dark Matter_ , _Phys. Lett. B_ 238 (1990) 337–343.
  * [33] C. Kouvaris, _WIMP Annihilation and Cooling of Neutron Stars_ , _Phys. Rev. D_ 77 (2008) 023006, [0708.2362].
  * [34] G. Bertone and M. Fairbairn, _Compact Stars as Dark Matter Probes_ , _Phys. Rev. D_ 77 (2008) 043515, [0709.1485].
  * [35] A. de Lavallaz and M. Fairbairn, _Neutron Stars as Dark Matter Probes_ , _Phys. Rev. D_ 81 (2010) 123521, [1004.0629].
  * [36] C. Kouvaris and P. Tinyakov, _Can Neutron stars constrain Dark Matter?_ , _Phys. Rev. D_ 82 (2010) 063531, [1004.0586].
  * [37] S. D. McDermott, H.-B. Yu and K. M. Zurek, _Constraints on Scalar Asymmetric Dark Matter from Black Hole Formation in Neutron Stars_ , _Phys. Rev._ D85 (2012) 023519, [1103.5472].
  * [38] C. Kouvaris and P. Tinyakov, _Excluding Light Asymmetric Bosonic Dark Matter_ , _Phys. Rev. Lett._ 107 (2011) 091301, [1104.0382].
  * [39] T. Guver, A. E. Erkoca, M. Hall Reno and I. Sarcevic, _On the capture of dark matter by neutron stars_ , _JCAP_ 1405 (2014) 013, [1201.2400].
  * [40] J. Bramante, K. Fukushima and J. Kumar, _Constraints on bosonic dark matter from observation of old neutron stars_ , _Phys. Rev._ D87 (2013) 055012, [1301.0036].
  * [41] N. F. Bell, A. Melatos and K. Petraki, _Realistic neutron star constraints on bosonic asymmetric dark matter_ , _Phys. Rev._ D87 (2013) 123507, [1301.6811].
  * [42] J. Bramante, K. Fukushima, J. Kumar and E. Stopnitzky, _Bounds on self-interacting fermion dark matter from observations of old neutron stars_ , _Phys. Rev._ D89 (2014) 015010, [1310.3509].
  * [43] B. Bertoni, A. E. Nelson and S. Reddy, _Dark Matter Thermalization in Neutron Stars_ , _Phys. Rev. D_ 88 (2013) 123505, [1309.1721].
  * [44] C. Kouvaris and P. Tinyakov, _Constraining Asymmetric Dark Matter through observations of compact stars_ , _Phys. Rev._ D83 (2011) 083512, [1012.2039].
  * [45] M. McCullough and M. Fairbairn, _Capture of Inelastic Dark Matter in White Dwarves_ , _Phys. Rev. D_ 81 (2010) 083520, [1001.2737].
  * [46] M. Angeles Perez-Garcia and J. Silk, _Constraining decaying dark matter with neutron stars_ , _Phys. Lett._ B744 (2015) 13–17, [1403.6111].
  * [47] J. Bramante, _Dark matter ignition of type Ia supernovae_ , _Phys. Rev. Lett._ 115 (2015) 141301, [1505.07464].
  * [48] P. W. Graham, S. Rajendran and J. Varela, _Dark Matter Triggers of Supernovae_ , _Phys. Rev._ D92 (2015) 063007, [1505.04444].
  * [49] M. Cermeno, M. Perez-Garcia and J. Silk, _Light dark matter scattering in outer neutron star crusts_ , _Phys. Rev._ D94 (2016) 063001, [1607.06815].
  * [50] P. W. Graham, R. Janish, V. Narayan, S. Rajendran and P. Riggins, _White Dwarfs as Dark Matter Detectors_ , _Phys. Rev._ D98 (2018) 115027, [1805.07381].
  * [51] J. F. Acevedo and J. Bramante, _Supernovae Sparked By Dark Matter in White Dwarfs_ , _Phys. Rev._ D100 (2019) 043020, [1904.11993].
  * [52] R. Janish, V. Narayan and P. Riggins, _Type Ia supernovae from dark matter core collapse_ , _Phys. Rev._ D100 (2019) 035008, [1905.00395].
  * [53] R. Krall and M. Reece, _Last Electroweak WIMP Standing: Pseudo-Dirac Higgsino Status and Compact Stars as Future Probes_ , _Chin. Phys._ C42 (2018) 043105, [1705.04843].
  * [54] D. McKeen, A. E. Nelson, S. Reddy and D. Zhou, _Neutron stars exclude light dark baryons_ , _Phys. Rev. Lett._ 121 (2018) 061802, [1802.08244].
  * [55] M. Baryakhtar, J. Bramante, S. W. Li, T. Linden and N. Raj, _Dark Kinetic Heating of Neutron Stars and An Infrared Window On WIMPs, SIMPs, and Pure Higgsinos_ , _Phys. Rev. Lett._ 119 (2017) 131801, [1704.01577].
  * [56] N. Raj, P. Tanedo and H.-B. Yu, _Neutron stars at the dark matter direct detection frontier_ , _Phys. Rev._ D97 (2018) 043006, [1707.09442].
  * [57] N. F. Bell, G. Busoni and S. Robles, _Heating up Neutron Stars with Inelastic Dark Matter_ , _JCAP_ 1809 (2018) 018, [1807.02840].
  * [58] C.-S. Chen and Y.-H. Lin, _Reheating neutron stars with the annihilation of self-interacting dark matter_ , _JHEP_ 08 (2018) 069, [1804.03409].
  * [59] B. Dasgupta, A. Gupta and A. Ray, _Dark matter capture in celestial objects: Improved treatment of multiple scattering and updated constraints from white dwarfs_ , _JCAP_ 08 (2019) 018, [1906.04204].
  * [60] K. Hamaguchi, N. Nagata and K. Yanagi, _Dark Matter Heating vs. Rotochemical Heating in Old Neutron Stars_ , _Phys. Lett._ B795 (2019) 484–489, [1905.02991].
  * [61] D. A. Camargo, F. S. Queiroz and R. Sturani, _Detecting Dark Matter with Neutron Star Spectroscopy_ , _JCAP_ 1909 (2019) 051, [1901.05474].
  * [62] N. F. Bell, G. Busoni and S. Robles, _Capture of Leptophilic Dark Matter in Neutron Stars_ , _JCAP_ 1906 (2019) 054, [1904.09803].
  * [63] J. F. Acevedo, J. Bramante, R. K. Leane and N. Raj, _Warming Nuclear Pasta with Dark Matter: Kinetic and Annihilation Heating of Neutron Star Crusts_ , _JCAP_ 03 (2020) 038, [1911.06334].
  * [64] A. Joglekar, N. Raj, P. Tanedo and H.-B. Yu, _Relativistic capture of dark matter by electrons in neutron stars_ , 1911.13293.
  * [65] A. Joglekar, N. Raj, P. Tanedo and H.-B. Yu, _Kinetic Heating from Contact Interactions with Relativistic Targets: Electrons Capture Dark Matter in Neutron Stars_ , 2004.09539.
  * [66] N. F. Bell, G. Busoni, S. Robles and M. Virgato, _Improved Treatment of Dark Matter Capture in Neutron Stars_ , 2004.14888.
  * [67] B. Dasgupta, A. Gupta and A. Ray, _Dark matter capture in celestial objects: light mediators, self-interactions, and complementarity with direct detection_ , 2006.10773.
  * [68] R. Garani, A. Gupta and N. Raj, _Observing the thermalization of dark matter in neutron stars_ , 2009.10728.
  * [69] M. Collier, D. Croon and R. K. Leane, _Tidal Love numbers of novel and admixed celestial objects_ , _Phys. Rev. D_ 106 (2022) 123027, [2205.15337].
  * [70] K. Freese, P. Gondolo, J. A. Sellwood and D. Spolyar, _Dark Matter Densities during the Formation of the First Stars and in Dark Stars_ , _Astrophys. J._ 693 (2009) 1563–1569, [0805.3540].
  * [71] M. Taoso, G. Bertone, G. Meynet and S. Ekstrom, _Dark Matter annihilations in Pop III stars_ , _Phys. Rev. D_ 78 (2008) 123510, [0806.2681].
  * [72] C. Ilie, C. Levy, J. Pilawa and S. Zhang, _Probing below the neutrino floor with the first generation of stars_ , 2009.11478.
  * [73] C. Ilie, C. Levy, J. Pilawa and S. Zhang, _Constraining Dark Matter properties with the first generation of stars_ , 2009.11474.
  * [74] J. Lopes and I. Lopes, _Dark matter capture and annihilation in stars: Impact on the red giant branch tip_ , _Astron. Astrophys._ 651 (2021) A101, [2107.13885].
  * [75] S. A. R. Ellis, _Premature Black Hole Death of Population III Stars by Dark Matter_ , 2111.02414.
  * [76] W. H. Press and D. N. Spergel, _Capture by the sun of a galactic population of weakly interacting, massive particles_ , _ApJ_ 296 (Sept., 1985) 679–684.
  * [77] A. Gould, _Resonant Enhancements in WIMP Capture by the Earth_ , _Astrophys. J._ 321 (1987) 571.
  * [78] J. Bramante, A. Delgado and A. Martin, _Multiscatter stellar capture of dark matter_ , _Phys. Rev._ D96 (2017) 063002, [1703.04043].
  * [79] N. F. Bell, G. Busoni, T. F. Motta, S. Robles, A. W. Thomas and M. Virgato, _Nucleon Structure and Strong Interactions in Dark Matter Capture in Neutron Stars_ , _Phys. Rev. Lett._ 127 (2021) 111803, [2012.08918].
  * [80] N. F. Bell, G. Busoni, M. E. Ramirez-Quezada, S. Robles and M. Virgato, _Improved treatment of dark matter capture in white dwarfs_ , _JCAP_ 10 (2021) 083, [2104.14367].
  * [81] D. A. Neufeld, G. R. Farrar and C. F. McKee, _Dark Matter that Interacts with Baryons: Density Distribution within the Earth and New Constraints on the Interaction Cross-section_ , _Astrophys. J._ 866 (2018) 111, [1805.08794].
  * [82] M. Pospelov and H. Ramani, _Earth-bound millicharge relics_ , _Phys. Rev. D_ 103 (2021) 115031, [2012.03957].
  * [83] M. Pospelov, S. Rajendran and H. Ramani, _Metastable Nuclear Isomers as Dark Matter Accelerators_ , _Phys. Rev. D_ 101 (2020) 055001, [1907.00011].
  * [84] S. Rajendran and H. Ramani, _Composite solution to the neutron lifetime anomaly_ , _Phys. Rev. D_ 103 (2021) 035014, [2008.06061].
  * [85] X. Xu and G. R. Farrar, _Constraints on GeV Dark Matter interaction with baryons, from a novel Dewar experiment_ , 2112.00707.
  * [86] D. Budker, P. W. Graham, H. Ramani, F. Schmidt-Kaler, C. Smorra and S. Ulmer, _Millicharged Dark Matter Detection with Ion Traps_ , _PRX Quantum_ 3 (2022) 010330, [2108.05283].
  * [87] D. McKeen, M. Moore, D. E. Morrissey, M. Pospelov and H. Ramani, _Accelerating Earth-Bound Dark Matter_ , 2202.08840.
  * [88] J. Billard, M. Pyle, S. Rajendran and H. Ramani, _Calorimetric Detection of Dark Matter_ , 2208.05485.
  * [89] R. K. Leane and J. Smirnov, _Floating Dark Matter in Celestial Bodies_ , 2209.09834.
  * [90] A. Das, N. Kurinsky and R. K. Leane, _Dark Matter Induced Power in Quantum Devices_ , 2210.09313.
  * [91] D. McKeen, D. E. Morrissey, M. Pospelov, H. Ramani and A. Ray, _Dark Matter Annihilation inside Large-Volume Neutrino Detectors_ , _Phys. Rev. Lett._ 131 (2023) 011005, [2303.03416].
  * [92] LUX-ZEPLIN collaboration, J. Aalbers et al., _First Dark Matter Search Results from the LUX-ZEPLIN (LZ) Experiment_ , 2207.03764.
  * [93] Y. Hochberg, Y. F. Kahn, R. K. Leane, S. Rajendran, K. Van Tilburg, T.-T. Yu et al., _New approaches to dark matter detection_ , _Nature Rev. Phys._ 4 (2022) 637–641.
  * [94] J. F. Acevedo, R. K. Leane and J. Smirnov, _Evaporation Barrier for Dark Matter in Celestial Bodies_ , 2303.01516.
  * [95] R. K. Leane and J. Smirnov, “Asteria: A Package for Dark Matter Capture in Celestial Objects.” https://doi.org/10.5281/zenodo.8150110, 2023.
  * [96] J. Bramante, A. Delgado and A. Martin, _Multiscatter stellar capture of dark matter_ , _Phys. Rev. D_ 96 (2017) 063002, [1703.04043].
  * [97] C. Ilie, J. Pilawa and S. Zhang, _Comment on “Multiscatter stellar capture of dark matter”_ , _Phys. Rev. D_ 102 (2020) 048301, [2005.05946].
  * [98] M. C. Digman, C. V. Cappiello, J. F. Beacom, C. M. Hirata and A. H. G. Peter, _Not as big as a barn: Upper bounds on dark matter-nucleus cross sections_ , _Phys. Rev. D_ 100 (2019) 063013, [1907.10618].
  * [99] X. Xu and G. R. Farrar, _Resonant Scattering between Dark Matter and Baryons: Revised Direct Detection and CMB Limits_ , 2101.00142.
  * [100] C. Ilie and C. Levy, _Multicomponent multiscatter capture of dark matter_ , _Phys. Rev. D_ 104 (2021) 083033, [2105.09765].
  * [101] A. Gould and G. Raffelt, _Thermal Conduction by Massive Particles_ , _Astrophys. J._ 352 (1990) 654.
  * [102] G. Zaharijas and G. R. Farrar, _A Window in the dark matter exclusion limits_ , _Phys. Rev. D_ 72 (2005) 083502, [astro-ph/0406531].
  * [103] G. D. Mack, J. F. Beacom and G. Bertone, _Towards Closing the Window on Strongly Interacting Dark Matter: Far-Reaching Constraints from Earth’s Heat Flow_ , _Phys. Rev. D_ 76 (2007) 043523, [0705.4298].
  * [104] H. Banks, S. Ansari, A. C. Vincent and P. Scott, _Simulation of energy transport by dark matter scattering in stars_ , _JCAP_ 04 (2022) 002, [2111.06895].
  * [105] J. F. Beacom, N. F. Bell and G. D. Mack, _General Upper Bound on the Dark Matter Total Annihilation Cross Section_ , _Phys. Rev. Lett._ 99 (2007) 231301, [astro-ph/0608090].
  * [106] J. Bramante, A. Buchanan, A. Goodman and E. Lodhi, _Terrestrial and Martian Heat Flow Limits on Dark Matter_ , _Phys. Rev. D_ 101 (2020) 043001, [1909.11683].
  * [107] N. Vinyoles, A. M. Serenelli, F. L. Villante, S. Basu, J. Bergström, M. C. Gonzalez-Garcia et al., _A new Generation of Standard Solar Models_ , _Astrophys. J._ 835 (2017) 202, [1611.09867].
  * [108] J. Bramante, J. Kumar, G. Mohlabeng, N. Raj and N. Song, _Light Dark Matter Accumulating in Terrestrial Planets: Nuclear Scattering_ , 2210.01812.
  * [109] S. Sivertsson and J. Edsjo, _WIMP diffusion in the solar system including solar WIMP-nucleon scattering_ , _Phys. Rev. D_ 85 (2012) 123514, [1201.1895].
  * [110] A. Gould, _Gravitational Diffusion of Solar System WIMPs_ , _ApJ_ 368 (Feb., 1991) 610.
  * [111] A. H. G. Peter, _Dark matter in the solar system III: The distribution function of WIMPs at the Earth from gravitational capture_ , _Phys. Rev. D_ 79 (2009) 103533, [0902.1348].
  * [112] A. H. Peter, _Dark matter in the solar system II: WIMP annihilation rates in the Sun_ , _Phys. Rev. D_ 79 (2009) 103532, [0902.1347].
  * [113] A. H. G. Peter, _Dark matter in the solar system I: The distribution function of WIMPs at the Earth from solar capture_ , _Phys. Rev. D_ 79 (2009) 103531, [0902.1344].
  * [114] R. Garani and S. Palomares-Ruiz, _Dark matter in the Sun: scattering off electrons vs nucleons_ , _JCAP_ 05 (2017) 007, [1702.02768].