# Hitting the Thermal Target for Leptophilic Dark Matter

Cari Cesarotti  0000-0001-5128-7919<EMAIL_ADDRESS>Center for Theoretical
Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA
Gordan Krnjaic  0000-0001-7420-9577<EMAIL_ADDRESS>Theoretical Physics
Department, Fermi National Accelerator Laboratory, Batavia, Illinois 60510
Department of Astronomy and Astrophysics, University of Chicago, Chicago, IL
60637  Kavli Institute for Cosmological Physics, University of Chicago,
Chicago, IL 60637

###### Abstract

We study future lepton collider prospects for testing predictive models of
leptophilic dark matter candidates with a thermal origin. We calculate
experimental milestones for testing the parameter space compatible with
freeze-out and the associated collider signals at past, present, and future
facilities. This analysis places new limits on such models by leveraging the
utility of lepton colliders. At $e^{+}e^{-}$ machines, we make projections
using precision $Z$-pole observables from
$e^{+}e^{-}\to\ell^{+}\ell^{-}+\not{E}$ signatures at LEP and future
projections for FCC-ee in these channels. Additionally, a muon collider could
also probe new thermal relic parameter space in this scenario via
$\mu^{+}\mu^{-}\to X+\not{E}$ where $X$ is any easy identifiable SM object.
Collectively, these processes can probe much all of the parameter space for
which DM direct annihilation to $\ell^{+}\ell^{-}$ yields the observed relic
density in Higgs-like models with mass-proportional couplings to charged
leptons.

††preprint: FERMILAB-PUB-24-0124-T††preprint: MIT-CTP/5703

##  I Introduction

Thermal freeze-out is the only dark matter (DM) production mechanism that is
insensitive to unknown high-energy physics above the scale of the DM mass.
Over the past several decades, this mechanism has primarily been associated
with the Weakly Interacting Massive Particle (WIMP) paradigm, which predicts
the DM mass to be roughly within 10 GeV – 10 TeV window [1]. While there is a
rich experimental program to search for WIMPs, no verified signal has yet been
identified, which motivates broadening the search effort outside this mass
window and to a wider range of interactions.

It is well known that high energy colliders are powerful probes of thermal
dark matter; indeed, given the interactions strengths that such machines can
probe, nearly every DM candidate testable at colliders was in chemical
equilibrium in the early universe. There are currently strong limits on
various dark matter models based on existing LEP [2], Tevatron [3], and LHC
[4] data. There are also multiple studies of DM sensitivity at proposed future
facilities including the ILC [5], FCC-ee [6], FCC-hh [7], CLIC [8], CCC [9],
CEPC [10], HL-LHC [11], and various muon collider run options [12, 13, 14, 15,
16].

Figure 1: Experimental milestone for $\chi\chi\to\ell^{+}\ell^{-}$ thermal
freeze-out (black curve) plotted alongside gray-shaded constraints from the
muon anomalous magnetic moment $a_{\mu}$ [17], electron recoil direct
detection limits from DAMIC-M [18] and XENON-1T [19, 20], existing collider
limits from LEP [21] and Belle II [22]. The green band indicates the favored
parameter space for resolving the $(g-2)_{\mu}$ anomaly. The four plots show
various parameter choices to illustrate a range of phenomenology: the ratio of
$m_{\varphi}/m_{\chi}=3,10$ and the coupling within the dark sector
$\alpha_{\chi}$ = 0.5, $10^{-2}$. The dashed curves are explained throughout
Section III.

In many of these studies, the most promising DM discovery channel is missing
energy in association with a single SM object $X$ (typically a single jet or
gauge boson). However, if DM couples preferentially to heavy flavor, there may
be additional signals that improve upon the mono-$X$ strategy, offering more
handles for both signal detection and background rejection. For example, DM
that couples mainly to $b$-quarks, motivated by the galactic center excess,
can be probed using $pp\to b\bar{b}\not{E}$ production [23, 24]. This strategy
can readily be adapted to leptophilic DM that couples primarily to muons and
taus. Production of DM is therefore accompanied by additional charged leptons,
whose presence may enhance signal sensitivity.

In this paper we consider the future collider discovery potential for models
of leptophilic dark matter with mass-proportional couplings [25], to Standard
Model (SM) leptons. In the predictive parameter space for this class of
models, the scalar particle $\varphi$ that mediates dark-visible interactions
is heavier than the DM $\chi$ ($m_{\varphi}>m_{\chi})$, so that, during
freezeout, the process

$\displaystyle\chi\chi\to\varphi^{*}\to\ell^{+}\ell^{-}$ (1)

directly annihilates DM into SM leptons to yield the observed relic density
$\Omega_{\chi}h^{2}=0.120\pm 0.001$ [26]. Since the cross section for this
process depends on the coupling between $\varphi$ and charged leptons, there
is a minimum value required for viable freeze-out. Thus, a thermal origin
through the direct annihilation in Eq. (1) predicts a minimum signal yield for
DM collider production via

$\displaystyle
e^{+}e^{-}\,,\,\mu^{+}\mu^{-}\to\ell^{+}\ell^{-}\varphi~{}~{},~{}~{}$ (2)

and there is a corresponding experimental milestone for decisively testing
models with this annihilation topology. For $m_{\varphi}>2m_{\chi}$, the
mediator decays predominantly to dark matter, so the processes in Eq. (2)
complement traditional mono-$X$ searches

$\displaystyle e^{+}e^{-}\,,\,\mu^{+}\mu^{-}\to\,X\varphi\,$ (3)

which are also present in this part of parameter space and included in our
analysis of muon collider reach. Our study finds that lepton colliders with
$e^{+}e^{-}$ or $\mu^{+}\mu^{-}$ beams are powerful probes this class of
models and can test nearly all the remaining parameter space that is
compatible with a thermal relic origin.

The outline of this work is as follows. In Section II we give a theoretical
overview for this scenario and calculate the relic abundance, in Section III
we present limits and future projections for various lepton collider searches,
in Section IV we review existing bounds from previous experiments, and finally
in Section V we offer some concluding remarks.

##  II Theory

###  II.1 Model Description

Our representative benchmark scenario consists of a Majorana fermion DM
candidate $\chi$ with mass $m_{\chi}$ coupled to a CP-even scalar mediator
$\varphi$ with mass $m_{\varphi}$ via

$\displaystyle{\cal L}_{\rm
int}=-\frac{g_{\chi}}{2}\varphi\chi\chi-\varphi\\!\sum_{\ell=e,\mu,\tau}\\!g_{\ell}\,\bar{\ell}\ell~{}~{},~{}~{}g_{\ell}=g_{e}\left(\frac{m_{\ell}}{m_{e}}\right),$
(4)

where $g_{\chi}$ and $g_{\ell}$ are respectively the Yukawa couplings to DM
and charged leptons $\ell$. Note that the mass proportionality in $g_{\ell}$
can naturally arise in a Type-III two-Higgs doublet model supplemented with
our singlet scalar $\varphi$ if the latter mass-mixes primarily with the heavy
CP-even eigenstate $H$ and has suppressed mixing with the SM-like 125 GeV
eigenstate $h$ (see Appendix A for a discussion and Ref. [27] for a review of
two Higgs doublet models.).

In order for thermal freeze-out to be predictive and testable with laboratory
measurements, the cross section for DM annihilation must depend on the
$\varphi$ coupling to the SM, i.e. annihilation cannot occur strictly in the
dark sector. This requirement is readily achieved by demanding
$m_{\varphi}>m_{\chi}$, so that the $\chi\chi\to\ell^{+}\ell^{-}$ topology is
kinematically allowed for all $\ell$ satisfying $m_{\chi}>m_{\ell}$. By
contrast, if $m_{\chi}>m_{\varphi}$, the dominant channel for annihilation
will be secluded annihilation via $\chi\chi\to\varphi\varphi$, which only
depends on $g_{\chi}$ and there is no minimum requirement on $g_{\ell}$ so
long as $\chi$ thermalizes with the SM in the early universe. Note that this
secluded channel yields an $s$-wave annihilation cross section, which is
excluded for thermal relics below $\sim$ 20 GeV (see Sec. IV below). Thus, for
the reminder of this work, we only consider the case where $\varphi$ is
heavier than $\chi$.

