data have found 11 new streams within 50 kpc of our location in the MW, these
streams are observed down to surface brightness of 23 $\rm mag/arcsec$ [402].
While currently DES is the state-of-the-art, with the most sensitive wide-area
view, with the advent of the Vera C. Rubin Observatory we are expected to see
of the order of a hundred streams out to the virial radius of the MW halo
[402]. Moreover, spectroscopic programs like the S5 survey and instruments
like the Maunakea Spectroscopic Explorer telescope will allow us to follow up
on some of the streams with spectroscopic data to observe the disruptions of
the streams in velocity space [403, 404, 405].

###  4.3 X-ray and Weak Lensing Observations of Clusters, Groups and Large
Ellipticals

The overall distribution of DM in halos can be probed in a variety of ways.
While rotation curves of galaxies can be measured where spectroscopic data is
available, the gas distribution in X-ray and gravitational lensing have been
used to probe the overall distribution and density profile of the DM out to
the viral radius for larger systems. Current and future generation surveys
like the Dark Energy Survey (DES; [406])999http://www.darkenergysurvey.org,
Hyper Suprime Camera (HSC; [407])101010http://hsc.mtk.nao.ac.jp/ssp/ survey
and LSST will provide precision measurements of the DM halo density profiles
through weak lensing of galaxies in the back ground of large, statistical
samples of cluster and group mass halos.

While the formation of a core is the most dramatic effect on the density
profile of DM halos in the presence of self-interactions, mass conservation
implies that the density at distances larger than the core becomes steeper
than the density profile in the CDM scenario. While the latter is a relatively
smaller effect, the scale at which this feature is expected, especially for
clusters, is less affected by baryonic physics than the core. Further, the
projected density profiles of clusters on these scales can be measured
accurately through weak lensing measurements. Fig. 12 from [95] shows the
stacked projected density profile, $\Delta\Sigma$ of cluster mass objects from
$N$-body simulations of SIDM with different cross sections. The red curve in
the figure corresponds to the current best constraints on the cross section
from the Bullet cluster [107] of $2$ cm${}^{2}/$g. The bottom panel shows the
ratio of the expected weak lensing signal for various SIDM cross sections and
CDM. The lightest grey bands show the measurement error bars for a similar
sample of objects from the DES Year 1 data. The progressively darker bands
show the projected error bars from DES Year 3 and LSST, taking into account
their deeper coverage. This observable, therefore, can potentially provide
much stronger constraints on the SIDM cross section at cluster scales,
especially once data from LSST becomes available.

Figure 12: The stacked lensing signal, $\Delta\Sigma$, of cluster mass DM
halos from simulations. The scales plotted are just just beyond the scale
radius of the DM halo. The density profiles steepen in this region. The
lightest grey band is the current measurement error from DES Y1 scaled to Y3
and LSST sky coverage as they get darker respectively. (This plot has been
recreated based on Ref. [95] with additional estimates for LSST).

Apart from the overall density profile, as mentioned previously DM self-
interactions can alter the shape of the DM halo, making the regions of higher
density more spherically symmetric. Constraints on scattering cross sections
from shapes of halos have been obtained using the X-Ray measurements of the
shapes of the hot halo gas of NGC720. The shape and the twisting of X-Ray
isophotes carry information about the three dimensional shape of the
potential. Measurements from NGC720 constrain the halo ellipticity to be $\sim
0.35\pm 0.03$, consistent with elliptical halos from CDM, this values is
consistent with $\sigma/m_{\text{DM}}<1$ cm2/g [165].

While the shape of DM halos was previously thought to be a promising probe of
the nature of DM, subsequent studies have found that introducing baryonic
physics can significantly complicate the inferences of cross sections from
halo shapes. For example, on galaxy scales [121] show that the asymmetric
potentials from central disks can make the shape of the DM halos aspherical
near the center on short timescales of a few Gyr depending on the size of the
disk. Note also that on cluster mass scales, while the DM can remain spherical
near the center for cross sections of $\sigma/m_{\text{DM}}<1$ cm${}^{2}/$g,
as shown in Fig. 7 from [164], the stars and the gas within the halo remain
aspherical for a longer time and their shapes are much slower to respond to
the change in the DM potential.

Weak lensing measurements of shapes of DM halos, stacked on the galaxy light
profiles or using three point functions can measure the shape of dark halo
directly as a function of radius. [408, 409, 410] measure shapes of Luminous
Red Galaxies and cluster mass halos using weak lensing and distribution of
satellite galaxies and find that the ellipticities are consistent with CDM,
[410] for example place strong constraints on the axis ratios of RedMaPPer
clusters in SDSS at $0.573\pm 0.039$, implying that the halo is overall quite
aspherical.

###  4.4 Major and Minor Mergers in Groups and Clusters (and associated
observables)

Just like terrestrial particle accelerators allow us to study the interactions
of fundamental particles, merging galaxy clusters can be thought of as “cosmic
accelerators” probing the properties of DM. Indeed, the Bullet Cluster is
arguably among the most sensitive and certainly the most famous probe of DM
self-interactions [169]. What makes this system so unique is that the gas in
the smaller of the two merging clusters exhibits a shock wave, which enables
us to conclude with certainty that the system is being observed shortly after
core passage. Combining this information with a lensing reconstruction of the
system’s DM profile [411] enables us to infer the velocity and impact
parameter of the collision [412, 413, 414]. Typical estimates yield a velocity
of around $v_{\mathrm{c}}\sim 4000\,\mathrm{km\,s^{-1}}$ at core passsage and
a projected (i.e. integrated) DM density of $\Sigma\sim
0.3\,\mathrm{g\,cm^{-2}}$ along the path of the subcluster.

While the frequent claim that the Bullet Cluster demonstrates DM to be
collisionless is clearly exaggerated, it is true that it places strong bounds
on the DM self-interaction cross section. Given the projected DM density of
the main cluster halo along the path of the subcluster, the main halo is
expected to become optically thick for DM self-interactions with
$\sigma/m_{\text{DM}}\gtrsim 3\,\rm{cm^{2}\,g^{-1}}$. This means that cross
sections larger than this will lead to almost complete destruction of the
subcluster’s DM halo, as can be seen in Fig. 13 for the case of a simulated
Bullet Cluster with $\sigma/m_{\text{DM}}=10\,\rm{cm^{2}\,g^{-1}}$.

The most stringent SIDM constraint from the Bullet Cluster comes from the
measurement of the mass-to-light ratio ($M/L$) within $150\,\mathrm{kpc}$ of
the subcluster centre, which is found to be consistent both with that of the
main cluster and with the typical range of $M/L$ ratios for clusters [169].
Ref. [167] argues that this agreement places an upper bound on the integrated
mass loss by the subcluster and that it cannot have lost more than 20–30
percent of its mass within the $150\,\mathrm{kpc}$ aperture. The resulting
bound on isotropic, hard-sphere scattering is approximately
$\sigma/m_{\text{DM}}\lesssim 0.7\,\mathrm{cm^{2}\,g^{-1}}$. It should be
noted that this result is subject to some caveats, most notably that there is
significant intrinsic scatter in $M/L$ of isolated clusters. Ref. [415] find a
55% scatter in $M/L$ within the virial radius of clusters, which would
presumably increase if only looking within a small aperture. The consistency
in $M/L$ between the main cluster and subcluster therefore does not exclude
the subcluster having lost a significant fraction of its pre-collision DM
mass, if it started with high $M/L$.

An independent constraint arises from the fact that the peak of the
subcluster’s projected density inferred from weak lensing measurements
coincides with the galaxy centroid from visible light to within $1\sigma$.
This observation rules out the presence of a substantial drag force, which
would lead to the deceleration of the DM halo and hence a separation between
DM and galaxies [167]. The observational results is that the DM is lagging
behind the galaxies by $25\pm 29\,\mathrm{kpc}$ [416, 154]. The first
numerical simulations of the system with SIDM translated this observation into
the constraint $\sigma/m_{\text{DM}}<1.25\,\,\mathrm{cm^{2}\,g^{-1}}$ [154].
However, more recent work has found this constraint to be over-stated, with a
cross section of this magnitude producing offsets of only
$10-15\,\mathrm{kpc}$ [91, 107]. This difference is driven by the methods used
to find the “positions” of the subcluster’s DM and galaxy components, with the
latter papers using methods more closely related to what was done
observationally.

Figure 13: Press release image of the Bullet Cluster (top) compared with
simulations of a Bullet Cluster-like system [92] with CDM (bottom-left) and
SIDM with $\sigma/m_{\text{DM}}=10\,\rm{cm^{2}\,g^{-1}}$ (bottom-right). Blue
shows a map of the DM (inferred from gravitational lensing in the observed
case), while pink shows the X-ray luminosity. The observed optical image has
been overlaid onto the simulation images. The simulation images were made from
the simulation snapshots where the separation between the two peaks in the
galaxy distribution matched that in the observed case. Press-release image
credit: X-ray: NASA/CXC/CfA/[417]; Optical: NASA/STScI;
Magellan/U.Arizona/[416]; Lensing Map: NASA/STScI; ESO WFI;
Magellan/U.Arizona/[416]. Videos of the simulations evolving are available at
https://www.youtube.com/watch?v=rLx_TXhTXbs and
https://www.youtube.com/watch?v=stvkDJVHPuQ.

In spite of the discovery of numerous merging galaxy clusters in recent years,
the Bullet Cluster remains the single most sensitive system, owing mostly to
the superior projected density of the main cluster. Progress can potentially
be made by the statistical combination of a large number of major and minor
mergers, with the first attempt at this producing an upper limit on the SIDM
cross section of $\sigma/m_{\text{DM}}<0.47\,\,\mathrm{cm^{2}\,g^{-1}}$ [112].
However, both the analytical model used to relate a DM–galaxy offset
(strictly, the ratio of the DM–galaxy offset to the gas–galaxy offset) to a
cross section, as well as the quality of the lensing maps used to infer DM
positions, have been called into question [173]. It remains to be seen whether
further discoveries of merging clusters can improve on the bounds obtained
from the Bullet Cluster.

One system that has received particular interest in recent years is the Abell
3827 cluster. The advantage of this system is that strong gravitational
lensing allows for a very precise reconstruction of the DM distribution. It
has thus been possible to identify four galaxies within the core of the main
cluster that all retain their own DM halos [172, 418]. Intriguingly, early
studies found a substantial separation between at least one of these DM halos
and the associated stars [419], which has been interpreted as positive
evidence for SIDM [420]. Unfortunately, more detailed observations have not
confirmed this separation and seem to imply that Abell 3827 is fully
consistent with the predictions for collisionless DM [421], in which the peaks
in stellar and DM density should be coincident [422].

Another controversial system is the ‘Train Wreck Cluster’ Abell 520, which has
been claimed to exhibit a dark core, i.e. a region of _enhanced_ $M/L$ ratio
close to the collision point [423, 424] (see however Ref. [425] for a
contradicting claim). While the most plausible explanation for such a dark
core may be a projection effect due to the presence of filaments [426], an
interpretation in terms of SIDM is also conceivable. It is, however, unclear
how elastic scattering of DM particles could lead to the formation of such a
dark core, since the scattered DM particles would typically have too much
kinetic energy to become gravitationally bound to each other. While it may be
possible in principle to accommodate such an observation in more complex SIDM
models (as discussed in Section 5), it remains challenging to bring the
required self-interaction cross section into agreement with the bounds from
the Bullet Cluster.

Finally, we mention the Musket Ball Cluster, which owes its name to the fact
that it is both older and slower than the Bullet Cluster [427]. This system
was initially considered to be particularly sensitive to DM self-interactions,
due to the fact that it is much further progressed than other merging
clusters. More detailed studies have however revealed, that it is typically
not advantageous to observe a system long after core passage, because the DM
halos quickly relax to their equilibrium position, so that the separation
between DM and stars decreases rather than increases with time [91, 149, 92].

To summarize, systems of merging clusters with the following properties would
make them ideal candidates for constraining DM self-interactions:

  1. 1.

The integrated density of DM along the axis of merger for the system should be
as large as possible, so that self-interactions have a larger effect.

  2. 2.

The system should allow for a precise reconstruction of the DM profile via
(strong) gravitational lensing so that any offset between the light profile
and the center-of-mass of the DM can be accurately measured.

  3. 3.

The system should be observed as shortly after core passage as possible, to
prevent the DM particles from relaxing back to their equilibirum positions.

  4. 4.

The collision should be observed face on to minimize projection effects.

Even if no single system is discovered that surpasses the Bullet Cluster in
terms of sensitivity to SIDM, progress can still be made through the
statistical combination of many colliding galaxy clusters. Expanding the
catalogue of such systems is therefore a highly promising avenue in the search
for SIDM.

###  4.5 Incidence of Multiple BCGs, their Stellar Kinematics and their
Spatial Separations

In Sec. 3.3, we discussed the counterintuitive result that galaxies do not
virialize if bound to the core of a cored DM halo. In simulations of equal-
mass mergers, Ref. [149] found that such galaxies oscillate about the center
of the halo with an amplitude equivalent to the core radius. The easiest
galaxies to definitively associate with cores are the biggest, the BCGs. If
the models of Ref. [149] are correct and hold for more typical mergers where
the masses of the two merging systems can be significantly different, then
BCGs should be offset in both position and velocity from the center of a
relaxed cluster. These studies from dark-matter-only simulations were followed
up in [428] to show that the effect of the BCG “wobble" in a cored profile is
present even in full hydrodynamic simulations.

Figure 14: Figure from [428] showing the evolution of the expected log-normal
median DM-BCG offset in observations. The magenta symbol is the median BCG
offset from the Local supercluster substructure survey measured in [429].

It is in this context that some key observational results about galaxy cluster
centers become even more interesting. Because measurement of the density
profile of clusters from weak lensing depends critically in accurate
determination of cluster centers, the weak lensing community is highly
invested in identifying the best observational tracer of halo centers. From
numerous studies, the main finding is that the position of the central galaxy
(usually, but not always, the BCG [430, 431, 432]) is by far the least noisy
and most accurate tracer of the cluster center [433, 434, 431]. Ref. [431]
finds that the median offset between the BCG and the X-ray center of massive
clusters is about 10 kpc. If true, and if the equal-mass merger prediction of
Ref. [149] holds more generally, the implied constraint on the SIDM hard-
sphere cross section is of order $0.1\leavevmode\nobreak\ $cm2/g. Ref. [429]
claim a detection of BCG offsets, which can be verified with the large and
well-studied cluster samples from current wide-field optical surveys. The
measurement from [429] is presented in Fig. 14, along with the simulation
results for the offset from [428]. The current constraints from the
measurement of BCG wobbles in the cores of galaxy clusters is
$\sigma/m_{\text{DM}}<0.39\ \rm cm^{2}/g$ at the 95% C.L. [428] and data from
future surveys like Euclid are expected to further improve these tight bounds.

There are several important lines of investigation that need to be pursued to
turn this promising probe of SIDM into a robust observable. First, more
simulations are required to determine how generic the galaxy sloshing in cores
is as a function of cluster assembly history (e.g., Ref. [49]). Second, it is
important that the position and velocity offsets be evaluated in the context
of both CDM and SIDM theory. This requires realistic mock observations to be
created from simulations. Third, from the observational side, the position and
velocity offsets for BCGs should be measured more precisely, and as a function
of halo mass. In the case of position offsets, enormous datasets are already
in hand from various imaging surveys, on account of the great interest in
clusters in the dark energy/wide-field survey community. Fourth, these offset
measurements should be performed not just for a single BCG, but for the
handful of brightest cluster galaxies, to trace the merger history of the
cluster. In the extreme case of an equal-mass merger with each component
initially seeded with its own BCG, we expect two sloshing BCGs in the merger
remnant sitting on opposite sides of the cluster center.

###  4.6 Dwarf Galaxies in the Local Group and beyond