Furthermore, since experimental limits on $g_{\ell}$ are generically stronger
than those on $g_{\chi}$, we will work in the regime where, in the
$m_{\varphi}>2m_{\chi}$ regime, the scalar will decay predominantly via
$\varphi\to\chi\chi$. Thus, we are interested in signatures involving heavy
lepton production in association with missing energy rather than $\varphi$
decaying back to SM particles as multi-lepton final states are well
constrained at colliders.

We also note that similarly predictive models with Higgs-mixed mediators,
whose couplings to all SM fermions is proportional to their massses [28].
However, in the parameter space where DM annihilates directly to SM particles,
this model faces stringent limits from various laboratory measurements (most
notably rare meson decays) which robustly exclude thermal freeze-out for all
DM masses. However, since most of these limits are based on the large top
Yukawa coupling, our leptophilic scenario remains viable across a wide range
of model parameters.

###  II.2 Cosmological Production

In the early universe for $T\gg m_{\chi}$, we assume the DM is in thermal
equilibrium with the SM. As the universe cools below $T\sim m_{\chi}$, the DM
population becomes Boltzmann suppressed and eventually the abundance of DM
freezes out. We calculate the $\chi$ relic density by solving the Boltzmann
equation

$\displaystyle\dot{n}_{\chi}+3Hn_{\chi}=-\langle\sigma
v\rangle[n_{\chi}^{2}-(n_{\chi}^{\rm eq})^{2}],$ (5)

where $H=1.66\sqrt{g_{\star}}T^{2}/M_{\rm Pl}$ is approximately the Hubble
expansion rate, $g_{\star}$ is the effective number of relativistic species,
$M_{\rm Pl}$ is the Planck mass, and $n_{\chi}^{\rm eq}$ is the equilibrium
$\chi$ number density. For $m_{\varphi}>m_{\chi}$, the thermally averaged
annihilation cross section satisfies

$\displaystyle\langle\sigma v\rangle=\sum_{\ell=e,\mu,\tau}\langle\sigma
v\rangle_{\chi\chi\to\ell\ell},$ (6)

where the sum is over all kinematically open channels ($m_{\chi}>m_{\ell}$).

To obtain the thermal average in Eq. (6), we numerically calculate the total
annihilation cross section as a function of Mandelstam variables and follow
the prescription from Ref. [29]:

$\displaystyle\langle\sigma
v\rangle_{\chi\chi\to\ell\ell}=\frac{1}{N}\int_{4m_{\chi}^{2}}^{\infty}\\!ds\,\sigma(s)(s-4m_{\chi}^{2})\sqrt{s}K_{1}\\!\left(\\!\frac{\sqrt{s}}{T}\right)\\!,\,$
(7)

where $N=8m^{4}_{\chi}TK^{2}_{2}\left(\frac{m_{\chi}}{T}\right)$ is a
normalization factor, and $K_{1,2}$ are modified Bessel functions of the first
and second kinds. For each $\chi\chi\to\ell^{+}\ell^{-}$ channel, the total
annihilation cross section is

$\displaystyle\sigma(s)=\frac{g_{\ell}^{2}g_{\chi}^{2}(s-4m_{\chi}^{2})^{3/2}}{16\pi[(s-m_{\varphi}^{2})^{2}+m_{\varphi}^{2}\Gamma_{\varphi}^{2}]}\sqrt{1-\frac{4m_{\chi}^{2}}{s}}~{},$
(8)

where $\Gamma_{\varphi}$ is the total $\varphi$ decay width, and we neglect
subleading corrections of order $m^{2}_{\ell}/s$.

Although our numerical results utilize the full thermally averaged cross
section in Eq. (7), it is useful to consider the form of $\sigma v$ in the
non-relativistic limit prior to performing the thermal average. In this
regime, each channel contributes

$\displaystyle\sigma
v_{\chi\chi\to\ell\ell}\approx\frac{g_{\chi}^{2}g_{\ell}^{2}m_{\chi}^{2}v^{2}}{8\pi(m_{\varphi}^{2}-4m_{\chi}^{2})^{2}}\propto
g_{\chi}^{2}g_{\ell}^{2}\left(\frac{m_{\chi}}{m_{\varphi}}\right)^{4}\frac{1}{m_{\chi}^{2}},$
(9)

where $v$ is the relative velocity between initial state $\chi$ particles.
Matching the convention in Refs. [30, 31], we define the dimensionless
variable

$\displaystyle
y\equiv\epsilon^{2}\alpha_{\chi}\left(\frac{m_{\chi}}{m_{\varphi}}\right)^{4},~{}~{}\epsilon\equiv
g_{e}/e~{},$ (10)

where $\alpha_{\chi}\equiv g_{\chi}^{2}/(4\pi)$. We normalize the parameter
$\epsilon$ in units of the electron charge $e$ in analogy with the familiar
kinetic-mixing parameter in models with dark photon mediators [32]. This
quantity is useful because, away from the fine tuned resonance region
$m_{\varphi}\approx 2m_{\chi}$, there is a one-to-one relationship between
$\sigma v$ and $y$ independently of the $m_{\varphi}/m_{\chi}$ or
$g_{\chi}/g_{\ell}$ ratios; here the relic density for each value of
$m_{\chi}$ is simply set by the product that defines $y$. In Fig. 1, the black
curves in the $y$ vs. $m_{\chi}$ plane represent the parameter points for
which $\chi\chi\to\ell\ell$ annihilation yields the observed dark matter
abundance through thermal freeze-out.

Importantly, the non-relativistic expansion in Eq. (9) also shows that this
process is $p$-wave, so the annihilation rate falls precipitously as the
universe expands and the temperature cools. As discussed below in Section IV,
this feature sharply suppresses the annihilation to SM particles during and
after the epoch of recombination.

##  III Collider Searches

In this section we present the lepton collider signatures that the model in
Section II predicts. Here we put new limits on this scenario using existing
data from LEP and Belle II, and we make future predictions for FCC-ee and a
future muon collier. The results are summarized in Fig. 1, with the excluded
regions shaded gray and predictions for future experiments given in dashed
lines.

###  III.1 Belle II

Since $\varphi$ has sizeable couplings to heavy flavor, it can be produced as
final state radiation at experiments where 2nd and 3rd generation dilepton
pairs are copiously produced, such as at $B-$factories via

$\displaystyle e^{+}e^{-}\to\mu^{+}\mu^{-}\varphi,~{}~{}(\varphi\to\chi\chi),$
(11)

where the dark matter $\chi$ in the final state produces missing energy. This
signal has an irreducible background from
$e^{+}e^{-}\to\mu^{+}\mu^{-}\bar{\nu}\nu$ from the missing energy in the final
state neutrinos, as well as reducible backgrounds from double bremsstrahlung
$e^{+}e^{-}\to\mu^{+}\mu^{-}\gamma\gamma$. The Belle II collaboration has
placed constraints on an $L_{\mu}-L_{\tau}$ gauge boson through this search
channel [22]. We recast these limits for our scalar scenario by using
MadGraph5_aMC@NLO [33] to calculate the cross section for
$e^{+}e^{-}\to\mu^{+}\mu^{-}\varphi$ production and rescale the limits for the
$L_{\mu}-L_{\tau}$ vector model with the corresponding cross section for
vector particles. We neglect any signal acceptance differences between scalars
and vector and vector production.