Dwarf satellite galaxies are observed around the MW and Andromeda [435, 436],
and around galaxies in the local Universe [437].  In hierarchical structure
formation scenarios, these dwarfs are hosted by subhalos bound to the parent
galaxy [32, 31].  There are about 60 known satellite galaxies of the MW [438],
many of which have measurements of stellar kinematics from large ground-based
telescopes and orbital properties (such as pericenter distances) constrained
with Gaia data. Several surveys are now able to detect satellites in MW-like
galaxies in the Local Volume and this number will grow in the coming years
[439, 437, 440, 441, 56, 442, 443, 444, 445, 446, 447, 448, 449, 450]. The
implications of these measurements for general SIDM models remains to be fully
explored.

The dispersion-dominated dwarf galaxies—known as dwarf Spheroidals (dSphs)—in
the MW span a large range in stellar sizes and luminosities, from ultra-faint
dwarfs such as Segue 1 with stellar mass of ${\cal O}(10^{2})M_{\odot}$ and
half-light radius of $30\rm\leavevmode\nobreak\ pc$ to Fornax with about
$10^{7}M_{\odot}$ of stellar mass within about a kpc. This is shown in the
left panel of Fig. 15. The dynamical masses of these dSphs are significantly
in excess of their stellar masses, making them some of the most dark-matter-
dominated systems observed. A simple way to estimate the dynamical masses of
the dSphs is to use the results of Refs. [451, 452], who argue that the mass
within the half-light radius can be related to the observed stellar dispersion
for a wide range of halo profiles and anisotropies (see also [453]). In the
right panel of Fig. 15, densities within the half-light radius estimated in
this way are shown revealing a striking anti-correlation of density with half-
light radius [131, 454].

Comparisons between simulated subhalos and observed dSph populations reveal
apparent discrepancies that have been the subject of intense scrutiny for the
past two decades. The most obvious mismatch is between the number of simulated
subhalos with virial mass larger than $M_{\rm vir}\gtrsim 10^{8}$ M⊙ (the mass
corresponding to the smallest virial temperature, $T_{\rm vir}\gtrsim 10^{4}$
K, that enables collisional cooling of neutral hydrogen gas and hence star
formation), and the much smaller number of observed Galactic satellites. This
“missing satellites" problem [20, 19] provided motivation to models like
‘warm’ dark matter (WDM), that invoke primordial velocity dispersions to
suppress gravitational collapse on scales below the thermal free streaming
length [282, 455]. Other proposed solutions invoke the astrophysics of galaxy
formation; in these scenarios, the vast majority of DM subhalos fail to host
galaxy formation, due to a combination of factors including mass thresholds,
orbital history and cosmic reionization [35, 456]. Incompleteness and biases
in observations have also been a significant factor in exacerbating the
problem [457, 65].

A satisfactory resolution of the luminosity function of the satellites in the
context of the mass function of subhalos is still missing. Recent work has
pointed out that there may be, in fact, no missing satellites problem [458,
68] and that the radial distribution of the satellites can be explained in the
context of CDM models. This requires, however, that satellites be hosted in
subhalos with masses below the HI cooling temperature [69, 459]—a puzzling
aspect that deserves further study. This conclusion seems to critically depend
on the radial distribution of satellites and the impact of the Large
Magellanic Cloud, with Ref. [460] finding no significant evidence for
populating satellites below the HI cooling limit. Another critical factor is
the limited resolution of N-body simulations, which can lead to artificial
destruction of what could be ultra-faint dSphs at the present time at radii
below 50 kpc [461]. These are clearly important issues to sort out and the
development of detailed models like those in Refs. [461, 460] will facilitate
this exciting science as more satellites are discovered by the Rubin
Observatory.

The key take-away message is that there are recipes for matching the
luminosity function of satellites (taking into account the observational
incompleteness) to the subhalo mass function. However, a more detailed
understanding of the physics responsible for this empirical matching is not
yet available and without that, we cannot claim to have fully understood the
“missing satellites” problem. Regardless of the eventual understanding, it is
clear now that there is limited room for models that do not form subhalos
below a certain mass scale (e.g., WDM, SIDM with additional hidden sector
interactions or fuzzy DM [151, 462, 150, 151]).

In the typical scenario where the luminosity function is abundance-matched to
the subhalo mass functions in CDM simulations to assign dwarf galaxies to
subhalos, the brighest dwarfs would reside in the most massive subhalos.
However, the measurements of the central density of dwarf spheroidals (dSphs)
are less than the central density of subhalos in CDM simulations (dark-matter-
only). Thus, the stellar kinematics of dSphs imply that the most massive
subhalos in the simulations are dark, which is a puzzle as they are "too big
to fail" to form stars [463, 464] (i.e., they would have retained their gas
and cooled to form stars).

Figure 15: Plot on the left from Ref. [438] shows the MW dSphs in the V-band
magnitude vs. half-light radius plane showing how deeper surveys have led to
more dwarfs being uncovered and how the dSphs properties compare against MW
globular clusters. Plot on the right from Ref. [131] shows the inferred total
mass density (almost all DM) within the half-light radius and the strong anti-
correlation between the mass density and half-light radius and the large
diversity in the properties of these dSphs. In the background are subhalo mass
profiles from a MW halo (Kauket from ELVIS) [465].

There are multiple resolutions to the too-big-to-fail problem for satellite
galaxies. The most promising among them is that the discrepancy is just an
artifact produced by comparing observations to the dark-matter-only
simulations. Simulations with baryons indicate that the subhalo density
profiles are much flatter than in dark-matter-only CDM simulations. A flatter
density profile goes a long way toward matching the central density, as
measured observationally, to the total halo mass [466, 467, 468]. Moreover,
including the disk of the MW is known to have an impact on the satellites that
have pericenters close to the disk [469, 467]. However, not all dSphs have
come close to the disk of the MW, so the validity of this solution remains to
be fully tested. All the key pieces of the too-big-to-fail problem in CDM and
SIDM models, with the exception of core collapse is visible in Fig. 16 where
the impact of the disk and pericenter distances in models with with
$\sigma/m_{\text{DM}}=0\leavevmode\nobreak\ \rm{and}\leavevmode\nobreak\ 1\rm
cm^{2}/g$ are shown.

The presence of elastic self-interactions changes the evolution of a satellite
within the host halo in many different ways as listed below and discussed in
more detail in Refs. [152, 134].

  * •

Dark ram-pressure stripping: the high velocity interactions of the subhalo
particles with the host halo particles can lead to mass being unbound from the
subhalo.

  * •

Deceleration of the subhalo due to self-interactions and enhanced dynamical
friction because of orbital decay.

  * •

Thermal evolution, including core collapse and expansion, of the core impacted
by the mass removal by tides.

  * •

Response of the DM and stellar density profiles to the removal of mass by
tidal stripping and ram-pressure stripping.

The different aspects of the orbital and thermal evolution of a subhalo are
sensitive to the velocity dependence of the self-interaction cross section
[153, 470, 127] and the presence of the disk [471, 472]. The thermalization of
the subhalo is controlled by the cross section at relative velocities of about
$4V_{\rm max,\ subhalo}/\sqrt{\pi}$, while the ram-pressure stripping is
controlled by the cross section at relative velocities of order $V_{\rm max,\
main\ halo}$[75, 155, 153]. The gravothermal evolution of the SIDM subhalos is
complicated by the fact that the subhalos’ interactions drive the thermal
evolution in the presence of mass loss and structural changes due to the tidal
forces of the host halo[103, 127].

It has been proposed that the too-big-to-fail problem can be resolved by cored
halos in self interacting DM produced by cross sections of around
$1\rm\leavevmode\nobreak\ cm^{2}/g$ [74, 473, 117], i.e, for cross sections
lower than the threshold needed to trigger the gravothermal collapse phase.
Refs. [474, 475, 51, 131] have shown that this is unlikely given the large
differences in the central densities of the dwarf galaxies. A major issue is
that the least luminous satellites lie in dense halos, even as the more
massive, classical satellite galaxies are low in central density [453, 476,
151]. Cross sections preferred by the more massive satellites dramatically
underpredict the densities of the ultrafaint dwarfs [151]. Perhaps, the best
match in terms of scatter to the MW satellite densities is achieved with a
cross section of $10\rm\leavevmode\nobreak\ cm^{2}/g$ [74] or even larger once
the gravothermal collapse mechanism is considered in velocity-dependent SIDM
models [51, 471, 470]. In this context, it is important to keep in mind that
the core collapse process only needs to increase the central densities by a
factor of few typically to achieve consistency with measured densities. For a
large enough cross sections, much of the evolution of the subhalo could be
spent in the core collapse phase with central densities of between few to ten
times times $\rho_{s}$ [103].

Figure 16:  Circular velocity profiles for the ten subhalos with the largest
peak velocity across time in the simulations with the disk and bulge
potentials (right) and the DM-only simulations (left). Note that these are not
necessarily the subhalos with the largest peak velocity today. The SIDM cross
section over mass is set to $1\rm cm^{2}/g$, which implies that there is no
possibility of the subhalos entering the core collapse phase. The data for the
classical MW satellites brighter than $2\times 10^{5}\,\mathrm{L}_{\odot}$
[452] are denoted with filled squares with point sizes proportional to the
logarithm of their stellar masses. The points and curves are color coded by
the inferred pericenter distances; for the MW satellites the pericenter
distances calculated in [477] using Gaia DR2 have been used. Gray symbols are
satellites whose inferred pericenter distances have errors of the same order
as their mean value, and are hence uninformative. The thin dashed (faded)
lines indicate unresolved regions in the simulations. The CDM halo has three
subhalos as dense as Draco, while the $\sigma/m_{\text{DM}}=1\,\rm cm^{2}/g$
halo does not have any. Note that for all models, with and without the disk,
there is a trend that subhalos with smaller pericenters (cyan lines) have
lower densities. The MW satellites do not seem to adhere to this trend. Figure
adapted from Ref. [277].

The discussion of the proposed solutions based on SIDM is incomplete without
including the impact of the MW disk, which could have more dramatic effects
for SIDM subhalos compared to CDM subhalos for two reasons. One, recent work
shows that the core collapse of SIDM subhalos driven by the gravothermal
catastrophe is accelerated

by mass loss induced by tidal forces [103, 276, 127]. Two, the SIDM subhalos
that do not undergo core collapse will have large cores that makes them more
prone to tidal effects leading to further core expansion [478, 479, 153]. The
combination of these two effects leads to an enhanced diversity of the inner
DM distribution of SIDM subhalos for constant cross section SIDM models [276,
127]. Notice that velocity-dependent SIDM models with large cross sections at
the scale of dwarfs without the presence of the tidal effects of the MW disk
already produce a diverse subhalo population [51]. It would then be expected
that such environmental effects will enhance even further the diversity
observed in such models. Recently, Ref.[471] modeled the combined effects of
the gravothermal collapse, with the tidal interactions due to the MW disk
using orbital information from the MW satellites to find that a velocity-
dependent model similar to the one in Ref. [51] is able to fit the diversity
of the MW satellites. The presence of evaporation can delay or halt core
collapse [127], which provides an independent handle on the velocity
dependence of cross section since evaporation is sensitive to the cross
section. For sharp velocity dependence that imply that the cross section is
$\lesssim 1\ \rm cm^{2}/g$ at velocities of $200\ \rm km/s$, evaporation is a
small effect [472].

The too-big-to-fail problem relates to the mean DM density within galaxies’
half-light radii. Complementary information can be obtained by measuring the
density within other radii to get a better sense of the full DM profile. This
is more challenging although a few different methods have been proposed:

  * •

The method of Ref. [480] to use multiple stellar populations with different
half-light radii to infer the dynamical mass within different radii and hence
infer slope of the DM density profile on those scales. Cored DM profiles are
better fits than NFW to both Sculptor and Fornax according to some studies
[481, 482, 480, 475, 476], but other studies argue the situation is not clear
[483, 484]. Careful mass modeling of classical dwarf galaxies indicate that
the inner slopes are heterogeneous, with some dwarfs preferentially in cusped
halos [475, 476, 485]

  * •

Using the projected separation of globular clusters (GCs) in Fornax and
arguing that dynamical friction would have sunk these GCs to the center of the
dSph if the DM subhalo density profile were cuspy [486, 487]. This leads to
the inference of a large ($\sim$ kpc) core in Fornax [488]; although this
issue remains controversial and recent studies have found consistency of a
cuspy Fornax halo with its GC population, e.g. see Ref.[489, 490].

  * •

The existence of cold substructure in the dwarf spheroidals is another cited
evidence for a cored halo due to the reduced tidal effects as compared to a
cuspy central density [491].

Figure 17:  Two potential ways to disentangle CDM + supernova feedback from
SIDM as mechanisms of cusp-core transformations in dwarf galaxies. The upper
panels from Ref. [492] show the slope $\alpha$ (within $0.5-1\%$ of the virial
radius) of the halo density profile as a function of stellar mass for a full
baryonic physics implementation (solid symbols; FIRE implementation, see Ref.
[493]) and the corresponding DM-only simulations (open symbols). CDM in the
upper left and SIDM ($\sigma/m_{\text{DM}}=1$ cm${}^{2}/$g) in the upper
right. Note that this assumes that all the subhalos are in the core expansion
phase, which seems to be disfavored currently (see text). For SIDM models
where some of the subhalos are in the core collapse phase, we would expect
greater diversity in the slopes, possibly spanning the whole range from slopes
from 0 to -2.  For small stellar-to-halo mass ratios, the supernova feedback
mechanism is (energetically) less efficient in producing sizeable cores, while
the SIDM mechanism (constant cross section $\sigma/m_{\text{DM}}\lesssim 1$
cm${}^{2}/$g) remains efficient at all masses. The bottom panels are from the
isolated controlled DM-only simulations of a dwarf-size halo from Ref. [494].
The halos in both panels result in a $\sim 0.5$ kpc core after $t\sim 4$ Gyr,
with the simulation on the right having a toy supernova feedback model, while
the one on the right has a standard SIDM implementation
($\sigma/m_{\text{DM}}=2$cm${}^{2}/$g). The plots show the radial 2D phase
space density of 2000 tracers originally set as an orbital family (narrow
energy and angular momentum distribution) at the end of the simulation.
Supernova feedback is an impulsive mechanism that irreversible perturbs the
orbital family creating distinct shells in phase space, while the SIDM case is
an adiabatic mechanism, which leaves the orbital family of traces largely
unperturbed.

Given the current situation, a key avenue of progress is to identify
observables in dwarf galaxies that are truly distinct between CDM and SIDM
(see also Sec. 3.8). A few promising examples can be mentioned:

  * •

A relevant observational trend reported recently is that of the correlation
between the orbital pericenter distances and the central density measured for
the MW dSphs [131, 454]. The central densities are larger for dSphs that have
closer pericenters, which makes it difficult to see how the MW disk is solely
responsible for the diversity in central densities observed across the dSphs.
On the other hand, this could fit well within the predictions of SIDM models
if the high density dSphs (such as Draco and Ursa Minor) were in the core
collapse phase accelerated by tidal effects, while low density dSphs (such as
Fornax and Sextans) were in the usual core expansion phase [495, 276].

  * •

The efficiency of supernova feedback as a cusp-core transformation mechanism
is strongly dependent of the energy released relative to the change of
potential energy in the halo from being cuspy to cored [496, 497, 498, 499].
The efficiency peaks at around $M_{h}\sim 10^{11}$M⊙, just around the upper
mass end of dwarf galaxies, and drops sharply towards higher and lower masses
[500, 501]. For smaller dwarf galaxies, the stellar-to-halo mass ratio is
expected to be too low to be efficient at forming large cores [492] (although
see Ref.[502]). On the other hand, in the SIDM scenario with
$\sigma/m_{\text{DM}}\sim 1$ cm${}^{2}/$g, cores of a size comparable to the
scale radius of the halo are expected across all masses [106] in the core
expansion phase. In the core collapse phase, we would expect much greater
diversity in the values of the slopes as discussed in this section. For the
lowest concentration halos or halos with larger $V_{\rm max}$ where the cross
section is small enough, we could have cored halos. For smaller mass halos or
for halos with high concentrations or more tidally-stripped halos, one would
expect slopes between 0 and -2 depending on the size of the core and the
radius at which we are able to measure the slopes. Unambiguously identifying
and characterizing sizeable cores for stellar masses $\lesssim 10^{6}$M⊙ is a
promising avenue to constrain physical mechanisms of core formation (see top
panel of Fig. 17).

  * •

Cores formed due to DM self-interactions and due to supernova feedback have
been shown to be widely degenerate (when it comes to the shape of the density
profile, see Ref.[503]). However, the supernova-driven gas blowouts that are
responsible for expanding the orbits of DM particles to create cores need to
be impulsive (relative to the orbital timescales), while self-interactions are
effectively an adiabatic process of energy redistribution leading to the
formation of cores [256, 494, 499]. This intrinsically different nature leads
to distinct signatures in the orbits of stars: shells in phase space are
formed in the impulsive SN-feedback scenario, while they are not present in
the SIDM case [494]; see bottom panel of Fig. 17. The presence (or absence) of
such shells in the stellar distribution of dwarf galaxies with confirmed cores
would strengthen (weaken) the supernova feedback scenario as the main
mechanism of core formation.

###  4.7 Rotation curves of spiral galaxies

Problems associated with explaining the DM halo density profiles persist in
higher mass galaxies, where they manifest themselves in several different but
related forms. Unlike the situation with dSphs where it has proven difficult
to measure unambiguously the presence of constant density cores, the situation
in rotation-supported galaxies is clearer. The evidence for cores in many
spiral galaxies is strong, particularly those of low surface brightness (e.g.
[41, 504, 505, 506, 468]). This evidence is sufficiently ubiquitous to lead
some to suggest a universal cored density profile with a well-defined relation
between core density and size [507, 508, 509]. NFW fits to the rotation curves
of individual galaxies of a wide range of type often perform poorly or require
concentrations significantly outside the $\Lambda$CDM mass–concentration
relation scatter (e.g. [510, 511]).

Fits to individual galaxies can suffer from systematic issues related to
inclination angle estimation, non-circular motion and presence of bars [512,
513, 514, 515]. The results discussed in the previous paragraph can be
extended in a statistical sense across the galaxy population by employing
dynamical scaling relations between the masses, sizes and characteristic
velocities of galaxies. Empirical methods for assigning galaxies to the
$\Lambda$CDM halo population from $N$-body simulations (e.g. subhalo abundance
matching [516, 517]) generate a larger ratio of dark to luminous matter than
observed across the rotation curves of galaxies. This is manifest in the fact
that the $\Lambda$CDM halo mass function with NFW (or similar) profile
normalises both the Tully–Fisher and Faber–Jackson relations too high [518,
519, 520], and overestimates velocities in the low-acceleration outer regions
of galaxies [521]. This may be seen as extending the Too Big To Fail problem
[463] to higher masses. There is evidence that the overprediction of dynamical
mass in $\Lambda$CDM is in fact more general than implied solely by these
relations [522]. On the other hand, relaxing assumptions about the density
profile of the dark-matter halos, allowing for cored halos, leads to better
consistency between the halo masses preferred by rotation curves and those
implied by abundance matching [523, 524, 468].

Several other aspects of galaxy phenomenology have proven challenging to
$\Lambda$CDM and may hold clues to the nature of DM. Scaling relations like
the baryonic Tully–Fisher [525] and mass discrepancy–acceleration [526, 527]
(or radial acceleration; [528]) relations and Fundamental Plane appear to have
smaller scatter than the $\Lambda$CDM halo mass–concentration relation and
galaxy–halo connection (e.g. $M_{*}-M_{\text{vir}}$), from which they should
receive contributions [529, 530, 521]. The diversity of rotation curve shapes
has also been seen as problematic for CDM and may point towards DM
interactions [45, 48, 46] or strong feedback [271]. However, recent work has
argued that models with cuspy dark matter halos (i.e., $\Lambda$CDM halos
without strong feedback) are likely to be more consistent with the observed
rotation curves if the gas disks in late-type spirals are like those simulated
in EAGLE simulations [274, 275]. The origin of this conclusion lies in the
discrepancy between recovered and actual rotation curves caused by systematic
effects related to non-circular motion, thickness of gas disk and departures
from dynamical equilibrium in a systematic way that makes cuspy halos more
likely to look like cored halos. If this conclusion holds true, this may be an
indication that core collapse due to self interactions is more common in dwarf
halos than discussed in the SIDM investigations thus far.

Other observations that are potentially surprising in the context of
$\Lambda$CDM include spatially and kinematically correlated planes of
satellites around the MW and nearby galaxies [531], a dearth of DM in the
central regions of high-mass ellipticals [532], and the strong point-by-point
correlation of kinematics and baryonic surface density in spirals [300, 528].

DM self-interactions may offer a solution to some of these problems. To
explain the smallest core densities in the most extreme LSBs, the cross
section needs to be $\gtrsim 3\text{cm}^{2}/\text{g}$ [50]. Such large cross
sections are disfavored from cluster measurements as discussed in Sec. 4.1.1.
These arguments therefore directly motivate the velocity-dependent cross
section model, which we have argued is generic in the context of dark sector
models. The velocity dependence allows for large effects on galaxy scales
while maintaining consistency with constraints from cluster lensing studies.
Small cross sections (well below $1\rm cm^{2}/g$) may also work but will
likely not be distinguishable from CDM.

The mass discrepancy–acceleration (or radial acceleration) relation has also
been explicitly studied in the context of SIDM [50] (see Fig. 18); the good
agreement with observations provides support for the redistribution of DM mass
effected by self-interactions.

Figure 18: The radial acceleration relation from SIDM fits to low-redshift
galaxies with $\sigma/m_{\text{DM}}=3$ cm2/g. The black solid and red dashed
lines show an empirical fitting function (and its $1\sigma$ scatter) that
describes the observational data well. The inset shows residuals with the fit
function removed. From [50].

Solutions based on baryonic feedback, modified gravity and other DM physics
are also possible however (see Sec. 3.8).  High-redshift studies of galaxy
kinematics (especially at low mass) may provide a way to disentangle SIDM from
the effects of baryons by measuring the time scale of core formation. This is
likely to differ between the different models, although generic predictions
for SIDM are lacking.

Currently, the majority of detailed studies of kinematics in low mass galaxies
are limited to the Local Volume, however there are a few cases where
spectroscopic measurements of strongly lensed galaxies up to z$\simeq$2 have
been made (e.g. [533, 534, 535, 536]). It is also possible to measure the
Tully–Fisher relation as a function of redshift ($z\lesssim 2$) to study the
relative evolution of galaxies’ stellar and DM components. Evidence for
variation in the Tully–Fisher normalisation is however controversial, with
some authors finding no significant trend (e.g. [537, 538, 539]) and other a
trend towards relatively more DM at high $z$ (e.g. [540, 541, 542, 543]) The
uncertainties in each of these methods are large due to the difficulty of
reliably measuring the kinematic tracers and small sample sizes but they hold
promise for testing the physics being implemented in hydrodynamic simulations.

Advancing this field will require a synergy of observational and modelling
techniques. On the one hand, empirical knowledge will benefit both from the
compilation of homogeneous datasets with high quality photometry and
kinematics (e.g. [544]), and from the increase in statistics and precision
afforded by current and next-generation IFU studies such as MaNGA [545].
Problems primarily based on photometry will benefit from DES, LSST, Euclid and
Nancy Grace Roman Space Telescope [546]111111https://roman.gsfc.nasa.gov/
data. Large future datasets will enable stringent cuts on quantities like
inclination and the degree of pressure support, reducing their confounding
effects. On the other hand, advances in modelling will sharpen the likelihood
function for dynamical observables for given theoretical assumptions. This
requires not only further modelling of a range of systems in SIDM, but also
increased resolution in $\Lambda$CDM hydrodynamic simulations—to push more of
the relevant physics above the grid scale and hence reduce uncertainty in the
effect of feedback—and an increase in the sophistication with which models and
data are compared.

###  4.8 Ultra-diffuse galaxies

Long suspected of dominating the low-luminosity galaxy population [547], LSBs
with central surface brightness fainter than $\mu=23\leavevmode\nobreak\
\hbox{mag}/\Box^{\prime\prime}$ are experiencing a second wave of discovery
after the initial discoveries in the 1990’s enabled by CCDs [548, 549, 550,
551, 552]. With new telescope and survey designs, LSBs with surface
brightnesses as faint as $\mu=31-32\leavevmode\nobreak\
\hbox{mag}/\Box^{\prime\prime}$ can be detected either as resolved stellar
populations [553] or in diffuse light [554]. This population is heterogeneous.
Many are “ordinary" dwarf galaxies, which are low surface brightness by virtue
of the size-luminosity relation [555, 556, 557]. However, a number appear to
be

bigger outliers in either the size-luminosity relation or in the stellar
mass–halo mass relation.

Ultra-diffuse galaxies (UDGs) are an outlier population of LSBs with large
sizes of $r_{e}>1.5$ kpc. Since 47 of these galaxies were found by the
Dragonfly telescope in the Coma cluster [558], large populations have been
observed in well-studied clusters, including Coma, Virgo, and Fornax [559,
560, 561, 116, 559, 562]. In recent years, UDGs have been observed in a
diversity of environments, including the field [563, 550, 564, 551, 552, 565].
There are differences in the properties of UDGs in field and cluster
environments–field galaxies tend to be gas-rich, blue, and star-forming, while
those in clusters tend to be redder and quenched [550, 561, 562, 566, 567].
Moreover, a number of cluster UDGs have unusual properties, such as being
surrounded by anomalously many globular clusters given their stellar mass, and
inhabiting halos far larger than indicated by the mean stellar mass–halo mass
relation [568, 569, 570, 571, 572, 573, 574]. Because of these differences in
field and cluster environments, and because of the odd properties of cluster
UDGs, an important question to answer is if these galaxies are true outliers
in galaxy formation, or if they are part of a continuum of LSBs. It is
important to understand if there are different formation mechanisms for the
field and cluster populations [575], and how those mechanisms arise in the
context of CDM and SIDM galaxy formation models.

In the context of CDM simulations, separate formation pathways are identified
for field and cluster UDGs, and SIDM may additionally play a role in both. In
the field, although high spin was originally posited as the driver of UDG
formation [576], the current conventional wisdom is that flybys and minor
mergers drive the large size of $M_{*}\sim 10^{8}M_{\odot}$ UDGs [577, 578].
Vigorous outflows can play a role too [579], and enhance the effect of flybys
and mergers. However, the presence of a large reservoir of gas in the field
UDGs makes the fly-by and concomitant tidal stripping unlikely explanations
for field UDGs [567, 580]. In addition, the rotation curves of these field
UDGs present challenges for all models of galaxy formation (whether based on
CDM or SIDM) making this one of the most interesting avenues to pursue in the
coming years [581]. Self-interactions of DM may play a role in UDG formation,
by enhancing the distribution of galaxy sizes and rotation curve shapes at
fixed stellar mass, as hinted at for the broader class of LSBs [48, 50], but
more work is clearly needed to flesh out the predictions of SIDM models for
these outliers.

Cluster UDGs are shaped by their environment. CDM simulations show that
although some cluster UDGs were previously UDGs in the field, many became UDGs
after infall [163, 582, 583, 584]. These simulations show that “typical” dwarf
galaxies can become UDGs if they fall into clusters early and experience
quenching (by ram-pressure stripping), fading of light as the stellar
population ages [575], and tide-driven expansion of the stellar population.
Intriguingly, Refs. [585, 163] argued that galaxies in DM halos with large
cores more easily become UDGs. Idealized simulation including DM self-
interactions in Ref. [586] have shown that it is possible to reproduce the
properties of NGC1052-DF2 and DF4 in the SIDM models.

The stellar component in dark-matter-dominated galaxies expands faster in
cored halos rather than cusped halos in response to tides, which has been
shown in the context of CDM [155, 587] as well as SIDM [75]. In the semi-
analytic model of Ref. [163], cored halos play the essential role in
populating the tail encompassing large UDG sizes, as shown in Fig. 19. This
model with cored halos at infall was able to reproduce the observed sizes and
stellar masses of UDGs, as well as their abundance as a function of host halo
mass. However, the assertion that cored halos are necessary has been
questioned in Ref. [588] who argue that it is possible for a large number of
UDGs to form in the field as LSBs and then be accreted into the cluster. In
this scenario, there are still be some UDGs created close to the center of the
cluster due to extreme tidal interactions in cuspy halos. More work on
discriminating the two scenarios is warranted, as are cosmological CDM and
SIDM simulations tailored to this problem.

Figure 19: Observed and simulated stellar masses and half light radii. From
[163].

###  4.9 Warping of stellar disks

We described in Sec. 3.6 how the drag force experienced in SIDM by halos
falling through a medium of DM can cause a separation between the centroids of
stars and DM in a galaxy. Both numerical [589] and simple analytic [590]
models show that the potential gradient that this establishes across the
stellar disk induces a U-shaped warp in the disk facing the direction of
infall (Fig. 20, left), and a longer-lasting disk thickening. While S-shaped
warps formed by tidal distortions of the stellar light profile are abundantly
observed in cluster environments, indicating that they are readily produced by
physics unrelated to DM self-interactions, such physics is unlikely to
generate prominent U-shaped warps. Such warps are only formed by a
differential force on the disk and its halo, e.g. due to SIDM drag.

The simulations of Ref. [589] indicate that the warping effect may be
observable on a galaxy-by-galaxy basis in current and next-generation
photometric data when SIDM cross sections are $\gtrsim
0.5-1\,\mathrm{cm^{2}\,g^{-1}}$. This has inspired work aimed at constraining
SIDM by means of this signal [590]. Here, the $r$-band images of $\sim$3000
edge-on disk field galaxies with $M_{*}>10^{9}M_{\odot}$ and $D<250$ Mpc from
the Nasa Sloan Atlas are reduced to yield estimates of their degree of
U-shaped warping, and a forward-model is developed to predict this in SIDM as
a function of interaction cross section and range, halo properties, and DM
environment. This is then combined with a simple model for astrophysical and
measurement noise to create a galaxy-by-galaxy likelihood for the warp
statistic, and the noise and other nuisance parameters are marginalised over
to obtain constraints on the SIDM cross section. Due to the requirement that
the system satisfies the fluid approximation and that drag rather than
evaporation is the dominant result of the self-interactions, only long-range
(light mediator) models can be constrained in this way, although a sample of
galaxies in higher-density multi-streaming environments (e.g. clusters) may in
the future be used to constrain shorter-range forces.

No correlation of the warp strength or direction with the SIDM expectation is
found in the analysis, leading to the constraints shown in Fig. 20 (right) on
the momentum transfer cross section at $300$ km/s. These are plotted as a
function of the median relative velocity of the galaxies and the background
medium because this is considerably uncertain. The two models shown correspond
to a) a constant velocity over the entire sample and b) velocities
proportional to those from the CosmicFlows-3 catalogue, respectively. At face
value these bounds are an order of magnitude tighter than those from other
SIDM effects, e.g. dwarf galaxy evaporation. The model does however make a
number of simplifying assumptions to render the inference tractable, for
example that the halos are fully spherical, thermalised and equilibrated, and
that all parts of the halo move together. Simulations of halos in environments
similar to those of galaxies with measurable morphology will be needed to
validate and refine this novel method for testing SIDM.

Figure 20: _Left:_ Cartoon of the formation of a U-shaped stellar disk warp in
SIDM. _Right:_ Constraints on the momentum transfer cross section at 300 km/s
for long-range (Rutherford) scattering from 3,213 measured warps in the nearby
Universe, as a function of the median galaxy velocity relative to the
background. Adapted from Ref. [590].