We note that, in principle, $B$-factories can be sensitive to $\varphi$
production as final state radiation in $e^{+}e^{-}\to\tau^{+}\tau^{-}\varphi$
channels (considered below for higher energy lepton colliers). However,
despite the large $\tau$ enhancement in this channel, to date, there has been
no dedicated search for such a process at these facilities for $\varphi$
decaying invisibly111 However, there are strong BABAR limits on
$e^{+}e^{-}\to\tau^{+}\tau^{-}$ followed by a $\varphi\to\mu^{+}\mu^{-}$ decay
[34]. We leave a detailed study of such signals for future work.

Figure 2: Representative Feynman diagram of a rare
$Z\to\tau^{+}\tau^{-}\varphi$ decay that generates our signal of interest at
LEP and FCC-ee. Here the $\varphi\to\chi\chi$ produces additional missing
energy relative to SM $Z\to\tau^{+}\tau^{-}$ decays.

###  III.2 LEP

It has long been established that the LEP and LEP-2 $e^{+}e^{-}$ colliders at
CERN [21] have impressive sensitivity to dark matter production through
mono-$\gamma$ and missing energy $\not{E}$ signatures [2]. While these
canonical signatures can also arise in our scenario from the couplings in Eq.
(4), they are not the dominant channel as the production cross section for
$e^{+}e^{-}\to\gamma\varphi\to\gamma\not{E}$ is suppressed by the small
electron coupling $g^{2}_{e}$. Since the $\varphi$ coupling to the $\tau$ is
enhanced by a factor of $(m_{\tau}/m_{e})^{2}\sim 10^{6}$, here we instead
consider the signature

$\displaystyle
e^{+}e^{-}\to\tau^{+}\tau^{-}\varphi~{}~{},~{}~{}(\varphi\to\chi\chi),$ (12)

where the invisibly decaying mediator is produced as final state radiation and
predicts deviations within well-measured SM uncertainties. Since SM $\tau$
decays always produce missing energy via $\nu$, our signal process will be
constrained by any previous measurement of $Z$ decays that result in $\tau$
production [35]. Up to small differences in geometric acceptance and final
state $\tau$ kinematic distributions, $Z\to\tau^{+}\tau^{-}$ production is
experimentally indistinguishable from $Z\to\tau^{+}\tau^{-}\varphi$ for
$\varphi$ decaying invisibly.

The one of the primary motivations for a high-energy electron-positron
collider is a run at the $Z$ pole where $\sqrt{s}\approx m_{Z}$, so the $Z$
production cross section is resonantly enhanced and with copious production,
precision measurements can be made. Therefore, LEP is mainly sensitive to
$\varphi$ production through modifications to the well-measured properties of
the $Z$. Of particular interest is the signal process depicted in Fig. 2

$\displaystyle Z\to\tau^{+}\tau^{-}\varphi,~{}(\varphi\to\chi\chi)~{},$ (13)

where the $\varphi$ decay to dark matter yields additional missing energy in
the final state. In Appendix B we calculate the partial width for this 3-body
process. From LEP, we obtain a limit on $\varphi$ production from the
measurement of the di-tau partial width [36]

$\displaystyle\Gamma({Z\to\tau^{+}\tau^{-}})=84.08\pm 0.22\,\rm MeV.$ (14)

From the uncertainty we extract a 2$\sigma$ limit on the
BR($Z\to\tau^{+}\tau^{-}\varphi$), we find the value of $g_{\tau}$ at each
value of $m_{\varphi}$ that satisfies

$\displaystyle\Gamma(Z\to\tau^{+}\tau^{-}\varphi)=2(0.22\,\rm MeV)~{},$ (15)

and translate this into a constraint on $y$ for a given choice of
$\alpha_{\chi}$ and $m_{\chi}/m_{\varphi}$ ratio.

We also include a highly simplified estimate for the acceptance of
$Z\rightarrow\tau\tau\varphi$ events due to analysis cuts on kinematics, which
is discussed in Appendix C. In Fig. 1 this bound is shaded in gray and labeled
LEP.

###  III.3 Future Circular Collider

The Future Circular Collider (FCC) at CERN is a proposed high-energy collider
to succeed the LHC at CERN. The first stage of the project is the FCC-ee, an
$e^{+}e^{-}$ collider to run at similar run configurations of LEP but with
much greater luminosity, which aims to take data in the early 2040s [37, 6,
38, 39, 40]. As part of its overall physics program, FCC-ee will run at the
$Z$-pole and produce over one trillion $Z$ bosons for precision electroweak
measurements (called the Tera-Z run), thereby exceeding the previous LEP
record of approximately 17 million [21]. This allows for a diverse physics
program that involves precision measurements of $Z$ decays (e.g. Refs. [41,
42])

To make sensitivity projects of the FCC-ee at the $Z$ pole, we assume
statistical uncertainties dominate the FCC-ee error budget and rescale the
uncertainty on the $Z$ width measurement by
$\sqrt{N_{\text{LEP}}/N_{\text{FCC-ee}}}$ for $N_{X}$ the total number of $Z$
bosons produced. With this estimation we predict sensitivities to $Z$-boson
branching ratios comparable to what is stated in Ref. [37]. In Fig. 1 we show
projections of the improvement of FCC-ee over LEP with the red dashed curve.

###  III.4 Future Muon Collider

Another possible option for a future high-energy collider is a multi-TeV muon
collider (MuC) [43, 44, 45, 46, 47, 48].

Due to their moderate mass and fundamental nature, muons can theoretically be
accelerated to collide at $\sqrt{\hat{s}}\sim 10\text{ TeV}$ energies without
insurmountable loses due to synchroton radiation. Additionally, as muons
couple predominantly to electroweak gauge bosons and do not produce enormous
hadronic backgrounds, a high-energy muon collider could provide a window into
both the 10-TeV frontier as well as precision electroweak measurements within
the next 30 years [49, 12, 50, 51, 52, 53, 54, 55]. While this technology
needs a full demonstrator program to prove feasibility, consider theoretical,
experimental, and accelerator work has already gone into its design (see Ref.
[56] for a recent summary), and its construction is already partially funded.

As a muon collider is still in its very early stages of planning, the run
configurations have not yet been fixed and a multitude of possible center-of-
mass energies are considered. The guiding principle for projections is that
the total luminosity must scale with $s$ such that we maintain sensitivity to
$s$-channel processes at high energies. However, due to the fact that muons
decay, higher-energy muon beams will have a higher luminosity from time
dilation effects.

In this paper, we utilize the high energy of the muon collider to probe
complementary parameter space to what can be done at the FCC-ee. We consider
the process $\mu^{+}\mu^{-}\to\gamma\varphi\to\gamma\not{E}$ with a
mono-$\gamma$ signature. Because of the kinematics of the process, the
$\gamma$ will have energy roughly $E_{\gamma}\sim\sqrt{s}/2$ depending on the
width of the $\varphi$ and energy resolution $\sigma_{E}$. The primary SM
background for this process is $\mu^{+}\mu^{-}\to\gamma\bar{\nu}\nu$.222We
also considered the 3$\gamma$ background, but as two photons, one with $E\sim$
TeV have to be missed, this background is comparatively negligible. We
simulate both the signal and background events in MadGraph5_aMC@NLO [33]. To
mitigate this background, we utilize the fully-constrained kinematics of the
2-body final state of the signal and cut on the energy of the visible photon:

$E_{\gamma}\in\left(E_{0}(1-2\sigma_{E}),E_{0}(1+2\sigma_{E})\right)$ (16)

where $E_{0}\equiv\frac{E_{\rm CM}^{2}-m_{\varphi}^{2}}{2E_{\rm CM}}$ and
$\sigma_{E}$ is the energy resolution we fix nominally at $\sigma_{E}=3$%. We
also take a conservative photon reconstruction efficiency of 98%, although we
expect this number for TeV photons might be higher. The search strategy in
this channel is generic as long as the mediator decays invisibly via
$\varphi\to\chi\chi(\not{E})$ and the other SM particles $X$ recoil in a
diametrically opposed direction to the missing energy and momentum in the
final state. A similar strategy with partially visible decays was considered
in Ref. [57].