###  4.10 Large-scale structure constraints

As we have discussed in Sec. 3.7, there are several particle physics models
that can modify the power spectrum of density fluctuations, thereby leaving
signatures on the temperature and polarization of the CMB, the clustering of
galaxies and the halo mass function. A necessary feature in these models is
the presence of light mediators. As discussed in Section 3.7 dark matter
models that tightly couple to dark radiation in the early universe can lead to
dark acoustic oscillations and diffusion damping in dark matter (like the the
BAO and Silk damping), modifying the linear power spectrum. In Ref. [244],
cosmological data from the lensing of the CMB, the clustering of galaxies on
linear scale and the BAO were analyzed to derive constraints on the fraction
of dark matter that can couple to dark radiation and the strength of the
coupling. Their results show that a fraction of partially interacting dark
matter, $f_{\rm int}>5\%$ is in severe tension with current observations of
large scale structure.

The damping scale of the power spectrum due to interactions between dark
matter and radiation is set by the horizon at the epoch of kinetic decoupling.
For a kinetic decoupling temperature, $T_{\rm kd}$, the suppression in the
power spectrum occurs approximately for modes $k>10(T_{\rm kd}/\mathrm{keV})\
\rm Mpc^{-1}$ [591, 592, 593, 8, 192]. SIDM coupling to a dark radiation
component like dark photons or sterile neutrinos can significantly delay the
decoupling of dark matter from the thermal bath leading to damping on scales
of order $0.1$ Mpc/h or smaller. An effective framework to describe these
interactions and the associated phenomenology has been developed by the ETHOS
collaboration [594, 179]. One of the central predictions is a potentially
observable cut-off in the halo mass function [233, 234, 235, 251, 245, 5]
arising from the damping of the power spectrum.

One of the ways to probe these scales is through the absorption spectrum of
Lyman-$\alpha$ forest, which measures the density of neutral hydrogen on large
scales. In Ref. [5] the Lyman-$\alpha$ forest power spectrum measurements was
used to put strong constraints on the damping scale (by recasting the warm
dark matter constraints) to set an upper limit on the decoupling temperature
of $T_{\rm kd}>1$ keV and a minimum halo mass of $M_{\rm
min}<10^{8}M_{\odot}/h$. The origin of this constraint is rooted firmly in the
physics of the dark sector, including the production mechanism for the dark
matter. The massive gauge boson that mediates the self-interaction in simple
SIDM models has to decay in order to avoid over-closing the universe and in
Ref. [5] decay of the gauge boson to light DM particles were assumed. Decays
to Standard Model particles have much stronger direct detection constraints
[595] and indirect search constraints [4]. In the future, we expect the
combined measurements of the Lyman-$\alpha$ power spectrum and the halo mass
function from faint galaxies will strongly constrain the space of allowed
parameter space for the models with dark radiation. Prospects for the future
in this direction are bright, with CMB-S4 [596]121212https://cmb-s4.org/
combined with LSST likely to improve constraints on $\Delta N_{\rm eff}$ in
these models by an order of magnitude [231]. Another promising avenue to
constrain simple SIDM models and those with additional light degrees of
freedom and interactions is to probe the linear and quasi-linear power
spectrum to higher $k$ values using 21-cm data. Future measurements offer the
possibility to constrain the power spectrum amplitude at the 10% level
integrated over $k=40-80\ \rm Mpc^{-1}$ [597], which can probe decoupling
temperatures well above a keV.

###  4.11 Summary of constraints on the self-interaction cross section

The various observational probes for SIDM can be thought of as a mapping
between the velocity dependence of the cross section and observable processes
at various energy scales within halos of different masses during hierarchical
structure formation. On the galaxy cluster scales, where the internal velocity
dispersions are of the order of a few 1000s of km/s, the most stringent
constraints currently come from measured core densities from strong lensing
[55]. The small sizes of the cores and the high central densities in these
clusters constrain cross section to
$\sigma/m_{\text{DM}}<0.13{\leavevmode\nobreak\ \rm cm}^{2}$/g at $95\%$
confidence in one study [353] and
$\sigma/m_{\text{DM}}<0.35{\leavevmode\nobreak\ \rm cm}^{2}$/g at $95\%$
confidence in another [354], as discussed in more detail in Sec. 4.1.1. Even
though there is some evidence of a non-zero cross section in these analyses,
it is not possible to conclude that we have a detection. In the future,
detailed simulation work combined with joint analyses of strong lensing, weak
lensing, X-ray and velocity dispersion data for a large number of clusters
could provide a concrete path to measuring a non-zero cross section at
relative velocities of about $1500\ \rm cm^{2}/g$. Similar analysis at the
group scale gives an upper limit of $1.1\ \rm cm^{2}/g$, with a mild
preference for a non-zero cross section ($0.5\pm 0.2\ \rm cm^{2}/g$; [354]). A
variety of other probes like the displacement between galaxies and dark matter
in merging systems, the steepening of the profiles beyond the dark matter core
and the existence of displaced BCGs will further help establish consistency
between different observations once the measurement systematics and the
interplay between baryonic and dark matter evolution are better constrained
and understood. Among these, constraining BCG wobbles, which are correlated
with the size of the constant density core, using offset BCGs seems highly
promising with current constraints at $\sigma/m_{\text{DM}}<0.39\ \rm
cm^{2}/g$ at the 95% C.L. [428].

At lower velocities, a wider range of cross sections is in principle still
allowed as long as the model is consistent with the constraints at high
velocities. As baryonic effects are more relevant for lower mass halos, the
combination of baryons and SIDM leads to a wide range of halo mass profiles
inside a region of about the NFW scale radius. In particular, for galaxies
like the MW and large elliptical galaxies, the large concentration of baryons
in the inner regions makes the density profiles of SIDM and CDM models
indistinguishable [119, 277, 598, 599]. We can expect constraints from
evaporation of DM in satellites of the MW (and similar galaxies) due to
scattering events between subhalo and host halo DM particles when the cross
section at $200\leavevmode\nobreak\ \rm km/s$ exceeds $5--10\rm cm^{2}/g$
(e.g., [153, 276, 470]) but the constraints have not been robustly inferred.
At lower velocities, fitting the rotation curves in SPARC database has shown
that $\sigma/m_{\text{DM}}>3\ \rm cm^{2}/g$ at about $100\rm km/s$ is needed
to adequately explain the low-surface brightness galaxies with the lowest
central DM densities [50]. Note that the predictions for the cored profiles of
SIDM halos do not change dramatically as the cross section is increased, so
cross sections even an order of magnitude larger cannot be easily
distinguished based on galaxy rotation curve data alone [117]. If the cross
sections are $\gtrsim 100\ \rm cm^{2}/g$, then core collapse will set in for
median concentrations and it may become difficult to explain the low density
DM cores of low-surface brightness galaxies.

Within the local universe, a self-consistent explanation for the distribution
and inner structure of the satellites of the Milky way is not yet fully known.
While the “missing satellites problem” is largely thought to be resolved given
the current census of satellites, the exact nature of the galaxy–subhalo
connection is still unclear and currently requires populating galaxies in
subhalos with $V_{\rm peak}\sim\leavevmode\nobreak\ 7$ km/s that are below the
collisional cooling limit. Recent work has shown that the minimum $V_{\rm
peak}$ will have to be smaller for moderate cross sections $\sigma/\rm{m}_{\rm
dm}\lesssim 5\ \rm cm^{2}/g$ [472]. While the “Too Big To Fail” problem is
resolved by moderate SIDM cross sections of $\sim$1 cm${}^{2}/$g [600, 277],
this underpredicts the central densities of ultra-faint satellites [601, 51,
276]. These conclusions are highly sensitive to the presence of the MW disk,
but recent work including the disk potential has demonstrated the validity of
these arguments [472]. The issue is not restricted to the ultra-faints
dwarfs—using analytic methods it has been found that if all the bright Milky
Way satellites are in the core expansion phase, it is difficult to fit their
measured stellar dispersions with a single cross section [47]. A velocity-
dependent cross section with $\sigma/m_{\text{DM}}>10\ \rm cm^{2}/g$ at
velocities relevant for the internal dynamics of satellites ($5-20\ \rm km/s$)
that drives the faint satellites to core-collapse is perhaps the best bet to
explain these observations [103, 276, 51, 471, 470, 472].

Putting together these constraints suggests that a significant velocity
dependence is necessary to explain data constraining DM mass profiles from
dwarf galaxy to clusters of galaxies. Compared to earlier work [55], the
motivations for a strong velocity dependence have increased over time. An
exciting possibility within models that feature large cross sections at low
velocities is that many of the halos and subhalos could be in the core
collapse phase of their evolution, with predictions for observables that are
different from CDM as discussed in this review. Examples of such sharp
velocity dependence are discussed in Refs. [471, 130, 470], where analytic
methods and N-body simulations are used to predict the central densities of
Milky Way satellites in the core expansion and core collapse phases. More
comprehensive analytic and numerical treatments including the effect of the MW
disk are required to make further progress. However, it seems clear that the
satellites of the MW can provide a definite test of the SIDM models that aim
to explain DM density profiles across the mass range from dwarf galaxies to
galaxy clusters.

##  5 Beyond the simple SIDM framework

In the basic setup of SIDM, a single species of DM particles interact through
the exchange of a light (or massless) scalar or vector mediator with a cross
sections per unit mass of $\mathcal{O}(1\leavevmode\nobreak\
\mathrm{cm}^{2}/\mathrm{g})$, as discussed in Sec. 2. We have, however,
largely neglected the implications of a hidden sector beyond the phenomenology
of the DM particle itself. There are also more complex hidden sectors that
predict interesting phenomenological signatures, not only in astrophysical
systems, but in cosmology and potentially in laboratory experiments as well.
While there are model-dependent constraints on any specific particle physics
realization of SIDM, many models may produce similar phenomenological
signatures and are thus subject to similar restrictions from data.

In this section, we review the consequences of non-minimal SIDM models and of
the hidden sectors themselves. We first discuss two broad generalizations of
the simple SIDM scenario that we have assumed thus far. In Sec. 5.1, we relax
the assumption that self interactions are elastic and allow a mechanism for
SIDM to dissipate energy. In Sec. 5.2, we discuss the possibility that SIDM as
a subcomponent of the overall DM content, allowing it to have very strong self
interactions. The hidden sector contains light mediators that allow DM
particles to self interact, and we discuss the consequences of the light
mediator in Sec. 5.3. Finally, in Sec. 5.4, we mention various complementary
constraints on DM sectors that interact with the SM.

###  5.1 SIDM with dissipation

If SIDM is not a simple single-state particle, it is possible for SIDM to
undergo inelastic processes between its ground state and an excited state. The
presence of exicited states can arise in models of composite SIDM, such as
atomic-like DM [192, 86] and nuclear-like DM [85], in which fundamental dark-
sector particles are bound under a dark-sector force. The composite nature of
DM permits a spectrum of internal energy levels, along with a naturally large
self-interaction cross section. Excited states also arise if DM is comprised
of different fundamental particles with a mass splitting and an off-diagonal
interaction [602, 603, 604, 605]. In either case, collisions may induce a
transition to the excited state, which may then relax back to the ground state
through collisional de-excitation or decay via the emission of some form of
radiation. This radiation may be a massless force carrier, such as a dark
photon, or a light mediating particle, with a mass smaller than the energy
splitting between the ground and excited states.

Assuming there is sufficient kinetic energy available, collisions within DM
halos can excite DM particles. If the decay width of the excited state is
large and if the halo is optically thin to the emitted radiation, the halo
loses energy and could eventually collapse [86]. As discussed in Sec. 3.2, the
gravothermal collapse process with elastic SIDM interactions can be
accelerated by environmental effects [103]; for an isolated halo, the collapse
process can be accelerated through energy dissipation [118], assuming DM
particles are not ejected from the halo [606]. The observation of uncollapsed
dwarf galaxies and LSBs can place restrictions on the acceptable amount of
halo cooling in those systems [118]. On the other hand, the gravothermal
collapse of DM halos with dissipation may be able to form the seeds for the
observed population of supermassive black holes [607].

###  5.2 Subcomponent SIDM

Some models of DM postulate that while most of the DM in the universe may be
(approximately) collisionless, there may be a small component that has very
large self-interactions, including dissipative interactions. This small
component could have a large impact on the observational predictions. For
example, even a very small amount of ultra-strong SIDM can cause a halo to
experience gravothermal collapse, as discussed in Sec. 3.2, and form seed
black holes [85, 102]. Alternatively, some models also propose that that up to
15% of the DM of the Universe exhibit a behavior similar to baryons (i.e.,
experience energy dissipation), including disk collapse [195, 194] and plasma
instabilities [608, 609, 610]. These effects have for example been studied in
the context of the formation of a dark core in major mergers such as Abell 520
[424, 611]. Dissipation could also lead to qualitatively new predictions for
the small-scale structure of DM, including asymmetric dark stars [612].

###  5.3 Phenomenology of light mediators

In many models of SIDM, the self interaction arises from the exchange of a
light mediator, such as a 1–100 MeV scalar or vector particle with a Yukawa
coupling to DM [81]. The mediating particles themselves can impact a wide
range of observations, allowing for additional probes of SIDM hidden sectors.

The production mechanism of SIDM in the early Universe can play an important
role throughout cosmic evolution (see, e.g., Ref. [613]). A common assumption
is that DM particles were at some point in thermal equilibrium with the bath
of CMB photons, through a common reheating mechanism. In a completely isolated
hidden sector, however, DM may be produced in its own thermal bath through a
separate reheating mechanism, resulting in isocurvature perturbations that are
constrained by CMB observations [614], and thus place constraints on self
interactions [615].

With a nontrivial dark sector, the SIDM relic abundance can be set through
mechanisms analogous to the standard thermal freeze-out scenario. For
instance, SIDM may undergo dark-sector freeze-out (in which SIDM particles
annihilate into lighter dark-sector particles) or the SIMP mechanism (in which
SIDM particles experience number-changing processes).

The mediator constitutes the thermal bath for DM and contributes an additional
$\Delta N_{\textrm{eff}}$ to the standard number of relativistic degrees of
freedom in the early Universe. Measurements of the CMB and of nuclei
abundances restrict the amount of $\Delta N_{\textrm{eff}}$ permitted at the
time of recombination and Big Bang Nucleosynthesis (BBN), respectively [1].
Moreover, for models of light DM (with masses in the MeV range), the DM itself
can contribute to $\Delta N_{\textrm{eff}}$ during BBN [616]. Since no
interaction with the SM is required in a minimal model of SIDM, the hidden-
sector bath can be slightly colder than the CMB and thus easily evade $\Delta
N_{\textrm{eff}}$ constraints [617, 85]. However, if the hidden sector does
communicate with the SM, $\Delta N_{\textrm{eff}}$ will set constraints on
such models.

The presence of a hidden thermal bath for SIDM also leads to the suppression
of small-scale fluctuations in the CMB, as described in Sec. 3 for the case of
dark photon mediators. Such signatures can be constrained by the CMB, as well
as large-scale structure observables, such as the Lyman-$\alpha$ forest [251,
5, 231] as discussed previously in Sec. 4.10. If models predict a strong
cutoff of the matter power spectrum, then additional powerful constraints at
small scales arise from the existence of the smallest DM halo observed [221].

###  5.4 Connections to the Standard Model

While not strictly required, there are good reasons to consider non-
gravitational interactions of SIDM with SM particles. These interactions allow
SIDM to be produced from the SM thermal bath, which can constrain the model
space stringently [4] and also provide simple viable models [613, 618, 619].
These interactions can impact the composition and cosmological evolution of
the Universe and may give rise to observational signatures in laboratory
experiments and astrophysical systems (see e.g., Refs. [620, 621, 622]).

While there are many different ways in which SIDM particles can couple to the
SM, most attention has been paid to two benchmark scenarios, in which the
couplings of SIDM particles are proportional to the electric charge and mass
of the various SM particles, respectively [623]. The former case arises from
the (kinetic) mixing between a new vector boson and the photon of the SM
[182], while the latter case arises from the mixing between a new scalar boson
and the Higgs boson of the SM [624]. The new exchange particle (often referred
to as the mediator) is typically considered to also be responsible for the
generation of strong self-interactions.

It may be possible to understand the cosmological abundance of SIDM particles
in terms of their annihilation and scattering rates [613]. If these rates are
large enough to bring SIDM particles into thermal equilibrium with the SM
plasma, the relic abundance can be set through the standard WIMP freeze-out
mechanism. Even if the interaction rates do not bring the two sectors into
equilibrium, the relic abundance may be set, for example, via the freeze-in
mechanism [3]. Moreover, many SIDM models require the presence of additional
(lighter) states in the dark sector, as discussed in Sec. 5.3. In order to
ensure that these additional particles do not constitute a large fraction of
DM or dark radiation, they typically need to be unstable, such that they decay
into SM particles. In order to evade constraints on $\Delta N_{\textrm{eff}}$
from BBN and the CMB and in order not to modify the predictions of BBN through
the photodisintegration of light elements, these decays should occur
sufficiently early in cosmological history, placing a lower bound on the
interaction strength between SIDM particles and SM states [625].

On the other hand, cosmological considerations also place an upper bound on
the interaction strength for DM annihilation into and DM scattering with SM
particles. DM annihilation can inject energy into the SM plasma, affecting the
optical depth of CMB photons post-recombination [626, 627, 628, 629, 630, 631,
4, 632, 633, 634]. CMB constraints on DM annihilation are strong for $s$-wave
annihilation, but significantly weaken for $p$-wave etc. due to the velocity
suppression of cold DM. Additionally, there are constraints on DM scattering
with ordinary matter, which can induce spectral distortions [215, 223] and a
suppression of anisotropies in the CMB [213, 208, 214, 216, 217, 218, 219,
220, 224, 635, 226] and the matter fluctuations [221, 57, 222, 225]. These
constraints are complementary to the more traditional search methods of direct
and indirect detection, described below.

Direct detection experiments place strong constraints on the interactions of
DM particles with the SM and are known to exclude many DM models that exploit
the freeze-out mechanism to explain the DM relic abundance.  Another concrete
way in which the direct detection of a SIDM particle can be motivated is the
requirement in some models that the mediators decay before the beginning of
BBN [636]. These experiments are particularly sensitive if the mediator of the
interaction is light, such that scattering with small momentum transfer is
enhanced [637, 638, 595]. In many models of SIDM, the DM particle mass is
below 1 GeV, for which the most promising strategy is to search for DM-
electron interactions [639].

A second source for observable signals of SIDM particles are annihilation
processes occurring in regions of enhanced DM density, such as MW satellites
or the Galactic Centre. The resulting fluxes of photons, neutrinos and/or
changed anti-particles can be searched for with indirect detection
experiments. These constraints are particularly strong if the annihilation
cross section in enhanced for small relative velocities (so-called Sommerfeld
enhancement [168, 640, 9, 641, 642, 10]). We note, however, that annihilation
signals may also be strongly suppressed or absent, for example if the
annihilation cross section is suppressed in the non-relativistic limit, or if
SIDM particles have a particle-antiparticle asymmetry similar to the one of
baryons [634].

Finally we note that SIDM particles may also be produced at accelerators,
although the resulting signatures are expected to look very similar to other
DM models, such as WIMPs. A notable exception occurs in the case of the SIDM
particles produced form a bound state (called WIMPonium or Darkonium), which
subsequently decays into SM particles [643, 644, 645, 646, 647, 7].

##  6 Future Prospects

###  6.1 Observational probes

Upcoming galaxy surveys like the LSST and Dark Energy Spectroscopic Instrument
(DESI; [648])131313http://desi.lbl.gov will allow us to observe the faint
satellite galaxies that inhabit the DM halo of the MW, allowing a follow-up
similar to the discovery of satellites and streams that DES, among other
surveys, has already made possible. Apart from observing the satellites of the
MW, there are several surveys that target the faint galaxies in the nearby
universe, among them the SAGA survey [437, 444] which searches for satellites
around MW analogs, ELVES [649] which provides a volume-limited sample of
classical satellites around MW analogs with confirmed distances, LBT-SONG
[650] which searches for satellites around intermediate-mass hosts beyond the
Local Group, and MADCASH [651] which observes Magellanic analogs in the Local
Volume. Apart from focussing on satellite systems, surveys like DELVE [652]—a
deep, multi-component survey that uses the Dark Energy Camera to study low-
surface brightness, DM dominated galaxies in the local universe—will observe
field dwarf galaxies. Similarly the Merian survey
[653]141414https://merian.sites.ucsc.edu/ uses the HSC to measure weak lensing
profiles of Small and Large Magellanic Cloud analogs in the local universe.
Spectroscopic follow-up of these galaxies will be enabled by Keck, the Thirty
Meter Telescope and proposed experiments like the Maunakea Spectroscopic
Explorer (MSE; [654]) allowing us to characterize the detailed properties of
these faint objects. Furthermore, surveys like the Widefield ASKAP L-band
Legacy All-sky Blind surveY (WALLABY) [655]—an HI survey in the local
universe—will illuminate the baryonic component of faint galaxies in the
nearby universe.

Within the next decade, the Rubin Observatory and Euclid are expected to
discover thousands of new galaxy-scale strong gravitational lenses [656]. In
order to exploit the future data to the maximum possible extent, new methods
have been developed to accelerate detections of such systems as well as
perturbers in the images [657, 658, 659, 660, 661, 662]. The vast increase in
sample sizes will provide much stronger constraints on DM scenarios. Similar
machine learning-based methods may also facilitate the analysis of other data
from upcoming surveys.

A dedicated follow-up effort for a well-chosen subsample of these lens systems
to obtain spectroscopic redshifts and (for Rubin-discovered systems) high-
resolution imaging will enable improved constraints on SIDM. On the one hand,
gravitational imaging analyses of extended arcs and rings are likely to lead
to the detection of several DM perturbers (either substructure or line-of-
sight halos) [360, 361, 362, 364]. Detailed modeling of the density profiles
of these perturbers [369, 370] will provide constraints on the SIDM cross
section at small velocities. On the other hands, a joint analysis of a large
population of quasar lenses [384, 378] could also be used to constrain the
SIDM cross section. In parallel, follow-up ALMA observations of lensed dusty
star-forming galaxies discovered in South Pole Telescope (SPT;
[663])151515https://pole.uchicago.edu/public/Home.html maps [664, 363] and
Very Long Baseline Interferometry (VLBI) observations [665, 666] of lensed
radio sources will also be used to discover small DM halos and subhalos. In
both cases, detailed studies of the density profile of these perturbers will
lead to constraints on SIDM. As we have discussed, the possibility of the
detection of a number of perturbers that are denser (or more compact) than
expected in CDM can be a clean path to the discovery of dark sector DM.
Finally, high-resolution imaging of galaxy-scale lenses could be used to
measure the small-scale power spectrum [392, 388, 379, 389, 390, 359]
independent of how of the density profile of line-of-sight halos and subhalos
are modeled, hence providing complementary constraints on SIDM.

The upcoming galaxy surveys across a wide range of halo masses, including both
satellite and field galaxies, will help constrain the cross section on
multiple velocity scales giving us a rich understanding of the varied
phenomenology that derives from the interplay of hierarchical structure
formation and the particle properties of DM. While observations discussed
above will help chart out the nearby Universe, we note also that data from the
James Webb Space Telescope [667]161616https://www.jwst.nasa.gov/ and SphereX
[668]171717https://www.jpl.nasa.gov/missions/spherex in combination with the
Rubin Observatory data will give us novel constraints on the nature of DM from
the first galaxies at high redshifts from the epoch of reionization.

In addition, constraints on SIDM cross sections from massive halos, cluster
and group profiles inferred via weak lensing [95] and core sizes [669, 354,
353] determined via strong lensing and stellar velocity are expected to
tighten significantly given the number of such massive clusters ($M\gtrsim
1\times 10^{14}M_{\odot}/h$) that will be detected by a combination of optical
surveys—such as Rubin and Euclid, as well as Y6 data from DES—and CMB surveys
such as the SPT, Atacama Cosmology Telescope (ACT;
[670])181818https://act.princeton.edu/, Simons Observatory
[671]191919https://simonsobservatory.org/ and CMB-S4. Given its sky coverage,
Rubin alone is expected to detect more that $50,000$ clusters in this mass
range, up from $\sim 5000$ clusters detected in the DES Y1 data. At slightly
lower mass scales ($5\times 10^{13}M_{\odot}/h\lesssim M\lesssim
10^{14}M_{\odot}/h$), the eROSITA mission
[672]202020https://www.mpe.mpg.de/eROSITA is expected to map out $\sim
100,000$ low-mass clusters and galaxy groups. On those parts of the sky where
eROSITA overlaps with optical surveys, weak lensing measurements around these
objects will improve constraints on SIDM cross sections at the relevant
velocity scale. Moreover, with an increasing sample of clusters and precise
weak lensing maps from telescopes like Roman and Euclid cross section
constraints from offsets between BCGs and potential centers of galaxy clusters
are also expected to improve by a factor of two [428]. Experiments like MSE
will allow complementary studies of positional offsets with velocity offsets.
A key advance in the future will be robustly combining weak lensing, strong
lensing, X-ray and stellar velocity dispersion data with a dedicated analysis
that includes a self-consistent halo model for SIDM. It is possible that we
may be able to measure cross sections of order $0.1\ \rm cm^{2}/g$ with the
next generation of observations and analysis tools. A necessary requirement
for this to become a reality is the availability of cosmological simulations
of galaxy clusters using cross sections of $0.1\ \rm cm^{2}/g$ and smaller. We
discuss simulations in more detail in the next subsection.

Finally, the upcoming surveys, several of which are optimized to answering
questions about the large-scale structure and evolution of the universe, will
provide unprecedented data sets to test models that predict modifications to
the power spectrum (Secs. 3.7, 4.10 and 5.3). Conservatively, CMB-S4 for
example is expected to measure the energy density in weakly-coupled light
particles, parametrized by $\Delta N_{\rm eff}$, to a precision of
$\sigma(N_{\rm eff})=0.02-0.03$ [596]. As we have discussed, $\Delta N_{\rm
eff}$ is a critical discriminant of SIDM models and dark sector physics in
general. These experiments will also be contemporaneous with galaxy surveys
from the Rubin Observatory, Euclid, DESI and the Roman Space Telescope, which
will measure cosmological parameters from the clustering of galaxies,
gravitational lensing and the abundance of clusters. Simons Observatory and
CMB-S4 will also improve CMB lensing measurements by an order of magnitude
compared to Planck [596]. Combined with measurements of Lyman-$\alpha$ forests
from spectrographs like Prime Focus Spectrograph on the Subaru telescope
[673], WEAVE on the William Herschel telescope [674], and in the future MOSAIC
on the Extremely Large Telescope [675], WFOS on TMT [676] and MSE [654], we
will have an opportunity to combine several different measurements to
constrain these models of dark matter with light mediators. Moreover, 21cm
surveys like HERA, SKA and LOFAR along with the James Webb Space Telescope
will allow us to push to the epoch of reionization and measure the small-scale
power spectrum [597]. Combined with the small-scale power spectrum, the mass
function of dark matter halos and progress in understanding galaxy formation
over a wide range of redshifts, we can be optimistic about making significant
strides in pinning down the microphysical nature of dark matter in the near
future.

###  6.2 Simulations and Modeling

With the wealth of data that will be available in the coming years, there is a
need for simultaneous advancement in simulation and modeling methods of
alternate DM models in order to use this data to constrain properties of DM.
For SIDM models of the type considered in this work, much of the progress has
been based on DM-only (i.e. without baryons) modified $N$-body simulations.
Different groups involved in the simulation efforts broadly agree on the types
of phenomena produced by the presence of elastic self-interactions, as well as
the size of these effects as a function of interaction cross section. On the
other hand, a number of powerful probes of SIDM presented in this work involve
modeling of small-scale structure, especially the evolution of substructure
inside a larger virialized object, where the different simulations techniques
have not been rigorously tested against each other. For collisionless
$\Lambda$CDM models, different implementations of $N$-body simulations have
been stringently calibrated against each other over the years, ensuring that
numerical artifacts from any specific implementation can be detected and
corrected [677, 678] (also see the AGORA code comparison project for
hydrodynamic simulation algorithms [679, 680, 681]). A code comparison project
of this type for different SIDM implementations will be hugely beneficial
toward the goal of predicting the effects of SIDM on small-scale structure
robustly (cf. [682]). Furthermore, since the small-scale structure of interest
often lie close to the resolution limits of these simulations, a thorough
understanding of whether the resolution affects the predictions of SIDM
effects, even within individual implementations, is needed.

In recent years, there has been a growing focus on the mechanism of disruption
of subhalos in collisionless CDM $N$-body simulations (see, e.g., [683, 684,
685, 466, 686]). It has been pointed out that much finer force resolution may
be needed to faithfully follow the tidal stripping of subhalos inside a more
massive host halo than what was estimated previously, and has been used in
much of the existing literature. Lack of sufficient force resolution results
in larger amounts of tidal stripping than is expected physically. A related
issue is that of subhalos being stripped to below the mass resolution of the
applied halo-finder and disappearing from the substructure catalog, even
though the central cusp retains its identity, and can continue to host a
luminous galaxy. While some of these issues can be overcome by additional
modeling, e.g. including an “orphan model" for galaxies living in the
substructure that fall below the halo finder detection threshold, these
numerical effects can have major implications on predictions of strong
lensing, satellite weak lensing signals and the phase space distributions of
satellites in the nearby Universe. An initial exploration of this effect was
done in [687], showing that accounting for disrupted orphans is important to
model satellite weak lensing in clusters and to understand the overall
distribution of subhalos. Additionally, it should be noted that tracking
orphans in SIDM simulations is non-trivial as interactions can often scatter
particles near the core of a subhalo, making methods like core-tracking
inherently difficult. A detailed understanding of the effects of artificial
disruption in SIDM $N$-body simulations is essential, especially as they can
be degenerate with certain physical signatures of DM self-interactions. Ram
pressure stripping of subhalos in SIDM simulations represent an additional
mass-loss mechanism, and the presence of a central core in SIDM subhalos
instead of a cusp can mean that the substructure can be completely destroyed.

Most algorithms that are currently employed to implement SIDM in $N$-body
simulations are designed explicitly to work in the limit of rare self-
interactions, i.e. when an individual simulation particle interacts via at
most one other particle through self-interactions in a given time-step. This
typically also corresponds to the regime where most particles do not scatter
in any particular time-step. On the other hand, certain SIDM models have been
proposed where it is necessary to model frequent self-interactions, i.e. where
a single particle can, in principle, interact with multiple other particles in
a single time step. The current simulation methods, therefore, need to be
updated to correctly predict the phenomenology of these models, a direction in
which there has already been some progress [114, 688]. The development of
these algorithms could also help robustly model the onset and subsequent
evolution of the core-collapse phase of SIDM in $N$-body simulations, where
interactions become progressively more frequent in the collapsing region and
its outskirts. Since the core-collapse mechanism is potentially interesting
for a whole range of observed phenomenology of dwarf satellites, advances in
simulation techniques will have large science impact. In this context, it is
also worth investigating whether $N$-body techniques are indeed the most
appropriate tools to model this regime, or whether gravothermal fluid
techniques in these regions can perform better.