The dashed blue curve in Fig. 1 is a sensitivity projection for
$\mu^{+}\mu^{-}\to X\varphi\to X\not{E}$ at a 3 TeV muon collider concept with
1 ab-1 of integrated luminosity.333We do not consider a lower energy center of
mass, like the Higgs pole, as higher-energy configurations have a broader
physics program and can even produce upwards of $10^{6}$ to $10^{8}$ Higgs.

##  IV Other Searches

In this section we describe the other constraints on this model from precision
measurements, cosmology, and direct detection experiments.

###  IV.1 Muon Anomalous Magnetic Moment

The Brookhaven E821 [58] and Fermilab Muon $g-2$ [17] experiments have
measured the anomalous magnetic moment of the muon, resulting in a world
average of

$\displaystyle a_{\mu}^{\rm exp}=116592059(22)\times 10^{-11}~{}~{}.$ (17)

The SM prediction based on R-ratio calculations is [59]

$\displaystyle a_{\mu}^{\rm SM}=116591810(43)\times 10^{-11}~{}~{},$ (18)

which differs from the observed value by

$\displaystyle\Delta a_{\mu}=(249\pm 48)\times 10^{-11},$ (19)

and thereby exceeds the nominal 5$\sigma$ benchmark for a new physics
discovery.

However, a lattice calculation from the BMW collaboration extracts a SM
prediction for $a_{\mu}$ that is closer to the observed value and in tension
with the R-ratio result. Additional lattice results will be necessary to
independently verify the BMW result and future R-ratio data will clarify the
tension between BABAR, KLOE, and CMD-3 measurements near $\sqrt{s}\approx$ 1
GeV [60, 59, 17]. Thus, in our analysis, we remain agnostic about whether the
discrepancy is due to new physics.

Since the leptophilic mediator $\varphi$ couples to muons, at one loop it
modifies $a_{\mu}$ by

$\displaystyle\Delta
a_{\mu}=\frac{g_{\mu}^{2}}{16\pi^{2}}\int_{0}^{1}dz\frac{(1-z)(1-z^{2})}{(1-z)^{2}+r^{2}z},$
(20)

where $r=m_{\varphi}/m_{\mu}$, which constrains the thermal relic parameter
space shown in Fig. 1. The gray region labelled “$\Delta a_{\mu}>5\sigma$” is
where model parameters make the tension with observation worse by 5 standard
deviations. For parameters within the green band labeled “$a_{\mu}$ favored,”
the $\varphi$ coupling to muons accounts for the discrepancy between the
measured world average and the $R$-ratio prediction within the Standard Model
assuming the difference is due to new physics.

###  IV.2 Fixed Target Experiments

For the sub-GeV DM mass range, proposed fixed target accelerators are powerful
laboratory probes of light weakly coupled particles [61]. The mediator
$\varphi$ can be radiatively produced in muon-nucleus fixed target
interactions via

$\displaystyle\mu N\to\mu N\varphi,~{}~{}(\varphi\to\chi\chi)~{},$ (21)

where $N$ is a target nucleus. For $m_{\varphi}>m_{\mu}$, the radiated
$\varphi$ typically acquires a large fraction of the incident beam energy,
which serves as a trigger for a dark matter search if other visible particles
produced in $\mu$-$N$ scattering can be vetoed. In Fig. 1 we show projections
for the proposed $M^{3}$ experiment [62, 57] based on this search strategy.
Comparable projections could also be computed for the CERN NA64$\mu$
experiment which performs a similar search at higher energy [63], but
translating the publicly available projections into the scalar model
considered here is beyond the scope of this work.

###  IV.3 Direct Detection

####  IV.3.1 Electron Recoils

Since the mediator $\varphi$ couples to electrons, this scenario is
constrained by direct detection experiments sensitive to scattering off
electrons. The non-relativistic cross section for $\chi e\to\chi e$ scattering
is

$\displaystyle\sigma_{\chi e}=\frac{g_{\chi}^{2}g_{e}^{2}\mu^{2}_{\chi
e}}{\pi[(\alpha
m_{e})^{2}+m_{\varphi}^{2}]^{2}}\to\frac{g_{e}^{2}g_{\chi}^{2}\mu_{e\chi}^{2}}{\pi
m_{\varphi}^{4}}\approx\left(\frac{4e^{2}m_{e}^{2}}{m_{\chi}^{4}}\right)y,~{}~{}$
(22)

where $\mu_{\chi e}$ is the $\chi$-$e$ reduced mass and in the last step we
took the the $m_{e}\ll m_{\varphi},m_{\chi}$ limit, where $\sigma_{\chi
e}\propto y$ and the direct detection cross section has the same parametric
dependence as the $\chi\chi\to e^{+}e^{-}$ annihilation cross section from Eq.
(9). Although the electron scattering sensitivities various experiments are
comparable in this mass range, in Fig. 1 we show the current best limits from
the DAMIC-M [18] and XENON-1T [19, 20] experiments.

####  IV.3.2 Nuclear Recoils

Although $\varphi$ has no tree-level coupling to quarks, at one-loop there is
induced mass-mixing between $\varphi$ and the SM Higgs, which yields the
interactions

$\displaystyle{\cal L}_{\rm int}\supset\theta y_{q}\varphi\bar{q}q\equiv
g_{q}\varphi\bar{q}q,$ (23)

where $q$ is any SM quark, $y_{q}=m_{q}/v_{h}$ is its SM Yukawa coupling to
the Higgs field, and $\theta$ is the $\varphi$-$h$ mixing angle. The leading
contribution to this mixing is given by integrating out the $\tau$ loop
$\theta\sim g_{\tau}y_{\tau}/(16\pi^{2})$. Based on this mixing, the elastic
spin-independent direct-detection cross section for $\chi$-$n$ scattering is

$\displaystyle\sigma_{\chi n}=\frac{g_{\chi}^{2}g^{2}_{n}\mu^{2}_{\chi n}}{\pi
m_{\varphi}^{4}},$ (24)

where $n$ is a target nucleon, $\mu_{\chi n}$ is the DM-nucleon reduced mass,
and the induced $\varphi$ nucleon coupling $g_{n}\propto\theta$ can be written
in terms of the mediator couplings to quarks in Eq. (23) following the
procedure in Ref. [64] – see also the appendix of Ref. [65] for a modern
treatment. Writing Eq. (24) in terms of $y$ from Eq. (10)

$\displaystyle y\approx 10^{-8}\left(\frac{\sigma_{\chi n}}{10^{-46}\,\rm
cm^{2}}\right)\left(\frac{m_{\chi}}{20\,\rm GeV}\right)~{}.$ (25)

Thus, for cross sections of order $\sigma_{\chi n}\approx 10^{-46}$ cm2 at the
maximum sensitivity of XENONnT [66] near $m_{\chi}\approx 20$ GeV, direct
detection does not meaningfully constrain the thermal relic values of $y$ as
shown in the black curve of Fig. 1. A similar argument can be used to show
that the two-loop QED scattering processes lead to similar results. Thus, we
omit nuclear recoil limits through loop induced elastic scattering in Fig. 1.
However, we note that if the Higgs-mixing parameter $\theta$ is enhanced
beyond the $\tau$ loop that we use as a reference value (e.g. from additional
tree level contributions), the model dependent constraint from this process
could be competitive.

###  IV.4 Cosmology

####  IV.4.1 BBN Dark Radiation

For $m_{\chi}\lesssim$ MeV, the $\chi$ freeze-out occurs around the epoch of
Big Bang Nucleosynthesis (BBN) and can modify light chemical yields relative
to SM predictions. In this light mass range, the DM has a large thermal
abundance at the onset of BBN, which modifies the Hubble expansion rate and
affects the timing of BBN events [67, 68]. Using the latest nuclear cross
sections and data on the deuterium fraction, BBN observables require
$m_{\chi}\gtrsim$ 5 MeV for Majorana DM candidates that transfer entropy to
photons upon freeze-out [69] and this defines the leftmost regions in the
panels of Fig. 1

####  IV.4.2 CMB Energy Injection

Thermal DM below the GeV scale faces stringent bounds from CMB energy
injection from out-of-equilibrium DM annihilation just after recombination
[70]. For an $s$-wave annihilation cross section, thermal WIMPs are excluded
below $\sim$ 20 GeV if their annihilation yields visible energy into the CMB
[26]. However, in our scenario the annihilation cross section in Eq. (9) is
$p$-wave and therefore suppressed by $v^{2}$ in the non-relativistic regime.
Thus, the CMB limit from Planck does not constrain the thermal relic parameter
space for this model.

##  V Discussion

In this paper we have studied a predictive leptophilic DM scenario in which
the thermal relic abundance arises from direct annihilation to SM particles
via $\chi\chi\to\ell^{+}\ell^{-}$ with mass-weighted couplings to charged
leptons. Since this scenario has no tree-level couplings to quarks, LHC and
nuclear-recoil direct detection signals are suppressed compared to traditional
WIMP DM models. We have set new limits on this model based on various probes,
including measurements of the muon anomalous magnetic moment, electron recoil
direct detection, $B$-factory production, and $Z$-pole measurements at LEP.
Collectively, these limits leave much of the thermal relic parameter space
unexplored, particularly for DM masses above the GeV scale where the relic
density arises predominantly through $\chi\chi\to\tau^{+}\tau^{-}$
annihilation in the early universe.

Future accelerators will play a key role in probing the remaining thermal
relic parameter space in this scenario. In the sub-GeV mass range, the
proposed fixed-target missing-momentum experiment $M^{3}$ can comprehensively
cover the full thermal target and complement ongoing efforts currently
underway with the NA64$\mu$ experiment. At higher mass, the proposed FCC-ee
$e^{+}e^{-}$ collider can probe large regions of viable parameter space with a
Tera-$Z$ run at the $Z$-pole, which can improve sensitivity to
$Z\to\tau^{+}\tau^{-}$ decays, currently limited by the luminosity of the
legacy LEP data set. And at the largest masses compatible with a viable
thermal history $m_{\chi}\sim$ TeV in this framework, a future TeV-scale muon
collider can probe nearly all of the remaining freeze-out parameter space with
mono-$\gamma/Z$ searches.

## Acknowledgements

We would like to thank André de Gouvêa, Rodolfo Capdevilla, Bertrand Echenard,
Maxim Pospelov, Jesse Thaler, Dario Buttazzo, Simone Pagan Griso, and Andrea
Wulzer for helpful conversations. C.C. is supported by the U.S. Department of
Energy (DOE) Office of High Energy Physics under Grant Contract No. DE-
SC0012567. Fermilab is operated by the Fermi Research Alliance, LLC under
Contract DE-AC02-07CH11359 with the U.S. Department of Energy. This material
is based partly on support from the Kavli Institute for Cosmological Physics
at the University of Chicago through an endowment from the Kavli Foundation
and its founder Fred Kavli.

##  Appendix A Realization in a Two Higgs Doublet Model

The mass weighted couplings in our scenario can arise in a type-III two Higgs
(2HDM) doublet model with an additional CP even scalar singlet $S$ that mixes
with the SM-like $h$ and heavy Higgs state $H$. In this realization of a 2HDM,
before electroweak symmetry breaking, the Yukawa interactions can be written

$\displaystyle{\cal
L}\supset\lambda_{u}H_{1}Q\bar{u}+\lambda_{d}QH^{\dagger}_{1}Q\bar{d}+\lambda_{\ell}H_{2}L\bar{e}~{},$
(26)

where $H_{1,2}$ are electroweak doublets with quantum numbers $2_{\pm 1/2}$
under $SU(2)_{L}\times U(1)_{Y}$, so $H_{1}$ gives mass to the quarks and
$H_{2}$ gives mass to charged leptons. Here we use two-component Weyl fermion
notation for the quark and lepton fields, where $Q,L$ are electroweak doublets
and $\bar{u},\bar{d},\bar{e}$ are $SU(2)_{L}$ singlet Dirac partners for the
up-type, down-type, and charged lepton fields, respectively.

Following the conventions in Ref. [27], after electroweak symmetry breaking,
the vacuum expectation values are $\langle H_{i}\rangle=v_{i}$. We define the
ratio $\tan\beta\equiv v_{2}/v_{1}$ and diagonalize the CP even mass
eigenstates with the rotation matrix

$\displaystyle\begin{pmatrix}\sqrt{2}\,{\rm Re}(H_{1})-v_{1}\\\ \sqrt{2}\,{\rm
Re}(H_{2})-v_{2}\\\
\end{pmatrix}=\begin{pmatrix}~{}~{}\cos\alpha&\sin\alpha\\\
-\sin\alpha&\cos\alpha\\\ \end{pmatrix}\begin{pmatrix}h\\\ H\end{pmatrix},$
(27)

where $h$ is the SM-like Higgs mass eigenstate and $H$ is the heavy CP-even
state. In the mass eigenbasis, the $h,H$ couplings to fermions and electroweak
gauge bosons (denoted by $V$) are given in Table 1

| $\lambda_{\rm 2HDM}/\lambda_{\rm SM}$
---|---
$hVV$ | $\sin(\beta-\alpha)$
$hQ\bar{u}$ | $\sin(\beta-\alpha)+\cos(\beta-\alpha)/\tan\beta$
$hQ\bar{d}$ | $\sin(\beta-\alpha)+\cos(\beta-\alpha)/\tan\beta$
$hL\bar{e}$ | $\sin(\beta-\alpha)-\tan\beta\cos(\beta-\alpha)$
$HVV$ | $\cos(\beta-\alpha)$
$HQ\bar{u}$ | $\cos(\beta-\alpha)-\sin(\beta-\alpha)/\tan\beta$
$HQ\bar{d}$ | $\cos(\beta-\alpha)-\sin(\beta-\alpha)/\tan\beta$
$HL\bar{e}$ | $\cos(\beta-\alpha)+\tan\beta\sin(\beta-\alpha)$
Table 1: Tree level couplings for CP even eigenstates $h$ and $H$ in a Type
III 2HDM normalized to their values in the SM. For large $\tan\beta\gg 1$ in
the alignment limit where $\sin(\beta-\alpha)\to 1$, $h$ recovers the SM like
couplings and $H$ becomes lepophilic with mass proportional Yukawa couplings.

To realize our scenario in Eq. (4) at low energy, we add a singlet scalar $S$
and introduce super-renormalizable Higgs-portal couplings

$\displaystyle
V(H_{1},H_{2},S)=S\left(\mu_{11}H_{1}^{\dagger}H_{1}+\mu_{12}H_{1}^{\dagger}H_{2}+\mu^{*}_{12}H_{2}^{\dagger}H_{1}+\mu_{22}H^{\dagger}_{2}H_{2}\right)~{},$
(28)

where $\mu_{ij}\ll v_{h}$ are dimensionful mixing parameters. Writing this
expression in terms of mass eigenstates in the electroweak breaking vacuum, we
have

$\displaystyle{\rm
Re}(H_{1})=\frac{h\cos\alpha+H\sin\alpha+v_{1}}{\sqrt{2}}~{}~{},~{}~{}{\rm
Re}(H_{2})=\frac{-h\sin\alpha+H\cos\alpha+v_{2}}{\sqrt{2}}~{},$ (29)