DM-only simulations, of course, make predictions for the evolution of DM halos
and subhalos. On the other hand, cosmological surveys generally measure the
positions, and in some cases, velocities of luminous galaxies. The connection
between the observed galaxies and the DM halos in which they live needs to be
effectively modeled in order for the simulations to be used in the
interpretation and analysis of observational data. In the context of CDM, the
galaxy–halo connection has been studied extensively (see [33] for a recent
review), and various effective parameterizations have been developed to
marginalize over the unknown galaxy formation physics, and its potential
impact back-reacting onto the dynamics and evolution of DM substructure. Such
studies are also needed in the context of SIDM to promote results from DM-only
SIDM simulations to analysis tools. An important question to answer in this
context is whether the same parameterizations developed in the context of CDM
can be applied to the SIDM context, or whether the galaxy–halo connection in
SIDM is qualitatively different, and needs a completely different approach.

Along with a detailed scrutiny of cosmological simulations of DM it is also
essential to understand the varied aspects of the interplay between baryonic
physics and self-interacting DM over a wide range of halo mass scales. Ideally
one would like to simultaneously explore the parameter space in SIDM along
with varying the feedback mechanisms in hydrodynamic simulations.

At the very least, when using simulations to investigate how observable
quantities are affected by the DM model one should compare this with the
differences seen with a fixed DM model when different implementations of
baryonic physics are used.

Another avenue for future code development is to include a larger number of
possible DM interactions. To date, the bulk of work on SIDM has considered
elastic scattering, typically with an isotropic differential cross section.
There has been a small amount of work that considered anisotropic scattering
cross sections [49, 95], or inelastic scattering [606], but these scenarios
have not been exhaustively explored and do not cover the full range of
phenomenology that DM with self-interactions could possess.

###  6.3 Promising tests of simple SIDM models

We have described a variety of tests of SIDM models in this review, ranging
from the least massive dark subhalos that could perturb lensed images to
galaxy clusters and large-scale structure. It is fascinating that the simple
modification of including elastic self-interactions to the CDM model, which is
well-motivated from the particle physics perspective, could impact structure
formation in so many different ways. It is also sobering to realize that
despite these striking predictions, the remarkable progress in simulating
galaxies, and the plethora of data that already exists, there are no
definitive findings in favor of or against SIDM. The reasons include the fact
that CDM models are an excellent fit to the majority of data, and where it may
have trouble either star-formation physics (with large modeling uncertainties)
is important or the data is not yet decisive. SIDM models predict deviations
from the CDM model predictions predominantly in the central regions of
galaxies, where baryonic physics is important. This is a feature of simple
SIDM models and not a fine-tuning, as we have reviewed in Sec. 3.

We have identified multiple avenues to test SIDM through astrophysical
observables. A key prediction of all DM models is the mass function of dark
subhalos and their internal density profiles. Strong lensing studies aimed at
testing these predictions will mature in the coming years, as we have
discussed. Observables associated with resolved and unresolved subhalos and
their power spectrum should provide a definitive test of a wide swathe of SIDM
models that predict halos and subhalos are mostly in the core-expansion phase
of their evolution. We have also discussed the recent progress in fleshing out
the predictions of SIDM models that predict some subhalos could be in the
core-collapse phase. These subhalos with their steep density profiles in the
inner regions (radii smaller than the NFW scale radius) provide a striking
target for future observations of strongly-lensed systems and stellar streams.
Much more work remains to be done to work out the predictions for these models
and the critical correlations between subhalo’s density profile, survival and
mass loss due to stripping and evaporation.

We also highlighted the possibility of throwing the full power of lensing
(strong and weak) reconstructions, X-ray observations, BCG offsets, mergers
(major and minor) and stellar kinematics in groups and clusters of galaxies to
put definitive constraints on the cross sections at relative velocities of
$1000-2000\ \rm km/s$.

It is only possible to put upper limits currently but the massive influx of
data in the coming years, along with the possibility of improvements in theory
and modeling, makes this an exciting avenue to search for a positive
detection.

The most concrete astrophysical motivations for SIDM have been the core-cusp
and diversity issues of rotation curves, and the too-big-to-fail problem of
the Milky Way satellites. We have reviewed these problems in detail along with
the constraints imposed by the growing census of Milky Way satellites. In the
near future, there will be an explosion of data on galaxies and satellites at
low redshift, which will provide critical tests of SIDM models. In this
context we have highlighted the unique opportunities provided by the ultra-
faint dwarf spheroidals of the Milky Way, other dwarfs in the local volume and
the ultra-diffuse galaxies in the field and clusters. It is safe to say that
there is not currently a concrete understanding of the dark matter or stellar
distribution of these objects in the context of either CDM or SIDM models. As
with the dark subhalos, the prevalence of core expansion vs core collapse
phases in the satellites will be critical. It seems clear that SIDM models
where interactions are not strong enough to push satellites into the core
collapse phase will be definitively tested in the near future.

SIDM models provide a viable, predictive and versatile framework for
interpreting astrophysical observables on galactic and sub-galactic scales.

The predictions of these models need to be developed further if we are to
maximize the science returns from the amazing data expected this decade and
next. SIDM models provide a unique and compelling opportunity to detect the
particle nature of DM through astrophysical observables, as we have
highlighted throughout this review. They also provide a well-defined and
predictive foil to the dominant CDM paradigm; with predictions for simple SIDM
models in hand, the tests of the CDM model become clearer. However, progress
will only be possible if there is a serious investment of resources in
simulating SIDM models. We have discussed many aspects of the required
simulations throughout this review. This investment in simulations is critical
given the incomplete progress in mapping out the small-scale structure of even
the simplest (single-component, elastic-scattering-only) SIDM models. It is
also vital because progress in understanding galaxy formation and the
phenomenology of DM self-interaction are intertwined. We need a concerted
effort on the theory side to take full advantage of the incredible data that
we will have in the near future and develop the astrophysical observables we
have reviewed into decisive probes of dark sector physics.

## Acknowledgements

Marisa March helped get this review and the Novel Probes Project off the
ground. We are very grateful for her energy, vision, and support in the early
stages of this work. We thank Alex Drlica-Wagner, Matthew Walker and Matthew
Buckley for early discussions.

KB acknowledges support from the NSF under Grant No. PHY-2112884. HD was
supported by a Junior Research Fellowship at St John’s College, Oxford, a
McWilliams Fellowship at Carnegie Mellon University and a Royal Society
University Research Fellowship (grant no. 211046). CD was partially supported
by the Department of Energy (DOE) Grant No. DE-SC0020223. BJ is supported in
part by the US Department of Energy Grant No. DE- SC0007901. FK was funded by
the Deutsche Forschungsgemeinschaft (DFG) through the Emmy Noether Grant No.
KA 4662/1-1. MK acknowledges support from the NSF under Grant No. PHY-1915005.
AHGP and FYCR were supported in part by the National Aeronautics and Space
Administration under award number 80NSSC18K1014. JZ acknowledges support by a
Grant of Excellence from the Icelandic Research fund (grant number 206930).

## References

  * [1] Planck, N. Aghanim et al.,  1807.06209arXiv1807.06209arXiv:1807.06209.
  * [2] LSST Dark Energy Science, C. Chang et al.,  1808.07335arXiv1808.07335arXiv:1808.07335.
  * [3] L. J. Hall, K. Jedamzik, J. March-Russell and S. M. West,  JHEP 03, 080 (2010), [0911.1120arXiv0911.1120arXiv:0911.1120].
  * [4] T. Bringmann, F. Kahlhoefer, K. Schmidt-Hoberg and P. Walia,  Phys. Rev. Lett. 118, 141802 (2017), [1612.00845arXiv1612.00845arXiv:1612.00845].
  * [5] R. Huo, M. Kaplinghat, Z. Pan and H.-B. Yu,  1709.09717arXiv1709.09717arXiv:1709.09717.
  * [6] T. Hambye, M. H. G. Tytgat, J. Vandecasteele and L. Vanderheyden,  Phys. Rev. D98, 075017 (2018), [1807.05022arXiv1807.05022arXiv:1807.05022].
  * [7] A. Aboubrahim, W.-Z. Feng, P. Nath and Z.-Y. Wang,  Phys. Rev. D 103, 075014 (2021), [2008.00529arXiv2008.00529arXiv:2008.00529].
  * [8] J. L. Feng, M. Kaplinghat, H. Tu and H.-B. Yu,  JCAP 07, 004 (2009), [0905.3039arXiv0905.3039arXiv:0905.3039].
  * [9] M. R. Buckley and P. J. Fox,  Phys. Rev. D81, 083522 (2010), [0911.3898arXiv0911.3898arXiv:0911.3898].
  * [10] A. Loeb and N. Weiner,  Phys. Rev. Lett. 106, 171302 (2011), [1011.6374arXiv1011.6374arXiv:1011.6374].
  * [11] S. Tulin, H.-B. Yu and K. M. Zurek,  Phys. Rev. Lett. 110, 111301 (2013), [1210.0900arXiv1210.0900arXiv:1210.0900].
  * [12] X. Chu, C. Garcia-Cely and H. Murayama,  Phys. Rev. Lett. 124, 041101 (2020), [1901.00075arXiv1901.00075arXiv:1901.00075].
  * [13] X. Chu, C. Garcia-Cely and H. Murayama,  Phys. Rev. Lett. 122, 071103 (2019), [1810.04709arXiv1810.04709arXiv:1810.04709].
  * [14] Y.-D. Tsai, R. McGehee and H. Murayama,  Phys. Rev. Lett. 128, 172001 (2022), [2008.08608arXiv2008.08608arXiv:2008.08608].
  * [15] W. J. G. de Blok and S. S. McGaugh,  MNRAS290, 533 (1997), [astro-ph/9704274arXivastro-ph/9704274arXiv:astro-ph/9704274].
  * [16] F. C. van den Bosch, B. E. Robertson, J. J. Dalcanton and W. J. G. de Blok,  AJ119, 1579 (2000), [astro-ph/9911372arXivastro-ph/9911372arXiv:astro-ph/9911372].
  * [17] W. J. G. de Blok, S. S. McGaugh, A. Bosma and V. C. Rubin,  ApJ552, L23 (2001), [astro-ph/0103102arXivastro-ph/0103102arXiv:astro-ph/0103102].
  * [18] R. Kuzio de Naray, S. S. McGaugh and W. J. G. de Blok,  ApJ676, 920 (2008), [arXiv:0712.0860arXivarXiv:0712.0860arXiv:arXiv:0712.0860].
  * [19] A. A. Klypin, A. V. Kravtsov, O. Valenzuela and F. Prada,  Astrophys. J. 522, 82 (1999), [astro-ph/9901240arXivastro-ph/9901240arXiv:astro-ph/9901240].
  * [20] B. Moore et al.,  Astrophys. J. 524, L19 (1999), [astro-ph/9907411arXivastro-ph/9907411arXiv:astro-ph/9907411].
  * [21] E. D. Carlson, M. E. Machacek and L. J. Hall,  Astrophys. J. 398, 43 (1992).
  * [22] F. Governato et al.,  Nature 463, 203 (2010), [0911.2237arXiv0911.2237arXiv:0911.2237].
  * [23] J. F. Navarro, C. S. Frenk and S. D. M. White,  Astrophys. J. 462, 563 (1996), [astro-ph/9508025arXivastro-ph/9508025arXiv:astro-ph/9508025].
  * [24] S.-H. Oh et al.,  AJ142, 24 (2011), [1011.2777arXiv1011.2777arXiv:1011.2777].
  * [25] D. N. Spergel and P. J. Steinhardt,  Phys. Rev. Lett. 84, 3760 (2000), [astro-ph/9909386arXivastro-ph/9909386arXiv:astro-ph/9909386].
  * [26] S. Hofmann, D. J. Schwarz and H. Stoecker,  Phys. Rev. D64, 083507 (2001), [astro-ph/0104173arXivastro-ph/0104173arXiv:astro-ph/0104173].
  * [27] A. M. Green, S. Hofmann and D. J. Schwarz,  Mon. Not. Roy. Astron. Soc. 353, L23 (2004), [astro-ph/0309621arXivastro-ph/0309621arXiv:astro-ph/0309621].
  * [28] S. Profumo, K. Sigurdson and M. Kamionkowski,  Phys. Rev. Lett. 97, 031301 (2006), [astro-ph/0603373arXivastro-ph/0603373arXiv:astro-ph/0603373].
  * [29] J. Diemand, M. Kuhlen and P. Madau,  ApJ649, 1 (2006), [astro-ph/0603250arXivastro-ph/0603250arXiv:astro-ph/0603250].
  * [30] L. Gao, S. D. M. White, A. Jenkins, F. Stoehr and V. Springel,  MNRAS355, 819 (2004), [astro-ph/0404589arXivastro-ph/0404589arXiv:astro-ph/0404589].
  * [31] J. Diemand, M. Kuhlen and P. Madau,  ApJ667, 859 (2007), [astro-ph/0703337arXivastro-ph/0703337arXiv:astro-ph/0703337].
  * [32] V. Springel et al.,  MNRAS391, 1685 (2008), [0809.0898arXiv0809.0898arXiv:0809.0898].
  * [33] R. H. Wechsler and J. L. Tinker,  Submitted to: Ann. Rev. Astron. Astrophys. (2018), [1804.03097arXiv1804.03097arXiv:1804.03097].
  * [34] R. Barkana and A. Loeb,  Astrophys. J. 523, 54 (1999), [astro-ph/9901114arXivastro-ph/9901114arXiv:astro-ph/9901114].
  * [35] J. S. Bullock, A. V. Kravtsov and D. H. Weinberg,  Astrophys. J. 539, 517 (2000), [astro-ph/0002214arXivastro-ph/0002214arXiv:astro-ph/0002214].
  * [36] A. J. Benson, C. S. Frenk, C. G. Lacey, C. M. Baugh and S. Cole,  Mon. Not. Roy. Astron. Soc. 333, 177 (2002), [astro-ph/0108218arXivastro-ph/0108218arXiv:astro-ph/0108218].
  * [37] R. S. Somerville,  ApJ572, L23 (2002), [astro-ph/0107507arXivastro-ph/0107507arXiv:astro-ph/0107507].
  * [38] T. Okamoto, L. Gao and T. Theuns,  MNRAS390, 920 (2008), [0806.0378arXiv0806.0378arXiv:0806.0378].
  * [39] P. Chen, J. H. Wise, M. L. Norman, H. Xu and B. W. O’Shea,  Astrophys. J. 795, 144 (2014), [1408.2523arXiv1408.2523arXiv:1408.2523].
  * [40] B. Smith, J. Wise, B. O’Shea, M. Norman and S. Khochfar,  Mon. Not. Roy. Astron. Soc. 452, 2822 (2015), [1504.07639arXiv1504.07639arXiv:1504.07639].
  * [41] W. J. G. de Blok, S. S. McGaugh, A. Bosma and V. C. Rubin,  ApJ552, L23 (2001), [astro-ph/0103102arXivastro-ph/0103102arXiv:astro-ph/0103102].
  * [42] J. D. Simon, A. D. Bolatto, A. Leroy, L. Blitz and E. L. Gates,  Astrophys. J. 621, 757 (2005), [astro-ph/0412035arXivastro-ph/0412035arXiv:astro-ph/0412035].
  * [43] J. J. Adams et al.,  Astrophys. J. 789, 63 (2014), [1405.4854arXiv1405.4854arXiv:1405.4854].
  * [44] R. Kuzio de Naray, G. D. Martinez, J. S. Bullock and M. Kaplinghat,  Astrophys. J. Lett. 710, L161 (2010), [0912.3518arXiv0912.3518arXiv:0912.3518].
  * [45] K. A. Oman et al.,  Mon. Not. Roy. Astron. Soc. 452, 3650 (2015), [1504.01437arXiv1504.01437arXiv:1504.01437].
  * [46] P. Creasey et al.,  MNRAS468, 2283 (2017), [1612.03903arXiv1612.03903arXiv:1612.03903].
  * [47] M. Valli and H.-B. Yu,  1711.03502arXiv1711.03502arXiv:1711.03502.
  * [48] A. Kamada, M. Kaplinghat, A. B. Pace and H.-B. Yu,  Phys. Rev. Lett. 119, 111102 (2017), [1611.02716arXiv1611.02716arXiv:1611.02716].
  * [49] A. Robertson et al.,  Mon. Not. Roy. Astron. Soc. 476, L20 (2018), [1711.09096arXiv1711.09096arXiv:1711.09096].
  * [50] T. Ren, A. Kwa, M. Kaplinghat and H.-B. Yu,  1808.05695arXiv1808.05695arXiv:1808.05695.
  * [51] J. Zavala, M. R. Lovell, M. Vogelsberger and J. D. Burger,  1904.09998arXiv1904.09998arXiv:1904.09998.
  * [52] A. B. Newman et al.,  ApJ765, 24 (2013), [1209.1391arXiv1209.1391arXiv:1209.1391].
  * [53] A. B. Newman, T. Treu, R. S. Ellis and D. J. Sand,  ApJ765, 25 (2013), [1209.1392arXiv1209.1392arXiv:1209.1392].
  * [54] A. B. Newman, R. S. Ellis and T. Treu,  ApJ814, 26 (2015), [1503.05282arXiv1503.05282arXiv:1503.05282].
  * [55] M. Kaplinghat, S. Tulin and H.-B. Yu,  Phys. Rev. Lett. 116, 041302 (2016), [1508.03339arXiv1508.03339arXiv:1508.03339].
  * [56] LSST Dark Matter Group, A. Drlica-Wagner et al.,  1902.01055arXiv1902.01055arXiv:1902.01055.
  * [57] DES, E. O. Nadler et al.,  Phys. Rev. Lett. 126, 091101 (2021), [2008.00022arXiv2008.00022arXiv:2008.00022].
  * [58] B. Willman et al.,  Astrophys. J. 626, L85 (2005), [astro-ph/0503552arXivastro-ph/0503552arXiv:astro-ph/0503552].
  * [59] DES, A. Drlica-Wagner et al.,  Astrophys. J. 813, 109 (2015), [1508.03622arXiv1508.03622arXiv:1508.03622].
  * [60] D. Kim, H. Jerjen, D. Mackey, G. S. Da Costa and A. P. Milone,  Astrophys. J. 804, L44 (2015), [1503.08268arXiv1503.08268arXiv:1503.08268].
  * [61] B. P. M. Laevens et al.,  Astrophys. J. 813, 44 (2015), [1507.07564arXiv1507.07564arXiv:1507.07564].
  * [62] N. Caldwell et al.,  Astrophys. J. 839, 20 (2017), [1612.06398arXiv1612.06398arXiv:1612.06398].
  * [63] D. Homma et al.,  Astrophys. J. 832, 21 (2016), [1609.04346arXiv1609.04346arXiv:1609.04346].
  * [64] G. Torrealba et al.,  Mon. Not. Roy. Astron. Soc. 475, 5085 (2018), [1801.07279arXiv1801.07279arXiv:1801.07279].
  * [65] E. J. Tollerud, J. S. Bullock, L. E. Strigari and B. Willman,  ApJ688, 277 (2008), [0806.4381arXiv0806.4381arXiv:0806.4381].
  * [66] S. E. Koposov et al.,  Astrophys. J. 696, 2179 (2009), [0901.2116arXiv0901.2116arXiv:0901.2116].
  * [67] J. R. Hargis, B. Willman and A. H. G. Peter,  ApJ795, L13 (2014), [1407.4470arXiv1407.4470arXiv:1407.4470].
  * [68] P. Jethwa, D. Erkal and V. Belokurov,  Mon. Not. Roy. Astron. Soc. 473, 2060 (2018), [1612.07834arXiv1612.07834arXiv:1612.07834].
  * [69] S. Y. Kim, A. H. G. Peter and J. R. Hargis,  1711.06267arXiv1711.06267arXiv:1711.06267.
  * [70] O. Newton, M. Cautun, A. Jenkins, C. S. Frenk and J. Helly,  ArXiv e-prints (2017), [1708.04247arXiv1708.04247arXiv:1708.04247].
  * [71] S. Bose, A. J. Deason and C. S. Frenk,  Astrophys. J. 863, 123 (2018), [1802.10096arXiv1802.10096arXiv:1802.10096].
  * [72] E. O. Nadler, Y.-Y. Mao, G. M. Green and R. H. Wechsler,  1809.05542arXiv1809.05542arXiv:1809.05542.
  * [73] A. S. Graus et al.,  ArXiv e-prints (2018), [1808.03654arXiv1808.03654arXiv:1808.03654].
  * [74] M. Vogelsberger, J. Zavala and A. Loeb,  Mon. Not. Roy. Astron. Soc. 423, 3740 (2012), [1201.5892arXiv1201.5892arXiv:1201.5892].
  * [75] G. A. Dooley, A. H. G. Peter, M. Vogelsberger, J. Zavala and A. Frebel,  Mon. Not. Roy. Astron. Soc. 461, 710 (2016), [1603.08919arXiv1603.08919arXiv:1603.08919].
  * [76] S. Tulin and H.-B. Yu,  Phys. Rept. 730, 1 (2018), [1705.02358arXiv1705.02358arXiv:1705.02358].
  * [77] J. S. Bullock and M. Boylan-Kolchin,  Ann. Rev. Astron. Astrophys. 55, 343 (2017), [1707.04256arXiv1707.04256arXiv:1707.04256].
  * [78] T. Baker et al.,  Rev. Mod. Phys. 93, 015003 (2021), [1908.03430arXiv1908.03430arXiv:1908.03430].
  * [79] H. W. Lin and A. Loeb,  JCAP 1603, 009 (2016), [1506.05471arXiv1506.05471arXiv:1506.05471].
  * [80] M. R. Buckley and E. T. Neil,  Phys. Rev. D87, 043510 (2013), [1209.6054arXiv1209.6054arXiv:1209.6054].
  * [81] S. Tulin, H.-B. Yu and K. M. Zurek,  Phys. Rev. D87, 115007 (2013), [1302.3898arXiv1302.3898arXiv:1302.3898].
  * [82] Y. Hochberg, E. Kuflik, T. Volansky and J. G. Wacker,  Phys. Rev. Lett. 113, 171301 (2014), [1402.5143arXiv1402.5143arXiv:1402.5143].
  * [83] Y. Hochberg, E. Kuflik and H. Murayama,  JHEP 05, 090 (2016), [1512.07917arXiv1512.07917arXiv:1512.07917].
  * [84] J. M. Cline, Z. Liu, G. Moore and W. Xue,  Phys. Rev. D89, 043514 (2014), [1311.6468arXiv1311.6468arXiv:1311.6468].
  * [85] K. K. Boddy, J. L. Feng, M. Kaplinghat and T. M. P. Tait,  Phys. Rev. D89, 115017 (2014), [1402.3629arXiv1402.3629arXiv:1402.3629].
  * [86] K. K. Boddy, M. Kaplinghat, A. Kwa and A. H. G. Peter,  Phys. Rev. D94, 123017 (2016), [1609.03592arXiv1609.03592arXiv:1609.03592].
  * [87] B. Colquhoun, S. Heeba, F. Kahlhoefer, L. Sagunski and S. Tulin,  Phys. Rev. D 103, 035006 (2021), [2011.04679arXiv2011.04679arXiv:2011.04679].
  * [88] L. Pitaevskii and E. Lifshitz,  Physical kinetics, 2012.
  * [89] S. Chapman, T. Cowling, D. Burnett and C. Cercignani,  The Mathematical Theory of Non-uniform Gases: An Account of the Kinetic Theory of Viscosity, Thermal Conduction and Diffusion in GasesCambridge Mathematical Library (Cambridge University Press, 1990).
  * [90] N. J. Outmezguine, K. K. Boddy, S. Gad-Nasr, M. Kaplinghat and L. Sagunski,  2204.06568arXiv2204.06568arXiv:2204.06568.
  * [91] F. Kahlhoefer, K. Schmidt-Hoberg, M. T. Frandsen and S. Sarkar,  Mon. Not. Roy. Astron. Soc. 437, 2865 (2014), [1308.3419arXiv1308.3419arXiv:1308.3419].
  * [92] A. Robertson, R. Massey and V. Eke,  Mon. Not. Roy. Astron. Soc. 467, 4719 (2017), [1612.03906arXiv1612.03906arXiv:1612.03906].
  * [93] P. Agrawal, F.-Y. Cyr-Racine, L. Randall and J. Scholtz,  JCAP 1705, 022 (2017), [1610.04611arXiv1610.04611arXiv:1610.04611].
  * [94] D. Yang and H.-B. Yu,  2205.03392arXiv2205.03392arXiv:2205.03392.
  * [95] A. Banerjee, S. Adhikari, N. Dalal, S. More and A. Kravtsov,  1906.12026arXiv1906.12026arXiv:1906.12026.
  * [96] V. Springel,  Annual Review of Astronomy and Astrophysics 48, 391 (2010), [https://doi.org/10.1146/annurev-astro-081309-130914arXivhttps://doi.org/10.1146/annurev-astro-081309-130914arXiv:https://doi.org/10.1146/annurev-astro-081309-130914].
  * [97] R. Teyssier,  Annual Review of Astronomy and Astrophysics 53, 325 (2015), [https://doi.org/10.1146/annurev-astro-082214-122309arXivhttps://doi.org/10.1146/annurev-astro-082214-122309arXiv:https://doi.org/10.1146/annurev-astro-082214-122309].
  * [98] J. I. Read,  Journal of Physics G Nuclear Physics 41, 063101 (2014), [1404.1938arXiv1404.1938arXiv:1404.1938].
  * [99] S. Balberg, S. L. Shapiro and S. Inagaki,  ApJ568, 475 (2002), [astro-ph/0110561arXivastro-ph/0110561arXiv:astro-ph/0110561].
  * [100] J. Koda and P. R. Shapiro,  MNRAS415, 1125 (2011), [1101.3097arXiv1101.3097arXiv:1101.3097].
  * [101] S. L. Shapiro and V. Paschalidis,  Phys. Rev. D89, 023506 (2014), [1402.0005arXiv1402.0005arXiv:1402.0005].
  * [102] J. Pollack, D. N. Spergel and P. J. Steinhardt,  ApJ804, 131 (2015), [1501.00017arXiv1501.00017arXiv:1501.00017].
  * [103] H. Nishikawa, K. K. Boddy and M. Kaplinghat,  arXiv e-prints (2019), [1901.00499arXiv1901.00499arXiv:1901.00499].
  * [104] S. Yang et al.,  2205.02957arXiv2205.02957arXiv:2205.02957.
  * [105] C. S. Kochanek and M. White,  ApJ543, 514 (2000), [astro-ph/0003483arXivastro-ph/0003483arXiv:astro-ph/0003483].
  * [106] M. Rocha et al.,  MNRAS430, 81 (2013), [1208.3025arXiv1208.3025arXiv:1208.3025].
  * [107] A. Robertson, R. Massey and V. Eke,  MNRAS465, 569 (2017), [1605.04307arXiv1605.04307arXiv:1605.04307].
  * [108] M. S. Fischer et al.,  Mon. Not. Roy. Astron. Soc. 505, 851 (2021), [2012.10277arXiv2012.10277arXiv:2012.10277].
  * [109] M. S. Fischer et al.,  Mon. Not. Roy. Astron. Soc. 510, 4080 (2022), [2109.10035arXiv2109.10035arXiv:2109.10035].
  * [110] M. S. Fischer et al.,  2205.02243arXiv2205.02243arXiv:2205.02243.
  * [111] L. Hernquist,  ApJ356, 359 (1990).
  * [112] D. Harvey, R. Massey, T. Kitching, A. Taylor and E. Tittley,  Science 347, 1462 (2015), [1503.07675arXiv1503.07675arXiv:1503.07675].
  * [113] M. Kesden and M. Kamionkowski,  Phys. Rev. D74, 083007 (2006), [astro-ph/0608095arXivastro-ph/0608095arXiv:astro-ph/0608095].
  * [114] J. Kummer, M. Brüggen, K. Dolag, F. Kahlhoefer and K. Schmidt-Hoberg,  1902.02330arXiv1902.02330arXiv:1902.02330.
  * [115] M. S. Fischer et al.,  arXiv e-prints , arXiv:2012.10277 (2020), [2012.10277arXiv2012.10277arXiv:2012.10277].
  * [116] J. Koda, M. Yagi, H. Yamanoi and Y. Komiyama,  Astrophys. J. 807, L2 (2015), [1506.01712arXiv1506.01712arXiv:1506.01712].
  * [117] O. D. Elbert et al.,  Mon. Not. Roy. Astron. Soc. 453, 29 (2015), [1412.1477arXiv1412.1477arXiv:1412.1477].
  * [118] R. Essig, H.-B. Yu, Y.-M. Zhong and S. D. Mcdermott,  1809.01144arXiv1809.01144arXiv:1809.01144.
  * [119] M. Kaplinghat, R. E. Keeley, T. Linden and H.-B. Yu,  Phys. Rev. Lett. 113, 021302 (2014), [1311.6524arXiv1311.6524arXiv:1311.6524].
  * [120] A. Robertson, R. Massey, V. Eke, J. Schaye and T. Theuns,  MNRAS501, 4610 (2021), [2009.07844arXiv2009.07844arXiv:2009.07844].
  * [121] O. Sameie et al.,  Mon. Not. Roy. Astron. Soc. 479, 359 (2018), [1801.09682arXiv1801.09682arXiv:1801.09682].
  * [122] P. Colín, V. Avila-Reese, O. Valenzuela and C. Firmani,  ApJ581, 777 (2002), [astro-ph/0205322arXivastro-ph/0205322arXiv:astro-ph/0205322].
  * [123] D. Lynden-Bell and R. Wood,  MNRAS138, 495 (1968).
  * [124] D. Lynden-Bell and P. P. Eggleton,  MNRAS191, 483 (1980).
  * [125] S. L. Shapiro and V. Paschalidis,  Phys. Rev. D 89, 023506 (2014), [1402.0005arXiv1402.0005arXiv:1402.0005].
  * [126] S. L. Shapiro,  Phys. Rev. D 98, 023021 (2018), [1809.02618arXiv1809.02618arXiv:1809.02618].
  * [127] Z. C. Zeng et al.,  2110.00259arXiv2110.00259arXiv:2110.00259.
  * [128] O. Sameie, H.-B. Yu, L. V. Sales, M. Vogelsberger and J. Zavala,  1904.07872arXiv1904.07872arXiv:1904.07872.
  * [129] F. Kahlhoefer, M. Kaplinghat, T. R. Slatyer and C.-L. Wu,  1904.10539arXiv1904.10539arXiv:1904.10539.
  * [130] H. C. Turner, M. R. Lovell, J. Zavala and M. Vogelsberger,  Mon. Not. Roy. Astron. Soc. 505, 5327 (2021), [2010.02924arXiv2010.02924arXiv:2010.02924].
  * [131] M. Kaplinghat, M. Valli and H.-B. Yu,  1904.04939arXiv1904.04939arXiv:1904.04939.
  * [132] O. D. Elbert et al.,  ApJ853, 109 (2018), [1609.08626arXiv1609.08626arXiv:1609.08626].
  * [133] D. Yang and H.-B. Yu,  Phys. Rev. D 104, 103031 (2021), [2102.02375arXiv2102.02375arXiv:2102.02375].
  * [134] F. Jiang et al.,  2206.12425arXiv2206.12425arXiv:2206.12425.
  * [135] W.-X. Feng, H.-B. Yu and Y.-M. Zhong,  Astrophys. J. Lett. 914, L26 (2021), [2010.15132arXiv2010.15132arXiv:2010.15132].
  * [136] S. Balberg and S. L. Shapiro,  Phys. Rev. Lett. 88, 101301 (2002), [astro-ph/0111176arXivastro-ph/0111176arXiv:astro-ph/0111176].
  * [137] W.-X. Feng, H.-B. Yu and Y.-M. Zhong,  JCAP 05, 036 (2022), [2108.11967arXiv2108.11967arXiv:2108.11967].
  * [138] S. Chandrasekhar,  ApJ97, 255 (1943).
  * [139] S. D. M. White,  ApJ274, 53 (1983).
  * [140] M. D. Weinberg,  MNRAS239, 549 (1989).
  * [141] S. Tremaine and M. D. Weinberg,  MNRAS209, 729 (1984).
  * [142] W. A. Mulder,  A&A117, 9 (1983).
  * [143] M. D. Weinberg,  ApJ300, 93 (1986).
  * [144] M. Boylan-Kolchin, C.-P. Ma and E. Quataert,  Mon. Not. Roy. Astron. Soc. 383, 93 (2008), [0707.2960arXiv0707.2960arXiv:0707.2960].
  * [145] J. I. Read et al.,  MNRAS373, 1451 (2006), [astro-ph/0606636arXivastro-ph/0606636arXiv:astro-ph/0606636].
  * [146] T. Goerdt, J. I. Read, B. Moore and J. Stadel,  Astrophys. J. 725, 1707 (2010), [0806.1951arXiv0806.1951arXiv:0806.1951].
  * [147] J. A. Petts, A. Gualandris and J. I. Read,  MNRAS454, 3778 (2015), [1509.07871arXiv1509.07871arXiv:1509.07871].
  * [148] J. A. Petts, J. I. Read and A. Gualandris,  MNRAS463, 858 (2016), [1607.04284arXiv1607.04284arXiv:1607.04284].
  * [149] S. Y. Kim, A. H. G. Peter and D. Wittman,  Mon. Not. Roy. Astron. Soc. 469, 1414 (2017), [1608.08630arXiv1608.08630arXiv:1608.08630].
  * [150] E. O. Nadler et al.,  Phys. Rev. Lett.126, 091101 (2021), [2008.00022arXiv2008.00022arXiv:2008.00022].
  * [151] S. Y. Kim and A. H. G. Peter,  arXiv e-prints , arXiv:2106.09050 (2021), [2106.09050arXiv2106.09050arXiv:2106.09050].
  * [152] J. Kummer, F. Kahlhoefer and K. Schmidt-Hoberg,  Mon. Not. Roy. Astron. Soc. 474, 388 (2018), [1706.04794arXiv1706.04794arXiv:1706.04794].
  * [153] E. O. Nadler, A. Banerjee, S. Adhikari, Y.-Y. Mao and R. H. Wechsler,  ApJ896, 112 (2020), [2001.08754arXiv2001.08754arXiv:2001.08754].
  * [154] S. W. Randall, M. Markevitch, D. Clowe, A. H. Gonzalez and M. Bradac,  Astrophys. J. 679, 1173 (2008), [0704.0261arXiv0704.0261arXiv:0704.0261].
  * [155] J. Penarrubia et al.,  Mon. Not. Roy. Astron. Soc. 406, 1290 (2010), [1002.3376arXiv1002.3376arXiv:1002.3376].
  * [156] O. Y. Gnedin and J. P. Ostriker,  Astrophys. J. 561, 61 (2001), [astro-ph/0010436arXivastro-ph/0010436arXiv:astro-ph/0010436].
  * [157] A. M. Brooks, M. Kuhlen, A. Zolotov and D. Hooper,  Astrophys. J. 765, 22 (2013), [1209.5394arXiv1209.5394arXiv:1209.5394].
  * [158] E. D’Onghia, V. Springel, L. Hernquist and D. Keres,  Astrophys. J. 709, 1138 (2010), [0907.3482arXiv0907.3482arXiv:0907.3482].
  * [159] T. Sawala et al.,  MNRAS467, 4383 (2017), [1609.01718arXiv1609.01718arXiv:1609.01718].
  * [160] S. Garrison-Kimmel et al.,  Mon. Not. Roy. Astron. Soc. 471, 1709 (2017), [1701.03792arXiv1701.03792arXiv:1701.03792].
  * [161] A. M. Nierenberg et al.,  Mon. Not. Roy. Astron. Soc. 471, 2224 (2017), [1701.05188arXiv1701.05188arXiv:1701.05188].
  * [162] A. J. Deason, Y.-Y. Mao and R. H. Wechsler,  Astrophys. J. 821, 5 (2016), [1601.07905arXiv1601.07905arXiv:1601.07905].
  * [163] T. Carleton et al.,  MNRAS485, 382 (2019), [1805.06896arXiv1805.06896arXiv:1805.06896].
  * [164] A. Robertson et al.,  MNRAS488, 3646 (2019), [1810.05649arXiv1810.05649arXiv:1810.05649].
  * [165] A. H. G. Peter, M. Rocha, J. S. Bullock and M. Kaplinghat,  Mon. Not. Roy. Astron. Soc. 430, 105 (2013), [1208.3026arXiv1208.3026arXiv:1208.3026].
  * [166] T. Brinckmann, J. Zavala, D. Rapetti, S. H. Hansen and M. Vogelsberger,  MNRAS474, 746 (2018), [1705.00623arXiv1705.00623arXiv:1705.00623].
  * [167] M. Markevitch et al.,  Astrophys. J. 606, 819 (2004), [astro-ph/0309303arXivastro-ph/0309303arXiv:astro-ph/0309303].
  * [168] L. Ackerman, M. R. Buckley, S. M. Carroll and M. Kamionkowski,  Phys. Rev. D79, 023519 (2009), [0810.5126arXiv0810.5126arXiv:0810.5126],  [,277(2008)].
  * [169] D. Clowe, A. Gonzalez and M. Markevitch,  Astrophys. J. 604, 596 (2004), [astro-ph/0312273arXivastro-ph/0312273arXiv:astro-ph/0312273].
  * [170] D. Harvey, A. Robertson, R. Massey and J.-P. Kneib,  Mon. Not. Roy. Astron. Soc. 464, 3991 (2017), [1610.05327arXiv1610.05327arXiv:1610.05327].
  * [171] P. Taylor et al.,  Mon. Not. Roy. Astron. Soc. 468, 5004 (2017), [1701.04412arXiv1701.04412arXiv:1701.04412].
  * [172] L. L. R. Williams and P. Saha,  Mon. Not. Roy. Astron. Soc. 415, 448 (2011), [1102.3943arXiv1102.3943arXiv:1102.3943].
  * [173] D. Wittman, N. Golovich and W. A. Dawson,  1701.05877arXiv1701.05877arXiv:1701.05877.
  * [174] B. Diemer and A. V. Kravtsov,  Astrophys. J. 789, 1 (2014), [1401.1216arXiv1401.1216arXiv:1401.1216].
  * [175] S. Adhikari, N. Dalal and R. T. Chamberlain,  JCAP 1411, 019 (2014), [1409.4482arXiv1409.4482arXiv:1409.4482].
  * [176] S. More et al.,  Astrophys. J. 825, 39 (2016), [1601.06063arXiv1601.06063arXiv:1601.06063].
  * [177] E. Baxter et al.,  Astrophys. J. 841, 18 (2017), [1702.01722arXiv1702.01722arXiv:1702.01722].
  * [178] Dark Energy Survey, C. Chang et al.,  1710.06808arXiv1710.06808arXiv:1710.06808.
  * [179] F.-Y. Cyr-Racine et al.,  Phys. Rev. D93, 123527 (2016), [1512.05344arXiv1512.05344arXiv:1512.05344].
  * [180] J. Silk,  Astrophys. J. 151, 459 (1968).
  * [181] B. Holdom,  Physics Letters B 166, 196 (1986).
  * [182] B. Holdom,  Physics Letters B 178, 65 (1986).
  * [183] A. D. Dolgov, S. L. Dubovsky, G. I. Rubtsov and I. I. Tkachev,  Phys. Rev. D88, 117701 (2013), [1310.2376arXiv1310.2376arXiv:1310.2376].
  * [184] S. D. McDermott, H.-B. Yu and K. M. Zurek,  Phys. Rev. D 83, 063509 (2011), [1011.2907arXiv1011.2907arXiv:1011.2907].
  * [185] J. L. Feng, M. Kaplinghat and H.-B. Yu,  Phys. Rev. Lett. 104, 151301 (2010), [0911.0422arXiv0911.0422arXiv:0911.0422].
  * [186] P. Ko and Y. Tang,  Phys. Lett. B762, 462 (2016), [1608.01083arXiv1608.01083arXiv:1608.01083].
  * [187] P. Ko, N. Nagata and Y. Tang,  Phys. Lett. B773, 513 (2017), [1706.05605arXiv1706.05605arXiv:1706.05605].
  * [188] H. Goldberg and L. J. Hall,  Phys. Lett. B174, 151 (1986).
  * [189] R. Foot and S. Mitra,  Phys. Rev. D66, 061301 (2002), [hep-ph/0204256arXivhep-ph/0204256arXiv:hep-ph/0204256].
  * [190] D. E. Kaplan, G. Z. Krnjaic, K. R. Rehermann and C. M. Wells,  J. Cosmology Astropart. Phys1005, 021 (2010), [0909.0753arXiv0909.0753arXiv:0909.0753].
  * [191] D. E. Kaplan, G. Z. Krnjaic, K. R. Rehermann and C. M. Wells,  J. Cosmology Astropart. Phys1110, 011 (2011), [1105.2073arXiv1105.2073arXiv:1105.2073].
  * [192] F.-Y. Cyr-Racine and K. Sigurdson,  Phys. Rev. D87, 103515 (2013), [1209.5752arXiv1209.5752arXiv:1209.5752].
  * [193] J. M. Cline, Z. Liu and W. Xue,  Phys. Rev. D 85, 101302 (2012), [1201.4858arXiv1201.4858arXiv:1201.4858].
  * [194] J. Fan, A. Katz, L. Randall and M. Reece,  Phys. Rev. Lett. 110, 211302 (2013), [1303.3271arXiv1303.3271arXiv:1303.3271].
  * [195] J. Fan, A. Katz, L. Randall and M. Reece,  Phys. Dark Univ. 2, 139 (2013), [1303.1521arXiv1303.1521arXiv:1303.1521].
  * [196] J. M. Cline, Z. Liu, G. Moore and W. Xue,  Phys. Rev. D 90, 015023 (2014), [1312.3325arXiv1312.3325arXiv:1312.3325].
  * [197] J. Choquette and J. M. Cline,  1509.05764arXiv1509.05764arXiv:1509.05764.
  * [198] P. Agrawal, F.-Y. Cyr-Racine, L. Randall and J. Scholtz,  JCAP 1708, 021 (2017), [1702.05482arXiv1702.05482arXiv:1702.05482].
  * [199] C. Boehm, P. Fayet and R. Schaeffer,  Phys. Lett. B518, 8 (2001), [astro-ph/0012504arXivastro-ph/0012504arXiv:astro-ph/0012504].
  * [200] C. Boehm, A. Riazuelo, S. H. Hansen and R. Schaeffer,  Phys. Rev. D 66, 083505 (2002), [astro-ph/0112522arXivastro-ph/0112522arXiv:astro-ph/0112522].
  * [201] C. Boehm and R. Schaeffer,  Astron. Astrophys. 438, 419 (2005), [astro-ph/0410591arXivastro-ph/0410591arXiv:astro-ph/0410591].
  * [202] G. Mangano, A. Melchiorri, P. Serra, A. Cooray and M. Kamionkowski,  Phys. Rev. D 74, 043517 (2006), [astro-ph/0606190arXivastro-ph/0606190arXiv:astro-ph/0606190].
  * [203] P. Serra, F. Zalamea, A. Cooray, G. Mangano and A. Melchiorri,  Phys. Rev. D 81, 043507 (2010), [0911.4411arXiv0911.4411arXiv:0911.4411].
  * [204] R. J. Wilkinson, C. Boehm and J. Lesgourgues,  JCAP 1405, 011 (2014), [1401.7597arXiv1401.7597arXiv:1401.7597].
  * [205] J. A. Schewtschenko, R. J. Wilkinson, C. M. Baugh, C. Boehm and S. Pascoli,  MNRAS449, 3587 (2015), [1412.4905arXiv1412.4905arXiv:1412.4905].
  * [206] M. Escudero, O. Mena, A. C. Vincent, R. J. Wilkinson and C. Boehm,  J. Cosmology Astropart. Phys1509, 034 (2015), [1505.06735arXiv1505.06735arXiv:1505.06735].
  * [207] K. Sigurdson and M. Kamionkowski,  Phys. Rev. Lett. 92, 171302 (2004), [astro-ph/0311486arXivastro-ph/0311486arXiv:astro-ph/0311486].
  * [208] K. Sigurdson, M. Doran, A. Kurylov, R. R. Caldwell and M. Kamionkowski,  Phys. Rev. D70, 083501 (2004), [astro-ph/0406355arXivastro-ph/0406355arXiv:astro-ph/0406355],  [Erratum: Phys. Rev.D73,089903(2006)].
  * [209] R. J. Wilkinson, J. Lesgourgues and C. Boehm,  J. Cosmology Astropart. Phys1404, 026 (2014), [1309.7588arXiv1309.7588arXiv:1309.7588].
  * [210] A. Kamada, N. Yoshida, K. Kohri and T. Takahashi,  JCAP 1303, 008 (2013), [1301.2744arXiv1301.2744arXiv:1301.2744].
  * [211] C. Boehm, J. A. Schewtschenko, R. J. Wilkinson, C. M. Baugh and S. Pascoli,  Mon. Not. Roy. Astron. Soc. 445, L31 (2014), [1404.7012arXiv1404.7012arXiv:1404.7012].
  * [212] A. Kamada, K. Kohri, T. Takahashi and N. Yoshida,  Phys. Rev. D95, 023502 (2017), [1604.07926arXiv1604.07926arXiv:1604.07926].
  * [213] X.-l. Chen, S. Hannestad and R. J. Scherrer,  Phys. Rev. D 65, 123515 (2002), [astro-ph/0202496arXivastro-ph/0202496arXiv:astro-ph/0202496].
  * [214] C. Dvorkin, K. Blum and M. Kamionkowski,  Phys. Rev. D 89, 023519 (2014), [1311.2937arXiv1311.2937arXiv:1311.2937].
  * [215] Y. Ali-Haïmoud, J. Chluba and M. Kamionkowski,  Phys. Rev. Lett. 115, 071304 (2015), [1506.04745arXiv1506.04745arXiv:1506.04745].
  * [216] V. Gluscevic and K. K. Boddy,  Phys. Rev. Lett. 121, 081301 (2018), [1712.07133arXiv1712.07133arXiv:1712.07133].
  * [217] K. K. Boddy and V. Gluscevic,  Phys. Rev. D 98, 083510 (2018), [1801.08609arXiv1801.08609arXiv:1801.08609].
  * [218] W. L. Xu, C. Dvorkin and A. Chael,  Phys. Rev. D 97, 103530 (2018), [1802.06788arXiv1802.06788arXiv:1802.06788].
  * [219] T. R. Slatyer and C.-L. Wu,  Phys. Rev. D 98, 023013 (2018), [1803.09734arXiv1803.09734arXiv:1803.09734].
  * [220] K. K. Boddy et al.,  Phys. Rev. D 98, 123506 (2018), [1808.00001arXiv1808.00001arXiv:1808.00001].
  * [221] E. O. Nadler, V. Gluscevic, K. K. Boddy and R. H. Wechsler,  Astrophys. J. 878, L32 (2019), [1904.10000arXiv1904.10000arXiv:1904.10000],  [Astrophys. J. Lett.878,32(2019)].
  * [222] K. Maamari, V. Gluscevic, K. K. Boddy, E. O. Nadler and R. H. Wechsler,  Astrophys. J. Lett. 907, L46 (2021), [2010.02936arXiv2010.02936arXiv:2010.02936].
  * [223] Y. Ali-Haïmoud,  Phys. Rev. D 103, 043541 (2021), [2101.04070arXiv2101.04070arXiv:2101.04070].