so the mixed potential involving $S$ becomes

$\displaystyle V(h,H,S)$ $\displaystyle=$
$\displaystyle\frac{S}{2}\biggl{[}\mu_{11}\left(h\cos\alpha+H\sin\alpha+v_{1}\right)^{2}+2{\rm
Re(\mu_{12})}\left(h\cos\alpha+H\sin\alpha+v_{1}\right)\left(-h\sin\alpha+H\cos\alpha+v_{2}\right)$
(30) $\displaystyle\hskip
85.35826pt+\mu_{22}\left(-h\sin\alpha+H\cos\alpha+v_{2}\right)^{2}\biggr{]}~{},$

and the mass-squared mixing coefficients for $S$-$H$ and $S$-$h$ mixing can
respectively be extracted via

$\displaystyle\mu_{SH}\equiv\frac{\partial^{2}V}{\partial S\partial
H}\biggr{|}_{h,H,S=0}=\sin\alpha(v_{1}\mu_{11}+v_{2}\mu_{12})+\cos\alpha(v_{1}\mu_{12}+v_{2}\mu_{22})$
$\displaystyle\mu_{Sh}\equiv\frac{\partial^{2}V}{\partial S\partial
h}\biggr{|}_{h,H,S=0}=\cos\alpha(v_{1}\mu_{11}+v_{2}\mu_{12})-\sin\alpha(v_{1}\mu_{12}+v_{2}\mu_{22})$
(31)

For parameter choices that satisfy $\mu_{Sh}\ll\mu_{SH}\ll v_{h}$ in the
$\tan\beta\gg 1$ limits, $S$ will mix predominantly with the $H$ state whose
couplings are leptophilic and proportional to charged lepton masses.
Diagonalizing this additional mixing between $S$-$H$ thus yields the
$S\bar{\ell}\ell$ interactions in Eq. (4). Note that, as in the case of one SM
Higgs doublet, a additional singlet scalar can mix through the super-
renormalizable portal wihtout acquiring a vacuum expectation value provided
that a suitably chosen linear $S$ term is included in the potential (see [71]
for a discussion).

##  Appendix B $Z\to\ell^{+}\ell^{-}\varphi$ decays

The amplitude for radiative $\varphi$ production in
$Z\to\ell^{+}(p_{1})\ell^{-}(p_{2})\varphi(p_{3})$ decays depicted in Fig. 2
can be written

$\displaystyle{\cal
M}=\frac{g_{\tau}m_{Z}}{v_{h}}\bar{u}(p_{1})\\!\biggl{[}\frac{(\not{p_{1}}+\not{p_{3}})+m_{\tau}}{m_{13}^{2}-m_{\ell}^{2}}\gamma^{\mu}(g_{V}+g_{A}\gamma^{5})+\gamma^{\mu}(g_{V}\\!+\\!g_{A}\gamma^{5})\frac{(\not{p_{2}}+\not{p_{3}})+m_{\tau}}{m_{12}^{2}-m_{\ell}^{2}}\biggr{]}v(p_{2}),$
(32)

where $p_{i}$ are the final state momenta with
$m_{ij}^{2}\equiv(p_{i}+p_{j})^{2}$, $v_{h}=246$ GeV is the SM Higgs vacuum
expectation value, and $g_{V,A}$ are respectively the tau vector and axial-
vector charges, which satisfy

$\displaystyle
g_{V}=-\frac{1}{2}(1-4\sin^{2}\theta_{W})~{}~{},~{}~{}g_{A}=-\frac{1}{2}~{}.$
(33)

The partial width can be written

$\displaystyle\Gamma(Z\to\ell^{+}\ell^{-}\varphi)=\frac{1}{256\pi^{3}m^{3}_{Z}}\int
dm^{2}_{12}\int dm^{2}_{23}\,\overline{|{\cal M}|^{2}}~{},$ (34)

where the $m^{2}_{12}$ integration limits satisfy $(m^{2}_{12})_{\rm
min}=4m^{2}_{\ell}~{},~{}(m^{2}_{12})_{\rm max}=(m_{Z}-m_{\varphi})^{2}~{},$
and the $m^{2}_{23}$ limits are

$\displaystyle(m^{2}_{23})_{\rm min}^{\rm
max}=(E_{2}^{*}+E_{3}^{*})^{2}-\\!\biggl{(}\\!\sqrt{E^{*2}_{2}-m_{\ell}^{2}}\mp\sqrt{E^{*2}_{3}-m_{\varphi}^{2}}\biggr{)}^{\\!2},$
(35)

and we have defined

$\displaystyle
E^{*}_{2}=\frac{m_{12}}{2}~{},~{}E_{3}^{*}=\frac{m_{Z}^{2}-m_{12}^{2}-m_{\varphi}^{2}}{2m_{12}}~{}.$
(36)

##  Appendix C Estimates on Acceptance Efficiency for
$Z\rightarrow\tau\tau\varphi$

Figure 3: The truth-level $p_{T}$ and $\eta$ distributions of the $\tau$
particles produced from $Z$ decays as a function of the $\varphi$ mass.

Since$\tau$decays are too complex to simply reconstruct the Z peak, as a proxy
for the much more involved analysis we perform truth-level cuts on taus from
MadGraph-simulated events [33] to estimate the number of new physics events
that would resemble SM events close enough to be included in the decay width
measurement. As done in LEP and LHC [35, 72] analyses, we put a cut on the
transverse momentum of taus $p_{T}>15$ GeV and on the total energy squared
$s^{\prime}$ of the $\tau\tau$ system to be $s^{\prime}>4m_{\tau}^{2}$. While
a more rigorous study would include evaluating the cuts of individual $\tau$
decay channels specifically and simulating the full $\tau$ decay, comparing
truth-level $\tau$ information with and without the $\phi$ production gives us
an order-of-magnitude understanding of how many new-physics events can be
mistaken for SM events.

We show the $\tau$ distributions in $p_{T}$ and $\eta$ in Figure 3. Note that
as the mass of the $\varphi$ is increased to $\mathcal{O}(m_{Z}\gtrsim
10\text{ GeV})$, the $p_{T}$ of the taus is dramatically affected and most of
the new physics events will be rejected. However, the $\eta$ distribution of
the taus is largely unchanging.

## References

  * [1] G. Jungman, M. Kamionkowski, and K. Griest, “Supersymmetric dark matter,” Phys. Rept. 267 (1996) 195–373, arXiv:hep-ph/9506380.
  * [2] P. J. Fox, R. Harnik, J. Kopp, and Y. Tsai, “LEP Shines Light on Dark Matter,” Phys. Rev. D 84 (2011) 014028, arXiv:1103.0240 [hep-ph].
  * [3] Y. Bai, P. J. Fox, and R. Harnik, “The Tevatron at the Frontier of Dark Matter Direct Detection,” JHEP 12 (2010) 048, arXiv:1005.3797 [hep-ph].
  * [4] F. Kahlhoefer, “Review of LHC Dark Matter Searches,” Int. J. Mod. Phys. A 32 no. 13, (2017) 1730006, arXiv:1702.02430 [hep-ph].
  * [5] H. Dreiner, M. Huck, M. Krämer, D. Schmeier, and J. Tattersall, “Illuminating Dark Matter at the ILC,” Phys. Rev. D 87 no. 7, (2013) 075015, arXiv:1211.2254 [hep-ph].
  * [6] G. Bernardi et al., “The Future Circular Collider: a Summary for the US 2021 Snowmass Process,” arXiv:2203.06520 [hep-ex].
  * [7] T. Han, S. Mukhopadhyay, and X. Wang, “Electroweak Dark Matter at Future Hadron Colliders,” Phys. Rev. D 98 no. 3, (2018) 035026, arXiv:1805.00015 [hep-ph].
  * [8] CLICdp Collaboration, J.-J. Blaising, P. Roloff, A. Sailer, and U. Schnoor, “Physics performance for Dark Matter searches at $\sqrt{s}=$ 3 TeV at CLIC using mono-photons and polarised beams,” arXiv:2103.06006 [hep-ex].
  * [9] M. Bai et al., “C3: A ”Cool” Route to the Higgs Boson and Beyond,” in Snowmass 2021.  10, 2021.  arXiv:2110.15800 [hep-ex].
  * [10] Z. Liu, Y.-H. Xu, and Y. Zhang, “Probing dark matter particles at CEPC,” JHEP 06 (2019) 009, arXiv:1903.12114 [hep-ph].
  * [11] M. Cepeda et al., “Report from Working Group 2: Higgs Physics at the HL-LHC and HE-LHC,” CERN Yellow Rep. Monogr. 7 (2019) 221–584, arXiv:1902.00134 [hep-ph].
  * [12] H. Al Ali et al., “The muon Smasher’s guide,” Rept. Prog. Phys. 85 no. 8, (2022) 084201, arXiv:2103.14043 [hep-ph].
  * [13] T. Han, Z. Liu, L.-T. Wang, and X. Wang, “WIMPs at High Energy Muon Colliders,” Phys. Rev. D 103 no. 7, (2021) 075004, arXiv:2009.11287 [hep-ph].
  * [14] T. Han, Z. Liu, L.-T. Wang, and X. Wang, “WIMP Dark Matter at High Energy Muon Colliders $-$A White Paper for Snowmass 2021,” in Snowmass 2021.  3, 2022.  arXiv:2203.07351 [hep-ph].
  * [15] S. Bottaro, D. Buttazzo, M. Costa, R. Franceschini, P. Panci, D. Redigolo, and L. Vittorio, “Closing the window on WIMP Dark Matter,” Eur. Phys. J. C 82 no. 1, (2022) 31, arXiv:2107.09688 [hep-ph].
  * [16] S. Bottaro, D. Buttazzo, M. Costa, R. Franceschini, P. Panci, D. Redigolo, and L. Vittorio, “The last complex WIMPs standing,” Eur. Phys. J. C 82 no. 11, (2022) 992, arXiv:2205.04486 [hep-ph].
  * [17] Muon g-2 Collaboration, D. P. Aguillard et al., “Measurement of the Positive Muon Anomalous Magnetic Moment to 0.20 ppm,” Phys. Rev. Lett. 131 no. 16, (2023) 161802, arXiv:2308.06230 [hep-ex].
  * [18] DAMIC-M Collaboration, I. Arnquist et al., “First Constraints from DAMIC-M on Sub-GeV Dark-Matter Particles Interacting with Electrons,” Phys. Rev. Lett. 130 no. 17, (2023) 171003, arXiv:2302.02372 [hep-ex].
  * [19] XENON Collaboration, E. Aprile et al., “Light Dark Matter Search with Ionization Signals in XENON1T,” Phys. Rev. Lett. 123 no. 25, (2019) 251801, arXiv:1907.11485 [hep-ex].
  * [20] XENON Collaboration, E. Aprile et al., “Emission of single and few electrons in XENON1T and limits on light dark matter,” Phys. Rev. D 106 no. 2, (2022) 022001, arXiv:2112.12116 [hep-ex].
  * [21] J. Drees, “Review of final LEP results, or, A Tribute to LEP,” Int. J. Mod. Phys. A 17 (2002) 3259–3283, arXiv:hep-ex/0110077.
  * [22] Belle-II Collaboration, I. Adachi et al., “Search for an Invisible Z’ in a Final State with Two Muons and Missing Energy at Belle II,” Phys. Rev. Lett. 130 no. 23, (2023) 231801, arXiv:2212.03066 [hep-ex].
  * [23] T. Lin, E. W. Kolb, and L.-T. Wang, “Probing dark matter couplings to top and bottom quarks at the LHC,” Phys. Rev. D 88 no. 6, (2013) 063510, arXiv:1303.6638 [hep-ph].
  * [24] E. Izaguirre, G. Krnjaic, and B. Shuve, “The Galactic Center Excess from the Bottom Up,” Phys. Rev. D 90 no. 5, (2014) 055002, arXiv:1404.2018 [hep-ph].
  * [25] V. Cirigliano, B. Grinstein, G. Isidori, and M. B. Wise, “Minimal flavor violation in the lepton sector,” Nucl. Phys. B 728 (2005) 121–134, arXiv:hep-ph/0507001.
  * [26] Planck Collaboration, N. Aghanim et al., “Planck 2018 results. VI. Cosmological parameters,” Astron. Astrophys. 641 (2020) A6, arXiv:1807.06209 [astro-ph.CO]. [Erratum: Astron.Astrophys. 652, C4 (2021)].
  * [27] N. Craig, J. Galloway, and S. Thomas, “Searching for Signs of the Second Higgs Doublet,” arXiv:1305.2424 [hep-ph].
  * [28] G. Krnjaic, “Probing Light Thermal Dark-Matter With a Higgs Portal Mediator,” Phys. Rev. D 94 no. 7, (2016) 073009, arXiv:1512.04119 [hep-ph].
  * [29] P. Gondolo and G. Gelmini, “Cosmic abundances of stable particles: Improved analysis,” Nucl. Phys. B 360 (1991) 145–179.
  * [30] E. Izaguirre, G. Krnjaic, P. Schuster, and N. Toro, “Analyzing the Discovery Potential for Light Dark Matter,” Phys. Rev. Lett. 115 no. 25, (2015) 251301, arXiv:1505.00011 [hep-ph].
  * [31] A. Berlin, N. Blinov, G. Krnjaic, P. Schuster, and N. Toro, “Dark Matter, Millicharges, Axion and Scalar Particles, Gauge Bosons, and Other New Physics with LDMX,” Phys. Rev. D 99 no. 7, (2019) 075001, arXiv:1807.01730 [hep-ph].
  * [32] M. Fabbrichesi, E. Gabrielli, and G. Lanfranchi, “The Dark Photon,” arXiv:2005.01515 [hep-ph].
  * [33] J. Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, O. Mattelaer, H. S. Shao, T. Stelzer, P. Torrielli, and M. Zaro, “The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations,” JHEP 07 (2014) 079, arXiv:1405.0301 [hep-ph].
  * [34] BaBar Collaboration, J. P. Lees et al., “Search for a Dark Leptophilic Scalar in $e^{+}e^{-}$ Collisions,” Phys. Rev. Lett. 125 no. 18, (2020) 181801, arXiv:2005.01885 [hep-ex].
  * [35] ALEPH, DELPHI, L3, OPAL, SLD, LEP Electroweak Working Group, SLD Electroweak Group, SLD Heavy Flavour Group Collaboration, S. Schael et al., “Precision electroweak measurements on the $Z$ resonance,” Phys. Rept. 427 (2006) 257–454, arXiv:hep-ex/0509008.
  * [36] Particle Data Group Collaboration, P. A. Zyla et al., “Review of Particle Physics,” PTEP 2020 no. 8, (2020) 083C01.
  * [37] I. Agapov et al., “Future Circular Lepton Collider FCC-ee: Overview and Status,” in Snowmass 2021.  3, 2022.  arXiv:2203.08310 [physics.acc-ph].
  * [38] FCC Collaboration, A. Abada et al., “FCC Physics Opportunities: Future Circular Collider Conceptual Design Report Volume 1,” Eur. Phys. J. C 79 no. 6, (2019) 474.
  * [39] FCC Collaboration, A. Abada et al., “FCC-hh: The Hadron Collider: Future Circular Collider Conceptual Design Report Volume 3,” Eur. Phys. J. ST 228 no. 4, (2019) 755–1107.
  * [40] FCC Collaboration, A. Abada et al., “FCC-ee: The Lepton Collider: Future Circular Collider Conceptual Design Report Volume 2,” Eur. Phys. J. ST 228 no. 2, (2019) 261–623.
  * [41] J. F. Kamenik, A. Korajac, M. Szewc, M. Tammaro, and J. Zupan, “Flavor-violating Higgs and Z boson decays at a future circular lepton collider,” Phys. Rev. D 109 no. 1, (2024) L011301, arXiv:2306.17520 [hep-ph].
  * [42] T. S. M. Ho, X.-H. Jiang, T. H. Kwok, L. Li, and T. Liu, “Testing Lepton Flavor Universality at Future $Z$ Factories,” arXiv:2212.02433 [hep-ph].
  * [43] G. I. Budker, “Accelerators and colliding beams,” Conf. Proc. C 690827 (1969) 33–39.
  * [44] V. V. Parkhomchuk and A. N. Skrinsky, “Ionization cooling: Physics and applications,” AIP Conference Proceedings 352 no. 1, (01, 1996) 7–9, https://pubs.aip.org/aip/acp/article-pdf/352/1/7/11394715/7_1_online.pdf. https://doi.org/10.1063/1.49355.
  * [45] D. Neuffer, “Principles and applications of muon cooling,” Part.Accel. . https://www.osti.gov/biblio/1156195.
  * [46] MICE Collaboration, T. A. Mohayai, “First Demonstration of Ionization Cooling in MICE,” in 9th International Particle Accelerator Conference.  6, 2018.  arXiv:1806.01807 [physics.acc-ph].
  * [47] J.-P. Delahaye et al., “A Staged Muon Accelerator Facility For Neutrino and Collider Physics,” in 5th International Particle Accelerator Conference, p. WEZA02.  6, 2014.  arXiv:1502.01647.
  * [48] M. Boscolo, J.-P. Delahaye, and M. Palmer, “The future prospects of muon colliders and neutrino factories,” Rev. Accel. Sci. Tech. 10 no. 01, (2019) 189–214, arXiv:1808.01858 [physics.acc-ph].
  * [49] V. D. Barger, M. S. Berger, J. F. Gunion, and T. Han, “Precision W boson and top quark mass determinations at a muon collider,” Phys. Rev. D 56 (1997) 1714–1722, arXiv:hep-ph/9702334.
  * [50] T. Han, Y. Ma, and K. Xie, “High energy leptonic collisions and electroweak parton distribution functions,” Phys. Rev. D 103 no. 3, (2021) L031301, arXiv:2007.14300 [hep-ph].
  * [51] T. Han, W. Kilian, N. Kreher, Y. Ma, J. Reuter, T. Striegl, and K. Xie, “Precision test of the muon-Higgs coupling at a high-energy muon collider,” JHEP 12 (2021) 162, arXiv:2108.05362 [hep-ph].
  * [52] D. Buttazzo, D. Redigolo, F. Sala, and A. Tesi, “Fusing Vectors into Scalars at High Energy Lepton Colliders,” JHEP 11 (2018) 144, arXiv:1807.04743 [hep-ph].
  * [53] D. Buttazzo, R. Franceschini, and A. Wulzer, “Two Paths Towards Precision at a Very High Energy Lepton Collider,” JHEP 05 (2021) 219, arXiv:2012.11555 [hep-ph].
  * [54] M. Forslund and P. Meade, “Precision Higgs width and couplings with a high energy muon collider,” JHEP 01 (2024) 182, arXiv:2308.02633 [hep-ph].
  * [55] M. Ruhdorfer, E. Salvioni, and A. Wulzer, “Invisible Higgs boson decay from forward muons at a muon collider,” Phys. Rev. D 107 no. 9, (2023) 095038, arXiv:2303.14202 [hep-ph].
  * [56] C. Accettura et al., “Towards a muon collider,” Eur. Phys. J. C 83 no. 9, (2023) 864, arXiv:2303.08533 [physics.acc-ph]. [Erratum: Eur.Phys.J.C 84, 36 (2024)].
  * [57] R. Capdevilla, D. Curtin, Y. Kahn, and G. Krnjaic, “Systematically testing singlet models for (g - 2),” JHEP 04 (2022) 129, arXiv:2112.08377 [hep-ph].
  * [58] Muon g-2 Collaboration, G. W. Bennett et al., “Final Report of the Muon E821 Anomalous Magnetic Moment Measurement at BNL,” Phys. Rev. D 73 (2006) 072003, arXiv:hep-ex/0602035.
  * [59] T. Aoyama et al., “The anomalous magnetic moment of the muon in the Standard Model,” Phys. Rept. 887 (2020) 1–166, arXiv:2006.04822 [hep-ph].
  * [60] CMD-3 Collaboration, F. V. Ignatov et al., “Measurement of the $e^{+}e^{-}\to\pi^{+}\pi^{-}$ cross section from threshold to 1.2 GeV with the CMD-3 detector,” arXiv:2302.08834 [hep-ex].
  * [61] G. Krnjaic et al., “A Snowmass Whitepaper: Dark Matter Production at Intensity-Frontier Experiments,” arXiv:2207.00597 [hep-ph].
  * [62] Y. Kahn, G. Krnjaic, N. Tran, and A. Whitbeck, “M3: a new muon missing momentum experiment to probe (g -2) and dark matter at Fermilab,” JHEP 09 (2018) 153, arXiv:1804.03144 [hep-ph].
  * [63] Y. M. Andreev et al., “Exploration of the Muon g-2 and Light Dark Matter explanations in NA64 with the CERN SPS high energy muon beam,” arXiv:2401.01708 [hep-ex].
  * [64] M. A. Shifman, A. I. Vainshtein, and V. I. Zakharov, “Remarks on Higgs Boson Interactions with Nucleons,” Phys. Lett. B 78 (1978) 443–446.
  * [65] M. Cirelli, E. Del Nobile, and P. Panci, “Tools for model-independent bounds in direct dark matter searches,” JCAP 10 (2013) 019, arXiv:1307.5955 [hep-ph].
  * [66] XENON Collaboration, E. Aprile et al., “First Dark Matter Search with Nuclear Recoils from the XENONnT Experiment,” Phys. Rev. Lett. 131 no. 4, (2023) 041003, arXiv:2303.14729 [hep-ex].
  * [67] C. Boehm, M. J. Dolan, and C. McCabe, “A Lower Bound on the Mass of Cold Thermal Dark Matter from Planck,” JCAP 08 (2013) 041, arXiv:1303.6270 [hep-ph].
  * [68] K. M. Nollett and G. Steigman, “BBN And The CMB Constrain Light, Electromagnetically Coupled WIMPs,” Phys. Rev. D 89 no. 8, (2014) 083508, arXiv:1312.5725 [astro-ph.CO].
  * [69] G. Krnjaic and S. D. McDermott, “Implications of BBN Bounds for Cosmic Ray Upscattered Dark Matter,” Phys. Rev. D 101 no. 12, (2020) 123022, arXiv:1908.00007 [hep-ph].
  * [70] D. P. Finkbeiner, S. Galli, T. Lin, and T. R. Slatyer, “Searching for dark matter in the CMB: A compact parametrization of energy injection from new physics,” Physical Review D 85 no. 4, (Feb, 2012) . https://doi.org/10.1103%2Fphysrevd.85.043522.
  * [71] B. Batell, D. McKeen, and M. Pospelov, “Singlet Neighbors of the Higgs Boson,” JHEP 10 (2012) 104, arXiv:1207.6252 [hep-ph].
  * [72] CMS Collaboration, A. M. Sirunyan et al., “Measurement of the $\mathrm{Z}\gamma^{*}\to\tau\tau$ cross section in pp collisions at $\sqrt{s}=$ 13 TeV and validation of $\tau$ lepton analysis techniques,” Eur. Phys. J. C 78 no. 9, (2018) 708, arXiv:1801.03535 [hep-ex].