capbtabboxtable[][] L>l< C>c<

# Dark matter - electron interactions in materials beyond the dark photon
model

Riccardo Catena  Daniel Cole  Timon Emken  Marek Matas  Nicola Spaldin  Walter
Tarantino  and Einar Urdshals

###### Abstract

The search for sub-GeV dark matter (DM) particles via electronic transitions
in underground detectors attracted much theoretical and experimental interest
in the past few years. A still open question in this field is whether
experimental results can in general be interpreted in a framework where the
response of detector materials to an external DM probe is described by a
single ionisation or crystal form factor, as expected for the so-called dark
photon model. Here, ionisation and crystal form factors are examples of
material response functions: interaction-specific integrals of the initial and
final state electron wave functions. In this work, we address this question
through a systematic classification of the material response functions induced
by a wide range of models for spin-0, spin-1/2 and spin-1 DM. We find several
examples for which an accurate description of the electronic transition rate
at DM direct detection experiments requires material response functions that
go beyond those expected for the dark photon model. This concretely
illustrates the limitations of a framework that is entirely based on the
standard ionisation and crystal form factors, and points towards the need for
the general response-function-based formalism we pushed forward recently [1,
2]. For the models that require non-standard atomic and crystal response
functions, we use the response functions of [1, 2] to calculate the DM-induced
electronic transition rate in atomic and crystal detectors, and to present 90%
confidence level exclusion limits on the strength of the DM-electron
interaction from the null results reported by XENON10, XENON1T, EDELWEISS and
SENSEI.

##  1 Introduction

The nature of the dark matter (DM) component of the universe, whose
gravitational effects are visible from local to cosmological scales, remains
elusive. The leading hypothesis in astroparticle physics is that DM is made of
new, yet undetected particles [3]. Direct detection experiments play a crucial
role in the quest for the DM particle, and search for DM scattering events in
low-background detectors located deep underground [4, 5, 6]. So far, this
class of experiments has focused on the search for the Weakly Interacting
Massive Particle (WIMP): a DM candidate of mass between few GeV and tens of
TeV that interacts with the known particles with weak scale interactions, and
predicted by theories which were not initially formulated in a DM context,
e.g. supersymmetry [7, 8]. The lack of an unambiguous WIMP detection has
recently motivated a systematic, theoretical and experimental exploration of
alternative paradigms [9].

An alternative scenario that could explain why WIMPs have so far escaped
detection is the one where DM is lighter than the nucleons bound to nuclei,
and thus too light to induce an observable nuclear recoil in a direct
detection experiment [10]. An interesting aspect of this scenario is that it
is testable, as a Milky Way DM particle in the keV-GeV mass range would have
enough kinetic energy to cause an observable electronic transition in a
detector material, provided it interacts strongly enough with the constituent
electrons [11]. Being a testable explanation for the lack of WIMP discovery,
the sub-GeV DM paradigm has recently attracted much attention [12, 9]. In
particular, a variety of experimental tests has been proposed or performed to
probe this idea. This includes the search for: atomic ionisations in noble gas
xenon and argon detectors [11, 13, 14, 1, 15, 16], electronic transitions in
semiconductor crystals [17, 18, 19, 20, 21, 22, 23, 24, 25, 2, 26, 27, 28, 29,
30] as well as in superconductors [31, 32] and 3D Dirac materials [33, 34,
35], electron ejections from graphene [36] and carbon nanotubes [37], and
collective phenomena such as phonons [38, 39] and magnons [40]. This is an
incomplete list, and we refer to [12, 9] for an extended review of the field.

A central element in the development and implementation of the experimental
tests mentioned above is the theoretical modelling of DM-electron
interactions. Typically, the interactions between DM particles and electrons
are described within the “dark photon model”, which extends the Standard Model
of particle physics by an additional $U(1)$ gauge group under which only the
DM particle is charged [41, 10]. In this model, the DM-electron interaction
arises from the “kinetic mixing” between the ordinary photon and the gauge
boson associated with the new $U(1)$ group, i.e. the dark photon. A crucial
prediction of the model is that the response of materials to DM-electron
interactions depends on a single, material-specific ionisation or crystal form
factor [18]. While the dark photon model has extensively been used in the
interpretation of results from DM direct detection experiments, it is also a
rather restrictive framework, as it a priori neglects DM-electron couplings
that are allowed by observations, and assumes that the interactions between DM
and electrons are necessarily mediated by a spin-1 particle. Furthermore, the
model is often complemented by the additional assumption that the DM particle
has either spin 0 or 1/2. However, there is no experimental evidence
supporting these restrictions and, therefore, a more general approach is in
many respects preferable.

In order to develop a description of DM-electron scattering in detector
materials that goes beyond the one of the dark photon model, we recently
developed a non-relativistic effective theory for the scattering of Milky Way
DM particles by the electrons bound to isolated atoms [1] or semiconductor
crystals [2]. This allowed us to systematically classify the interactions
between DM particles and electrons in terms of a finite set of Galilean and
rotational invariant quantum mechanical operators that are defined in the DM-
electron spin space. Virtually any DM particle theory can be mapped onto a
linear combination of such operators with momentum transfer dependent
coefficients. Importantly, this general approach to the theoretical modelling
of DM-electron interactions predicts that up to four (for isolated atoms) and
up to five (for crystals) electron wave function overlap integrals, or
“response functions”, are required in order to accurately describe the
scattering of DM particles by the electrons bound in detector materials. In
other words, the response of materials to an external DM perturbation can in
principle be significantly different from the one predicted by the dark photon
model. Recently, our formalism was adopted by the XENON collaboration to
analyse their latest S2-only data from XENON1T [42].

A still open question in this field is whether the novel material response
functions we identified in [1] and [2] are indeed important in the
interpretation of DM direct detection experiments, or they only appear in
marginal models. While in [1, 2] we provided relevant examples, such as the
anapole or magnetic dipole DM models, where our “response function” formalism
is crucial to accurately predict observable electronic transition rates, a
systematic classification of the DM-electron interaction models that can
generate a non-standard material response to an external DM perturbation is
still missing. Answering this question is important, as it would clarify
whether present and future DM direct detection experiments can indeed be
interpreted within a framework that entirely relies on a single ionisation or
crystal form factor, or a general analysis of the response to DM-electron
interactions of detector materials is crucial.

In this work, we answer the above question by calculating the rate of
electronic transitions induced by the non-relativistic scattering of DM
particles by electrons bounds to isolated atoms or silicon and germanium
crystals for a wide range of models for spin-0, spin-1/2 and spin-1 DM.
Remarkably, we find eight relativistic Lagrangians which in the non-
relativistic limit predict a non-standard material response to external DM-
electron interactions. This corroborates the results we reported in [1, 2] and
shows the importance of a “response function” formalism in the interpretation
of DM direct detection experiments. This work is organised as follows. In Sec.
2, we introduce the relativistic models for DM-electron interactions we are
interested in. In Sec. 3 we calculate the rate of DM-induced electronic
transitions in atoms or crystals predicted by the models of Sec. 2 within the
response function formalism of [1, 2]. Here, we also show that the response
functions for xenon, germanium and silicon can be expressed within a unified
formalism that allows for a direct comparison of different target materials.
Finally, in Sec. 4, we compare our theoretical predictions with the null
results reported by the XENON10 [43], XENON1T[15], EDELWEISS [23] and SENSEI
[24] DM direct detection experiments, setting 90% confidence level upper
bounds on the coupling constants of the models we consider in this work. We
summarise and conclude in Sec. 5, and refer to three appendices for the
details of our calculations, useful identities and complementary results.

##  2 General models for dark matter-electron interactions

This section introduces the framework we use to describe DM-electron
interactions in detector materials.

In the interpretation of data from direct detection experiments, DM-electron
interactions are typically described within the so-called dark photon model
[41] (for a review, see also [10] and references therein). This framework
extends the Standard Model of particle physics by an additional $U(1)$ gauge
group and by a DM candidate. Only the DM particle is charged under the new
gauge group, while a “kinetic mixing” between ordinary and dark photon
generates the interactions between the electrically charged particles in the
Standard Model and the DM candidate, which, within this framework, typically
consists of either scalar or fermionic particles. Here, by “dark photon” one
refers to the gauge boson associated with the new $U(1)$ group.

In this work, we extend the standard dark photon model in three ways. First,
we consider both spin-0 and spin-1 particles as “mediators” of the DM-electron
interaction, without restricting ourselves to vector mediators, as in the dark
photon model. Second, we focus on DM-mediator interaction vertices with a
general Lorentz structure. In contrast, in the dark photon model only
interactions arising from a kinetic mixing are considered. Finally, we allow
for DM to have spin 0, spin 1/2 or spin 1. These extensions are motivated by
the lack of an experimentally preferred model for sub-GeV DM, which, in our
view, should be compensated by a general approach to the theoretical modelling
of DM-electron interactions. The main implication of these extensions is that
in the models we consider here the response of detector materials to DM-
electron interactions cannot be described via the canonical “form factors”
introduced within the dark photon model. As we will see below, it must be
modelled within the “response function” formalism we developed by using
effective theory methods in [1], for isolated atoms, and in [2], for crystals.
From this point of view, the analysis presented here provide a crucial,
concrete application of the formalism we established in [1, 2], which is
therefore proven to be relevant and important in the theoretical modelling of
DM-electron interactions. Furthermore, it also shows how the effective
operators identified in [1, 2] can arise from the non-relativistic reduction
of relativistic theories for DM-electron interactions. Below, we introduce the
models we are interested in, focusing on scalar, fermionic and vector DM
separately.

Notice that the models we consider here do not provide a mechanism for the DM
particle and mediator mass generation, and should therefore be completed in
the high-energy (i.e. ultra-violet) limit. While this completion is expected
to introduce new, phenomenologically interesting particles, it goes beyond the
purposes of this work. Therefore, we refer to [7] and references therein for a
review of possible ultra-violet completions of the models introduced below.

###  2.1 Scalar dark matter

Let us start by introducing our models for spin-0 DM. Here, we denote by $S$
the complex scalar that describes the DM particle, and consider two cases for
the particle that mediates the DM-electron interactions, i.e. “the mediator”.
In the first case, we assume a spin-0 particle of mass $m_{\phi}$ described by
the real scalar field $\phi$. In the second case we assume a spin-1 particle
of mass $m_{G}$ described by the real vector $G_{\mu}$. In the case of a
scalar mediator, $\phi$, we introduce the Lagrangian

$\displaystyle\mathcal{L}_{S\phi e}$ $\displaystyle=$
$\displaystyle\partial_{\mu}S^{\dagger}\partial^{\mu}S-m_{S}^{2}S^{\dagger}S-\frac{\lambda_{S}}{2}(S^{\dagger}S)^{2}$
(2.1) $\displaystyle+$
$\displaystyle\frac{1}{2}\partial_{\mu}\phi\partial^{\mu}\phi-\frac{1}{2}m_{\phi}^{2}\phi^{2}-\frac{m_{\phi}\mu_{1}}{3}\phi^{3}-\frac{\mu_{2}}{4}\phi^{4}$
$\displaystyle+$ $\displaystyle i\overline{e}\not{D}e-m_{e}\overline{e}e$
$\displaystyle-$ $\displaystyle
g_{1}m_{S}S^{\dagger}S\phi-\frac{g_{2}}{2}S^{\dagger}S\phi^{2}-h_{1}\overline{e}e\phi-
ih_{2}\overline{e}\gamma^{5}e\phi\,,$

where $m_{S}$ is the DM mass, while $\lambda_{S}$, $\mu_{1}$, $\mu_{2}$,
$g_{1}$, $g_{2}$, $h_{1}$, and $h_{2}$ are coupling constants. In Eq. (2.1),
$m_{e}$ and $e$ are the electron mass and Dirac spinor, respectively. On the
other hand, when the interaction between DM particles and electrons is
mediated by a spin-1 particle, we assume

$\displaystyle\mathcal{L}_{SGe}$ $\displaystyle=$
$\displaystyle\partial_{\mu}S^{\dagger}\partial^{\mu}S-m_{S}^{2}S^{\dagger}S-\frac{\lambda_{S}}{2}(S^{\dagger}{S})^{2}$
(2.2) $\displaystyle-$
$\displaystyle\frac{1}{4}\mathcal{G}_{\mu\nu}\mathcal{G}^{\mu\nu}+\frac{1}{2}m_{G}^{2}G_{\mu}G^{\mu}-\frac{\lambda_{G}}{4}(G_{\mu}G^{\mu})^{2}$
$\displaystyle+$ $\displaystyle i\overline{e}\not{D}e-m_{e}\overline{e}e$
$\displaystyle-$
$\displaystyle\frac{g_{3}}{2}S^{\dagger}SG_{\mu}G^{\mu}-ig_{4}(S^{\dagger}\partial_{\mu}S-\partial_{\mu}S^{\dagger}S)G^{\mu}$
$\displaystyle-$ $\displaystyle
h_{3}(\overline{e}\gamma_{\mu}e)G^{\mu}-h_{4}(\overline{e}\gamma_{\mu}\gamma^{5}e)G^{\mu}\,,$

where $\lambda_{G}$, $g_{1}$, $g_{2}$, $h_{1}$, and $h_{2}$ are coupling
constants, while $\mathcal{G}_{\mu\nu}$ is the field strength tensor for
$G_{\mu}$.

###  2.2 Fermionic dark matter

We now focus on the case of spin-1/2 DM with mass $m_{\chi}$. When the
interactions between the DM particles and electrons are mediated by a scalar
particle of mass $m_{\phi}$, we assume the Lagrangian,

$\displaystyle\mathcal{L}_{\chi\phi e}$ $\displaystyle=$ $\displaystyle
i\overline{\chi}\not{D}\chi-m_{\chi}\overline{\chi}\chi$ (2.3)
$\displaystyle+$
$\displaystyle\frac{1}{2}\partial_{\mu}\phi\partial^{\mu}\phi-\frac{1}{2}m_{\phi}^{2}\phi^{2}-\frac{m_{\phi}\mu_{1}}{3}\phi^{3}-\frac{\mu_{2}}{4}\phi^{4}$
$\displaystyle+$ $\displaystyle i\overline{e}\not{D}e-m_{e}\overline{e}e$
$\displaystyle-$
$\displaystyle\lambda_{1}\phi\overline{\chi}\chi-i\lambda_{2}\phi\overline{\chi}\gamma^{5}\chi-
h_{1}\phi\overline{e}e-ih_{2}\phi\overline{e}\gamma^{5}e\,,$

where $\lambda_{1}$ and $\lambda_{2}$ are coupling constants, while $\chi$ is
the Dirac spinor that describes the DM particle. Alternatively, in the case of
a vector mediator, $G_{\mu}$, we assume the following Lagrangian

$\displaystyle\mathcal{L}_{\chi Ge}$ $\displaystyle=$ $\displaystyle
i\overline{\chi}\not{D}\chi-m_{\chi}\overline{\chi}\chi$ (2.4)
$\displaystyle-$
$\displaystyle\frac{1}{4}\mathcal{G}_{\mu\nu}\mathcal{G}^{\mu\nu}+\frac{1}{2}m_{G}^{2}G_{\mu}G^{\mu}-\frac{\lambda_{G}}{4}(G_{\mu}G^{\mu})^{2}$
$\displaystyle+$ $\displaystyle i\overline{e}\not{D}e-m_{e}\overline{e}e$
$\displaystyle-$ $\displaystyle\lambda_{3}\overline{\chi}\gamma^{\mu}\chi
G_{\mu}-\lambda_{4}\overline{\chi}\gamma^{\mu}\gamma^{5}\chi G_{\mu}$
$\displaystyle-$ $\displaystyle
h_{3}\overline{e}\gamma_{\mu}eG^{\mu}-h_{4}\overline{e}\gamma_{\mu}\gamma^{5}eG^{\mu}\,,$

where $h_{3}$ and $h_{4}$ are additional coupling constants, and the remaining
symbols were introduced above.

###  2.3 Vector dark matter

Finally, we focus on the case of spin-1 DM particles with mass $m_{X}$
described by a complex vector field, $X_{\mu}$. For DM-electron interactions
mediated by a scalar field, $\phi$, we assume the Lagrangian,

$\displaystyle\mathcal{L}_{X\phi e}$ $\displaystyle=$
$\displaystyle-\frac{1}{2}{\mathcal{X}}_{\mu\nu}^{\dagger}\mathcal{X}^{\mu\nu}+m_{X}^{2}X_{\mu}^{\dagger}X^{\mu}-\frac{\lambda_{X}}{2}(X_{\mu}^{\dagger}X^{\mu})^{2}$
(2.5) $\displaystyle+$
$\displaystyle\frac{1}{2}(\partial_{\mu}\phi)^{2}-\frac{1}{2}m_{\phi}^{2}\phi^{2}-\frac{m_{\phi}\mu_{1}}{3}\phi^{3}-\frac{\mu_{2}}{4}\phi^{4}$
$\displaystyle+$ $\displaystyle i\overline{e}\not{D}e-m_{e}\overline{e}e$
$\displaystyle-$ $\displaystyle b_{1}m_{X}\phi
X_{\mu}^{\dagger}X^{\mu}-\frac{b_{2}}{2}\phi^{2}X_{\mu}^{\dagger}X^{\mu}$
$\displaystyle-$ $\displaystyle
h_{1}\phi\overline{e}e-ih_{2}\phi\overline{e}\gamma^{5}e\,,$

where $\lambda_{X}$, $b_{1}$ and $b_{2}$ are new coupling constants, and
$\mathcal{X}_{\mu\nu}$ is the field strength tensor for $X_{\mu}$. On the
other hand, in the case of DM-electron interactions mediated by a spin-1
particle, we assume

$\displaystyle\mathcal{L}_{XGe}$ $\displaystyle=$
$\displaystyle-\frac{1}{2}\mathcal{X}^{\dagger}_{\mu\nu}\mathcal{X}^{\mu\nu}+m_{X}^{2}X^{\dagger}_{\mu}X^{\mu}-\frac{\lambda_{X}}{2}(X_{\mu}^{\dagger}X^{\mu})^{2}$
(2.6) $\displaystyle-$
$\displaystyle\frac{1}{4}\mathcal{G}_{\mu\nu}\mathcal{G}^{\mu\nu}+\frac{1}{2}m_{G}^{2}G_{\mu}G^{\mu}-\frac{\lambda_{G}}{4}(G_{\mu}G^{\mu})^{2}$
$\displaystyle+$ $\displaystyle i\overline{e}\not{D}e-m_{e}\overline{e}e$
$\displaystyle-$
$\displaystyle\frac{b_{3}}{2}G_{\mu}G^{\mu}(X^{\dagger}_{\nu}X^{\nu})-\frac{b_{4}}{2}(G^{\mu}G^{\nu})(X^{\dagger}_{\mu}X_{\nu})$
$\displaystyle-$
$\displaystyle\left[ib_{5}X_{\nu}^{\dagger}\partial_{\mu}X^{\nu}G^{\mu}+b_{6}X_{\mu}^{\dagger}\partial^{\mu}X_{\nu}G^{\nu}\right.$
$\displaystyle+$
$\displaystyle\left.b_{7}\varepsilon_{\mu\nu\rho\sigma}(X^{\dagger\mu}\partial^{\nu}X^{\rho})G^{\sigma}+h.c.\right]$
$\displaystyle-$ $\displaystyle
h_{3}G_{\mu}\overline{e}\gamma^{\mu}e-h_{4}G_{\mu}\overline{e}\gamma^{\mu}\gamma^{5}e\,,$

where $b_{3}$, $b_{4}$, $b_{5}$, $b_{6}$ and $b_{7}$ are coupling constants.
Contrary to all other coupling constants, which can without loss of generality
be assumed to be real, $b_{6}$ and $b_{7}$ are in general complex variables
[44].

##  3 Expected rate of electron transitions in detector materials

From the Lagrangians introduced in Sec. 2, we now calculate the expected rate
of electron transitions induced by DM scattering in isolated atoms and
crystals. We start by reviewing the results of [1, 2] (Sec. 3.1). We then
apply these results to the general models of Sec. 2 (Sec. 3.2).

###  3.1 General considerations

Let us consider the scattering of a DM particle of mass $m_{\rm DM}$, initial
velocity in the detector rest frame $\bm{v}$, and momentum $\bm{p}=m_{\rm
DM}\bm{v}$ by an electron bound to an isolated atom or to a crystal, and
described by the initial state $\ket{\mathbf{e}_{1}}$. Furthermore, let us
denote by $\ket{\mathbf{e}_{2}}$ the state of the outgoing electron, and by
$\bm{p}^{\prime}$ the momentum of the final state DM particle. The momentum
transfer in the scattering is thus $\bm{q}=\bm{p}-\bm{p}^{\prime}$, and the
rate for the $\ket{\mathbf{e}_{1}}$ to $\ket{\mathbf{e}_{2}}$ transition is
given by [1, 2]

$\displaystyle\mathscr{R}_{1\rightarrow 2}$ $\displaystyle=\frac{n_{\rm
DM}}{16m^{2}_{\rm DM}m^{2}_{e}}\,\int\frac{{\rm d}^{3}q}{(2\pi)^{3}}\int{\rm
d}^{3}\bm{v}\,f(\bm{v})(2\pi)\delta(E_{f}-E_{i})\overline{\left|\mathcal{M}_{1\rightarrow
2}\right|^{2}}\,,$ (3.1)

where $E_{i}$ and $E_{f}$ are the total initial and final energy,
respectively. They can be written as follows

$\displaystyle E_{i}$ $\displaystyle=m_{\rm
DM}+m_{e}+\frac{|\bm{p}|^{2}}{2m_{\rm DM}}+E_{1}$ $\displaystyle E_{f}$
$\displaystyle=m_{\rm DM}+m_{e}+\frac{|\bm{p}-\bm{q}|^{2}}{2m_{\rm
DM}}+E_{2}\,.$ (3.2)

Here, $E_{1}<0$ ($E_{2}>0$) is the energy eigenvalue of the state
$\ket{\mathbf{e}_{1}}$ ($\ket{\mathbf{e}_{2}}$). This convention for the signs
of $E_{1}$ and $E_{2}$ reflects our interest in transitions from bound atomic
(valence) to free (conduction) states (see Secs. 3.1.1 and 3.1.2). In Eq.
(3.1), $\mathcal{M}_{1\rightarrow 2}$ is the electron transition amplitude,
defined as the overlap integral,

$\displaystyle\mathcal{M}_{1\rightarrow 2}(\bm{q},\bm{v})$
$\displaystyle=\int\frac{{\rm
d}^{3}\bm{\ell}}{(2\pi)^{3}}\,\widetilde{\psi}_{2}^{*}(\bm{\ell}+\bm{q})\,\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el}(\bm{v},\bm{q},\bm{\ell}))\,\widetilde{\psi}_{1}(\bm{\ell})\,,$ (3.3)

where $\widetilde{\psi}_{1}$ ($\widetilde{\psi}_{2}$) is the Fourier transform
of $\psi_{1}$ ($\psi_{2}$), while $\psi_{1}$ ($\psi_{2}$) is the wave function
associated with the state $\ket{\mathbf{e}_{1}}$ ($\ket{\mathbf{e}_{2}}$), and
$\mathcal{M}$ is the free electron scattering amplitude. In the non-
relativistic, elastic scattering of a DM particle by free electrons,
$\mathcal{M}$ uniquely depends on the momentum transfer $\bm{q}$ and on the
transverse relative velocity $\bm{v}^{\perp}_{\rm el}=\bm{v}-\bm{q}/(2\mu_{\rm
DMe})-\bm{\ell}/m_{e}$, where $\mu_{\rm DMe}$ is the DM-electron reduced mass
and $\bm{\ell}$ is the initial state electron momentum. This simplification is
a consequence of Galilean invariance (that is, the invariance under constant
shifts of particle velocities) and three-dimensional momentum conservation,
which allow one to express $\mathcal{M}$ as a function of two kinematic
variables only (see [45, 46, 47, 1, 2] for further details). Notice that
Galilean invariance is fulfilled when both the DM particle and the electron
are non-relativistic. We can now expand $\mathcal{M}$ at linear order in
$\bm{\ell}/m_{e}$, which implies [1, 2],

$\displaystyle\mathcal{M}_{1\rightarrow 2}(\bm{q},\bm{v})$
$\displaystyle\simeq\left.\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})\right|_{\bm{\ell}=\mathbf{0}}f_{1\rightarrow
2}(\bm{q})+m_{e}\left.\nabla_{\bm{\ell}}\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})\right|_{\bm{\ell}=\mathbf{0}}\cdot\bm{f}_{1\rightarrow 2}(\bm{q})\,,$
(3.4)

where in the right-hand-side we introduced the scalar and vectorial wave
function overlap integrals,

$\displaystyle f_{1\rightarrow 2}(\bm{q})$
$\displaystyle\equiv\int\mathrm{d}^{3}\textbf{x}\,\psi_{2}^{*}(\mathbf{x})\,e^{i\bm{q}\cdot\mathbf{x}}\,\psi_{1}\left(\mathbf{x}\right)\,,$
(3.5) $\displaystyle\bm{f}_{1\rightarrow 2}(\bm{q})$
$\displaystyle\equiv\int\mathrm{d}^{3}\textbf{x}\,\psi_{2}^{*}(\mathbf{x})\,e^{i\bm{q}\cdot\mathbf{x}}\,\frac{i\nabla_{\mathbf{x}}}{m_{e}}\psi_{1}\left(\mathbf{x}\right)\,.$
(3.6)

In Eq. (3.1)$,n_{\rm DM}=\rho_{\rm DM}/m_{\rm DM}$ is the local DM number
density, $\rho_{\rm DM}=0.4\,\text{GeV cm}^{-3}$ the local DM mass density
[48], and $f(\bm{v})$ the local DM velocity distribution boosted to the
detector rest frame. For $f(\bm{v})$, we assume a truncated Maxwell-Boltzmann
distribution with local standard of rest speed $v_{0}=~{}238$ km s-1 [49],
galactic escape speed $v_{\mathrm{esc}}=544$ km s-1 [49] and Earth’s speed in
the galactic frame (where the mean DM particle velocity is zero) $v_{e}=250.5$
km s-1 [49].

The theoretical rate of DM-induced electron transitions in a detector
material, $\mathscr{R}_{\rm theory}$, is related to the rate for state
$\ket{\mathbf{e}_{1}}$ to state $\ket{\mathbf{e}_{2}}$ transitions,
$\mathscr{R}_{1\rightarrow 2}$, by a sum over all “relevant” initial and final
states. For detectors that can be approximated by a system of isolated atoms,
and that search for atomic ionisations induced by the scattering of DM
particles in the detector, the relevant initial states are bound states
characterised by the principal, angular and magnetic quantum numbers, $n$,
$\ell$ and $m$, respectively. The relevant final states are free electron
states characterised by an “asymptotic” momentum $k^{\prime}$, and by the
electron angular and magnetic quantum numbers, $\ell^{\prime}$ and
$m^{\prime}$, respectively. We refer to the appendices of [1] for an explicit
form of the corresponding wave functions. By replacing $1$ with the quantum
numbers $\\{n\ell m\\}$, and $2$ with
$\\{k^{\prime}\ell^{\prime}m^{\prime}\\}$ in Eq. (3.4), and inserting Eq.
(3.4) into Eq. (3.1), we obtain

$\displaystyle\mathscr{R}_{\rm theory}=\mathscr{R}_{\rm ion}$
$\displaystyle\equiv
2\sum_{n,\ell}\sum_{m=-\ell}^{\ell}\sum_{\ell^{\prime}=0}^{\infty}\sum_{m^{\prime}=-\ell^{\prime}}^{\ell^{\prime}}\int\mathrm{d}k^{\prime}\frac{Vk^{\prime
2}}{(2\pi)^{3}}\mathscr{R}_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}}\,,$ (3.7)

where the $(\ell,m)$ sum runs over the outermost occupied orbitals, e.g. the
4s, 4p, 4d, 5s, and 5p orbitals in the case of xenon atoms, which is of
interest for the present work. In Eq. (3.7) we assume that wave functions have
been normalised over a finite volume $V$, and the overall factor of 2 accounts
for the double occupation of each electronic state due to spin degeneracy. For
crystal detectors that search for electronic transitions from valence to
conduction bands in semiconductor targets, the relevant initial states are
occupied valence states labeled by a band index $i$ and a reciprocal lattice
momentum in the first Brillouin Zone (BZ), $\bm{k}$. Similarly, the relevant
final states are in this case unoccupied conduction states characterised by a
band index $i^{\prime}$ and a reciprocal lattice momentum in the first BZ,
$\bm{k}^{\prime}$. For crystal detectors, we thus find

$\displaystyle\mathscr{R}_{\rm theory}=\mathscr{R}_{\rm crystal}\equiv
2\sum_{ii^{\prime}}\int_{\rm BZ}\frac{V{\rm d}^{3}k}{(2\pi)^{3}}\int_{\rm
BZ}\frac{V{\rm d}^{3}k^{\prime}}{(2\pi)^{3}}\,\mathscr{R}_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}\,.$ (3.8)

We refer to [2] for an explicit form of the Bloch wave functions
$\psi_{i\bm{k}}$ and $\psi_{i^{\prime}\bm{k}^{\prime}}$ entering Eqs. (3.5)
and (3.6).

In general, we express the rate $\mathscr{R}_{\rm theory}$ within our
“response function” formalism [1, 2], where

$\displaystyle\mathscr{R}_{\rm theory}=\frac{n_{\rm DM}}{128\pi m_{\rm
DM}^{2}m_{e}^{2}}$ $\displaystyle\int\mathrm{d}(\ln\Delta
E)\int\mathrm{d}q\,q\,\widehat{\eta}\left(q,\Delta
E\right)\sum_{l=1}^{r}\mathrm{Re}\left[\mathcal{R}_{l}^{*}(q,v)\overline{\mathcal{W}}_{l}(q,\Delta
E)\right]\,.$ (3.9)

With an appropriate choice of “material response” functions
$\overline{\mathcal{W}}_{l}(q,\Delta E)$, and “DM response” functions
$\mathcal{R}_{l}(q,v)$, Eq. (3.9) applies to any detector material and
particle physics model. In particular, it applies to detectors that can be
approximated by a set of isolated atoms, such as XENON10 and XENON1T, and to
crystal detectors, such as SENSEI and EDELWEISS. The velocity operator in Eq.
(3.9), $\widehat{\eta}\left(q,\Delta E\right)$, depends on the momentum
transfer, $\bm{q}$, and on the deposited energy, $\Delta E$. It linearly acts
on the generic function of velocities, $g(\bm{v})$, as shown below

$\displaystyle\widehat{\eta}(q,\Delta E)\left[g(\bm{v})\right]$
$\displaystyle=\int_{|\bm{v}|\geq v_{\rm min}}{\rm
d}^{3}v\,g(\bm{v})\frac{f_{\chi}(\bm{v})}{v}\,,$ (3.10)

where

$\displaystyle v_{\rm min}=\frac{q}{2m_{\rm DM}}+\frac{\Delta E}{q}\,,$ (3.11)

is the minimum velocity required to transfer a momentum $\bm{q}$ when the
deposited energy is $\Delta E$, and $q\equiv|\bm{q}|$. Below, we first provide
explicit expressions for the $\overline{\mathcal{W}}_{l}(q,\Delta E)$,
$l=1,\dots r$ functions that characterise detectors consisting of isolated
atoms (Sec. 3.1.1) and crystals (Sec. 3.1.2). We then calculate the functions
$\mathcal{R}_{l}(q,v)$, $l=1,\dots r$, for all models in Sec. 2 (Sec. 3.2).

####  3.1.1 Isolated atoms

The material response functions $\overline{\mathcal{W}}_{l}$ are defined via
an integral over the azimuthal and polar angles of $\bm{q}$,

$\displaystyle\overline{\mathcal{W}}_{l}(q,\Delta E)\equiv$
$\displaystyle\int\mathrm{d}\Omega_{q}\,\mathcal{W}_{l}(\bm{q},\Delta E)\,,$
(3.12)

where, for detectors that can be modelled as a sample of isolated atoms, such
as XENON10 or XENON1T,

$\displaystyle\mathcal{W}_{l}(\bm{q},\Delta E)$
$\displaystyle=\frac{2}{\pi}\Delta
E\sum_{n,\ell,m}\sum_{\ell^{\prime},m^{\prime}}\int\frac{Vk^{\prime
2}\,\mathrm{d}k^{\prime}}{(2\pi)^{3}}\,\mathscr{B}_{l}\,\delta(\Delta
E-E_{k^{\prime}\ell^{\prime}m^{\prime}}+E_{n\ell m})\,,$ (3.13)

and

$\displaystyle\mathscr{B}_{1}=$ $\displaystyle\left|f_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}}\right|^{2}$ (3.14)
$\displaystyle\mathscr{B}_{2}=$
$\displaystyle\frac{\bm{q}}{m_{e}}\cdot(f_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}})(\bm{f}_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}})^{*}$ (3.15)
$\displaystyle\mathscr{B}_{3}=$ $\displaystyle\left|\bm{f}_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}}\right|^{2}$ (3.16)
$\displaystyle\mathscr{B}_{4}=$
$\displaystyle\left|\frac{\bm{q}}{m_{e}}\cdot\bm{f}_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}}\right|^{2}\,.$ (3.17)

$E_{n\ell m}$ ($E_{k^{\prime}\ell^{\prime}m^{\prime}}$) is the energy of the
initial (final) state electron and $f_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}}$ ($\bm{f}_{n\ell m\rightarrow
k^{\prime}\ell^{\prime}m^{\prime}}$) replaces $f_{1\rightarrow 2}$
($\bm{f}_{1\rightarrow 2}$). To compare the results reported here with those
in [1], we note that the $W_{l}^{n\ell}(k^{\prime},q)$ functions introduced in
[1] and the $\overline{\mathcal{W}}_{l}(q,\Delta E)$ functions defined here
are related by

$\displaystyle\overline{\mathcal{W}}_{l}(q,\Delta E)=2\Delta
E\sum_{n,\ell}\int\frac{{\rm
d}k^{\prime}}{k^{\prime}}\,W_{l}^{n\ell}(k^{\prime},q)\,\delta(\Delta
E-E_{k^{\prime}\ell^{\prime}m^{\prime}}+E_{n\ell m})\,.$ (3.18)

Since we are interested in the ionisation of isolated atoms, in Eq. (3.18)
$E_{k^{\prime}\ell^{\prime}m^{\prime}}=k^{\prime 2}/(2m_{e})$ and the initial
state energy, $E_{n\ell m}$, does not depend on $m$. Consequently, Eq. (3.18)
consistently gives a response function $\overline{\mathcal{W}}_{l}(q,\Delta
E)$ that is independent of quantum numbers. Also, requiring that $v_{\rm
min}(q,\Delta E)$ is less than the maximum speed a DM particle can have in the
galaxy, $v_{\rm max}=v_{e}+v_{\rm esc}$, implies that only certain regions in
the $(q,\Delta E)$ plane are kinematically allowed.

We conclude this subsection by providing an explicit relation between our
response function $\overline{\mathcal{W}}_{1}$ and the standard ionisation
form factor $f^{n\ell}_{\rm ion}$ used within the dark photon model, namely

$\displaystyle\overline{\mathcal{W}}_{1}(q,\Delta E)=2\Delta
E\sum_{n,\ell}\int\frac{{\rm d}k^{\prime}}{k^{\prime}}\,|f_{\rm
ion}^{n\ell}(k^{\prime},q)|^{2}\,\delta(\Delta
E-E_{k^{\prime}\ell^{\prime}m^{\prime}}+E_{n\ell m})\,.$ (3.19)

In the numerical applications, we compute the response functions
$W_{l}^{n\ell}(k^{\prime},q)$ by using the DarkART code [50], which implements
the equations we derived in [1].

####  3.1.2 Crystals

Similarly, in the case of crystal detectors such as SENSEI and EDELWEISS one
can use Eq. (3.12) with

$\displaystyle\mathcal{W}_{l}(\bm{q},\Delta E)$
$\displaystyle=\frac{2}{\pi}\Delta
E\sum_{ii^{\prime}}\int_{\textrm{BZ}}\frac{V\mathrm{d}^{3}k}{(2\pi)^{3}}\int_{\textrm{BZ}}\frac{V\mathrm{d}^{3}k^{\prime}}{(2\pi)^{3}}\,\mathscr{B}_{l}\,\delta(\Delta
E-E_{i^{\prime}\bm{k}^{\prime}}+E_{i\bm{k}})\,$ (3.20)

to calculate the response functions $\overline{\mathcal{W}}_{l}(q,\Delta E)$.
In this case, $r=5$ functions $\mathscr{B}_{l}$ are non zero [2],

$\displaystyle\mathscr{B}_{1}=$ $\displaystyle\left|f_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}\right|^{2}$ (3.21) $\displaystyle\mathscr{B}_{2}=$
$\displaystyle\frac{\bm{q}}{m_{e}}\cdot(f_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}^{\prime})(\bm{f}_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}})^{*}$ (3.22) $\displaystyle\mathscr{B}_{3}=$
$\displaystyle\left|\bm{f}_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}\right|^{2}$ (3.23) $\displaystyle\mathscr{B}_{4}=$
$\displaystyle\left|\frac{\bm{q}}{m_{e}}\cdot\bm{f}_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}\right|^{2}$ (3.24) $\displaystyle\mathscr{B}_{5}=$
$\displaystyle i\frac{\bm{q}}{m_{e}}\cdot\left[\bm{f}_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}\times\left(\bm{f}_{i\bm{k}\rightarrow
i^{\prime}\bm{k}^{\prime}}\right)^{*}\right]\,.$ (3.25)

Interestingly, the response function induced by $\mathscr{B}_{5}$ is
identically zero in the case of isolated atoms [1]. Here, $E_{i\bm{k}}$
($E_{i^{\prime}\bm{k}^{\prime}}$) is the energy of the $i\bm{k}$
($i^{\prime}\bm{k}^{\prime}$) initial state valence (final state conduction)
electron, while $f_{1\rightarrow 2}$ ($\bm{f}_{1\rightarrow 2}$) is replaced
by $f_{i\bm{k}}$ ($\bm{f}_{i^{\prime}\bm{k}^{\prime}}$), and $V=N_{\rm
cell}V_{\rm cell}$, where $V_{\rm cell}$ is the volume of a unit cell, and
$N_{\rm cell}$ the number of unit cells in the crystal. Note that the
functions $\overline{\mathcal{W}}_{l}$ defined via Eq. (3.20) and the
functions $\overline{W}_{l}$ introduced in [2] are related as follows

$\displaystyle\overline{\mathcal{W}}_{l}(q,\Delta E)=N_{\rm
cell}\,\overline{W}_{l}(q,\Delta E)\,.$ (3.26)

By absorbing $N_{\rm cell}$ in the definition of $\overline{\mathcal{W}}_{l}$,
Eq. (3.9) applies to both isolated atoms and crystals. Finally,
$\overline{\mathcal{W}}_{1}$ and the crystal form factor used in the dark
photon model, $f_{\text{crystal}}$, are related by

$\overline{\mathcal{W}}_{1}(q,\Delta E)=\frac{8\Delta E\alpha m_{e}^{2}N_{\rm
cell}}{q^{3}}\left|f_{\text{crystal}}\left(q,\Delta E\right)\right|^{2}\,,$
(3.27)

where $\alpha$ is the fine structure constant. In the numerical applications,
we compute the response functions $\overline{W}_{l}(q,\Delta E)$ by using the
QEdark-EFT code [51], which implements the equations we derived in [2].

####  3.1.3 Unified notation

Remarkably, the response functions $\mathcal{W}_{l}(\bm{q},\Delta E)$ for
isolated atoms, Eq. (3.13), and those for crystals, Eq. (3.20), have the very
same form. Specifically, they can both be expressed as follows

$\displaystyle\mathcal{W}_{l}(\bm{q},\Delta E)$
$\displaystyle=\frac{2}{\pi}\Delta
E\sum_{\\{1\\},\\{2\\}}\mathscr{B}_{l}\,\delta(\Delta E-E_{2}+E_{1})\,,$
(3.28)

where $\sum_{\\{1\\},\\{2\\}}$ represents the sums/integrals over all relevant
physical states in Eq. (3.13) and Eq. (3.20). In the case of isolated atoms,
Eq. (3.28) gives the “detector response” for a single atom, whereas in the
case of crystal detectors Eq. (3.28) gives the “detector response” for a
single crystal. In order to compare the two quantities, we introduce a
detector response per unit mass,

$\displaystyle\widetilde{\mathcal{W}}_{l}(q,\Delta
E)\equiv\overline{\mathcal{W}}_{l}(q,\Delta E)/\widetilde{m}\,,$ (3.29)

where for isolated atoms $\widetilde{m}=m_{\rm Xe}$, $m_{\rm Xe}$ being the
mass of a xenon atom, whereas for crystals, $\widetilde{m}=m_{\rm cell}N_{\rm
cell}$, $m_{\rm cell}$ and $N_{\rm cell}$ being the mass and number of unit
cells. Fig. 1 shows the ratios

$\displaystyle(\widetilde{\mathcal{W}}_{l}^{(\rm
Xe)}-\widetilde{\mathcal{W}}_{l}^{(\rm
Ge)})/(\widetilde{\mathcal{W}}_{l}^{(\rm
Xe)}+\widetilde{\mathcal{W}}_{l}^{(\rm Ge)})$ (3.30)

and

$\displaystyle(\widetilde{\mathcal{W}}_{l}^{(\rm
Xe)}-\widetilde{\mathcal{W}}_{l}^{(\rm
Si)})/(\widetilde{\mathcal{W}}_{l}^{(\rm
Xe)}+\widetilde{\mathcal{W}}_{l}^{(\rm Si)})$ (3.31)

in the $(q,\Delta E)$ plane for $l=1$ and $l=2$. Similarly, Fig. 2 shows the
same ratios now for $l=3$ and $l=4$. In both cases, a superscript specifies
the material for which $\widetilde{\mathcal{W}}_{l}(q,\Delta E)$ is evaluated.
Notice that only the regions above the black dashed lines are kinematically
allowed in the two figures for our choice of $v_{\rm esc}$ and $v_{e}$. For
larger $v_{\rm esc}$, for example, smaller $|\bm{q}|$ would be accessible. As
expected, for $\Delta E$ below the xenon ionisation threshold,
$\widetilde{\mathcal{W}}_{l}^{(\rm Ge)}$ and
$\widetilde{\mathcal{W}}_{l}^{(\rm Si)}$ dominate over
$\widetilde{\mathcal{W}}_{l}^{(\rm Xe)}$, showing that xenon and crystal
detectors are mainly sensitive to complementary regions in the $(q,\Delta E)$
plane. Despite this general trend, for deposited energies above about 30 eV
the germanium $3d$ bands contribute significantly to
$\widetilde{\mathcal{W}}_{l}^{(\rm Ge)}$ and the ratio
$\widetilde{\mathcal{W}}_{l}^{(\rm Ge)}/\widetilde{\mathcal{W}}_{l}^{(\rm
Xe)}$ can be of the order of 1 for kinematically allowed values of the
momentum transfer. While for $l=1$, $l=3$ and $l=4$ the response function
$\overline{\mathcal{W}}_{l}(q,\Delta E)$ is always real, for $l=2$, we replace
$\overline{\mathcal{W}}_{2}(q,\Delta E)$ with
$|\overline{\mathcal{W}}_{2}(q,\Delta E)|$ in Eqs. (3.30) and (3.31), as the
former takes complex values for crystals. In the case of isolated atoms, all
material response functions are real.

Figure 1: Relative differences between detector response functions per unit
mass for xenon and germanium (left), and for xenon and silicon (right)
detectors in the $(q,\Delta E)$ plane. The top (bottom) panels refer to $l=1$
($l=2$). Regions above the black dashed line are kinematically allowed.

Figure 2: The same as for Fig. 1, now for $l=3$ and $l=4$.

With the notation of Eq. (3.28), we can also write a general relation between
products of material and DM response functions, and the modulus squared of the
transition amplitude, namely

$\displaystyle\sum_{l=1}^{r}\mathrm{Re}\left[\mathcal{R}_{l}^{*}(q,v)\overline{\mathcal{W}}_{l}(q,\Delta
E)\right]$ $\displaystyle=\frac{2}{\pi}\Delta E\sum_{\\{1\\},\\{2\\}}\int{\rm
d}\Omega_{q}\,\overline{\overline{|\mathcal{M}_{1\rightarrow
2}(\bm{v},\bm{q})|^{2}}}\,\delta(\Delta E-E_{2}+E_{1})\,,$ (3.32)

where we used [2]

$\displaystyle\overline{\overline{|\mathcal{M}_{1\rightarrow
2}(\bm{v},\bm{q})|^{2}}}$ $\displaystyle\equiv\frac{1}{2\pi}\int{\rm
d}\phi\,{\overline{|\mathcal{M}_{1\rightarrow
2}(\bm{v},\bm{q})|^{2}}}_{\cos\theta=\xi}$
$\displaystyle=\sum_{l=1}^{r}\mathrm{Re}\left[\mathcal{R}_{l}^{*}(q,v)\mathscr{B}_{l}(\bm{q})\right]\,,$
(3.33)

$\theta$ and $\phi$ are polar and azimuthal angle of the velocity vector
$\bm{v}$ in a frame with z-axis in the direction of $\bm{q}$,

$\displaystyle\xi=\frac{q}{2m_{\rm DM}v}+\frac{\Delta E}{qv}\,,$ (3.34)

and, finally, a single bar denotes a sum (average) over final (initial)
DM/electron spin configurations. In the evaluation of the azimuthal average in
Eq. (3.33), we use the following identifies [2]

$\displaystyle\left.\frac{1}{2\pi}\int_{0}^{2\pi}\mathrm{d}\phi\,\bm{v}_{\rm
el}^{\perp}\cdot\bm{A}\right|_{\cos\theta=\xi}=\frac{\bm{q}}{m_{e}}\cdot\bm{A}\left(\frac{|\bm{q}|}{m_{e}}\right)^{-2}\frac{\bm{q}}{m_{e}}\cdot\bm{v}_{\rm
el}^{\perp}$
$\displaystyle\left.\frac{1}{2\pi}\int_{0}^{2\pi}\mathrm{d}\phi\,\bm{v}_{\rm
el}^{\perp}\cdot\left(\frac{\bm{q}}{m_{e}}\times\bm{A}\right)\right|_{\cos\theta=\xi}=0\,,$
(3.35)

where $\bm{A}$ is an arbitrary three-dimensional vector. We refer to [2] for a
derivation of these equations. In our numerical applications, we will be
interested in the non-relativistic limit of Eq. (3.33),

$\displaystyle\overline{\overline{|\mathcal{M}_{1\rightarrow
2}(\bm{v},\bm{q})|^{2}}}$
$\displaystyle\simeq\overline{\overline{|\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})|^{2}_{\bm{\ell}=\mathbf{0}}|f_{1\rightarrow 2}(\bm{q})|^{2}}}$
$\displaystyle+2m_{e}\mathrm{Re}\left[\overline{\overline{\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})_{\bm{\ell}=0}f_{1\rightarrow
2}(\bm{q})[\nabla_{\bm{\ell}}\mathcal{M}^{*}(\bm{q},\bm{v}^{\perp}_{\rm
el})]_{\bm{\ell}=0}\cdot\bm{f}^{*}_{1\rightarrow 2}(\bm{q})}}\right]$
$\displaystyle+m_{e}^{2}\overline{\overline{|\nabla_{\bm{\ell}}\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})]_{\bm{\ell}=0}\cdot\bm{f}_{1\rightarrow 2}(\bm{q})|^{2}}}\,.$ (3.36)

###  3.2 Tree-level application to specific models

In this section, we calculate the DM response functions $\mathcal{R}_{l}(q,v)$
associated with the models for DM-electron interactions of Sec. 2. We refer to
App. A for a detailed account of these calculations.

$\mathcal{O}_{1}=\mathds{1}_{\rm DM}\mathds{1}_{e}$ |
---|---
$\mathcal{O}_{4}=\bm{S}_{\rm DM}\cdot\bm{S}_{e}$ | $\mathcal{O}_{11}=i\bm{S}_{\rm DM}\cdot\frac{\bm{q}}{m_{e}}\mathds{1}_{e}$
$\mathcal{O}_{5}=i\bm{S}_{\rm DM}\cdot\left(\frac{\bm{q}}{m_{e}}\times\bm{v}^{\bot}\right)\mathds{1}_{e}$ | $\mathcal{O}_{14}=i\left(\bm{S}_{\rm DM}\cdot\frac{\bm{q}}{m_{e}}\right)\left(\bm{S}_{e}\cdot\bm{v}^{\bot}\right)$
$\mathcal{O}_{6}=\left(\bm{S}_{\rm DM}\cdot\frac{\bm{q}}{m_{e}}\right)\left(\bm{S}_{e}\cdot\frac{\bm{q}}{m_{e}}\right)$ | $\mathcal{O}_{15}=-\left(\bm{S}_{\rm DM}\cdot\frac{\bm{q}}{m_{e}}\right)\left[\left(\bm{S}_{e}\times\bm{v}^{\bot}\right)\cdot\frac{\bm{q}}{m_{e}}\right]$
$\mathcal{O}_{7}=\bm{S}_{e}\cdot\bm{v}^{\bot}\mathds{1}_{\rm DM}$ | $\mathcal{O}_{17}=i\frac{\bm{q}}{m_{e}}\cdot\bm{\mathcal{S}}\cdot\bm{v}^{\bot}\mathds{1}_{e}$
$\mathcal{O}_{8}=\bm{S}_{\rm DM}\cdot\bm{v}^{\bot}\mathds{1}_{e}$ | $\mathcal{O}_{18}=i\frac{\bm{q}}{m_{e}}\cdot\bm{\mathcal{S}}\cdot\bm{S}_{e}$
$\mathcal{O}_{9}=i\bm{S}_{\rm DM}\cdot\left(\bm{S}_{e}\times\frac{\bm{q}}{m_{e}}\right)$ | $\mathcal{O}_{19}=\frac{\bm{q}}{m_{e}}\cdot\bm{\mathcal{S}}\cdot\frac{\bm{q}}{m_{e}}\mathds{1}_{e}$
$\mathcal{O}_{10}=i\bm{S}_{e}\cdot\frac{\bm{q}}{m_{e}}\mathds{1}_{\rm DM}$ | $\mathcal{O}_{20}=\left(\bm{S}_{e}\times\frac{\bm{q}}{m_{e}}\right)\cdot\bm{\mathcal{S}}\cdot\frac{\bm{q}}{m_{e}}$
Table 1: List of operators that can contribute to the transition amplitude for
the models considered here. We denote by $\mathds{1}_{\rm DM}$
($\mathds{1}_{e}$) the identity in the DM (electron) spin space. Similarly,
$\bm{S}_{\rm DM}$ ($\bm{S}_{e}$) is the DM (electron) spin operator. For
spin-0 DM, $\langle\bm{S}_{\rm DM}\rangle=0$. For spin-1/2 DM
$\langle\bm{S}_{\rm
DM}\rangle=\delta^{r^{\prime}r}\xi^{s^{\prime}}\bm{\sigma}\xi^{s}/2$, where
$r$ ($s$) labels the initial electron (DM) spin configuration, $r^{\prime}$
($s^{\prime}$) identifies the final electron (DM) spin configuration and
$\bm{=}(\sigma_{1},\sigma_{2},\sigma_{3})$ consists of the three Pauli
matrices. For spin-1 DM, $\langle\bm{S}_{\rm
DM}\rangle=-i\delta^{r^{\prime}r}\bm{e}^{\prime}_{s^{\prime}}\times\bm{e}_{s}$,
where $\bm{e}_{s}$ and $\bm{e}^{\prime}_{s^{\prime}}$ are three-dimensional
polarisation vectors. Finally,
$\langle\bm{\mathcal{S}}\rangle=\delta^{r^{\prime}r}(e_{si}e^{\prime}_{s^{\prime}j}+e_{sj}e^{\prime}_{s^{\prime}i})/2$.
See App. B for further details.

####  3.2.1 Scalar dark matter

In the case of scalar DM with a scalar mediator, the free amplitude for DM-
electron scattering is

$\displaystyle
i\mathcal{M}=ig_{1}m_{S}\left[\frac{i}{(k^{\prime}-k)^{2}-m_{\phi}^{2}}\right]\left[\overline{u}^{r^{\prime}}(k^{\prime})i(h_{1}+ih_{2}\gamma^{5})u^{r}(k)\right]\,,$
(3.37)

where $u^{r}(k)$, $r=1,2$ and $\overline{u}^{r^{\prime}}(k^{\prime})$,
$r^{\prime}=1,2$ are free Dirac spinors for the initial and final state
electron, while $k$ and $k^{\prime}$ are four-momenta. In the non-relativistic
limit, Eq. (3.37) can be written as

$\displaystyle i\mathcal{M}\simeq
i2m_{S}m_{e}\left(\frac{g_{1}h_{1}}{|\bm{q}|^{2}+m_{\phi}^{2}}\langle\hat{\mathcal{O}}_{1}\rangle-\frac{g_{1}h_{2}}{|\bm{q}|^{2}+m_{\phi}^{2}}\langle\hat{\mathcal{O}}_{10}\rangle\right)\,,$
(3.38)

where the operators $\hat{\mathcal{O}}_{1}$ and $\hat{\mathcal{O}}_{10}$ are
given in Tab. 1, and angle brackets in general denote an expectation value
between initial, $s/r$, and final, $s^{\prime}/r^{\prime}$, DM/electron spin
configurations. For example,
$\langle\hat{\mathcal{O}}_{1}\rangle=\delta^{r^{\prime}r}\delta^{s^{\prime}s}$.
For scalar DM $\delta^{s^{\prime}s}\rightarrow 1$, while $s,s^{\prime}=1,2$
($s,s^{\prime}=1,2,3$) for spin-1/2 (spin-1) DM. By inserting Eq. (3.38) into
Eq. (3.4), we find the modulus squared transition amplitude

$\displaystyle\overline{|M_{1\rightarrow
2}|^{2}}=\left(c_{1}^{2}+\frac{1}{4}c_{10}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\right)\mathscr{B}_{1}\,,$
(3.39)

where the coupling constants $c_{1}$ and $c_{10}$ are related to the coupling
constants and masses of the underlying model as in Tab. 2. Notice the
$1/\Gamma(u,t)$ factor in the definition of $c_{1}$ and $c_{10}$, where

$\displaystyle\Gamma(u,t)^{-1}\equiv 4t(1+u^{2})^{-1}\,,$ (3.40)

arising from the mediator propagator. Comparing Eq. (3.39) with the second
line in Eq. (3.33), we find

$\displaystyle\mathcal{R}_{1}(q,v)=\left(c_{1}^{2}+\frac{1}{4}c_{10}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\right)\,,$
(3.41)

while
$\mathcal{R}_{2}(q,v)=\mathcal{R}_{3}(q,v)=\mathcal{R}_{4}(q,v)=\mathcal{R}_{5}(q,v)=0$.
The case of spin 0 DM with a scalar mediator is thus not characterised by new
features, since the observable electron transition rate can in this case be
computed by using the response function $\overline{\mathcal{W}}_{1}$, as in
the dark photon model.

As a first non trivial example where the ionisation and crystal form factors
commonly used in the case of the dark photon model provide a poor description
of the rate of electron transitions induced by DM scattering, let us consider
the case of scalar DM with a vector mediator. For this model, the free
amplitude for DM-electron scattering can be written as follows

$\displaystyle
i\mathcal{M}=ig_{4}P_{\mu}\left[\frac{-ig^{\mu\nu}}{(k^{\prime}-k)^{2}-m_{G}^{2}}\right]\left[\overline{u}^{r^{\prime}}(k^{\prime})i(h_{3}\gamma_{\nu}+ih_{4}\gamma_{\nu}\gamma^{5})u^{r}(k)\right]\,,$
(3.42)

where $P=p+p^{\prime}$, and $p$ ($p^{\prime}$) is the initial (final) DM four-
momentum. In the non-relativistic limit,

$\displaystyle i\mathcal{M}\simeq
i2m_{S}m_{e}\left(-\frac{2g_{4}h_{3}}{|\bm{q}|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{1}\rangle+\frac{4g_{4}h_{4}}{|\bm{q}|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{7}\rangle\right)\,,$
(3.43)

where the operators $\hat{\mathcal{O}}_{1}$ and $\hat{\mathcal{O}}_{7}$ are
listed in Tab. 1. For the azimuthal average of the squared modulus of the
transition amplitude $\mathcal{M}_{1\rightarrow 2}$ (defined in Eq. (3.33) and
used in Eq. (3.32)), we now find

$\displaystyle\overline{\overline{|\mathcal{M}_{1\rightarrow
2}(\bm{v},\bm{q})|^{2}}}=\left(c_{1}^{2}+\frac{1}{4}c_{7}^{2}\,|\bm{v}^{\perp}_{\rm
el}|^{2}\right)\mathscr{B}_{1}-\frac{1}{2}c_{7}^{2}\,\left(\frac{\bm{q}}{m_{e}}\cdot\bm{v}^{\perp}_{\rm
el}\right)\frac{m_{e}^{2}}{|\bm{q}|^{2}}\,\mathscr{B}_{2}+\frac{1}{4}c_{7}^{2}\,\mathscr{B}_{3}\,,$
(3.44)

which implies

$\displaystyle\mathcal{R}_{1}(q,v)$
$\displaystyle=\left(c_{1}^{2}+\frac{1}{4}c_{7}^{2}\,|\bm{v}^{\perp}_{\rm
el}|^{2}\right)$ $\displaystyle\mathcal{R}_{2}(q,v)$
$\displaystyle=-\frac{1}{2}c_{7}^{2}\left(\frac{\bm{q}}{m_{e}}\cdot\bm{v}^{\perp}_{\rm
el}\right)\frac{m_{e}^{2}}{|\bm{q}|^{2}}$ $\displaystyle\mathcal{R}_{3}(q,v)$
$\displaystyle=\frac{1}{4}c_{7}^{2}\,,$ (3.45)

with $\mathcal{R}_{4}(q,v)=\mathcal{R}_{5}(q,v)=0$. Tab. 2 relates $c_{1}$ and
$c_{7}$ to the parameters of the underling DM model. As in the previous
example, in Tab. 2 we include a $1/\Gamma(u,t)$ factor in the definition of
$c_{1}$ and $c_{7}$.

Scalar Mediator, | Vector Mediator,
---|---
$u=|\bm{q}|/m_{\phi}$, $t=m_{S}m_{e}/m^{2}_{\phi}$ |  $u=|\bm{q}|/m_{G}$, $t=m_{S}m_{e}/m^{2}_{G}$
$\displaystyle c_{1}=\frac{1}{2}\frac{g_{1}h_{1}}{\Gamma(u,t)}$ |  |  |  |  | $\displaystyle c_{10}=-\frac{1}{2}\frac{g_{1}h_{2}}{\Gamma(u,t)}$ |  | $\displaystyle c_{1}=-\frac{g_{4}h_{3}}{\Gamma(u,t)}$ |  |  |  |  | $\displaystyle c_{7}=\frac{2g_{4}h_{4}}{\Gamma(u,t)}$ |
Table 2: Relation between the coupling constants in the Lagrangians of Eq.
(2.1) and Eq. (2.2), and the non-relativistic coupling constants in the DM
response functions for, respectively, scalar DM with a scalar mediator, Eq.
(3.41), and scalar DM with a vector mediator, Eq. (3.45). The dimensionless
function $\Gamma(u,t)^{-1}\equiv 4t(1+u^{2})^{-1}$ accounts for the dependence
of the non-relativistic coupling constants on the momentum transfer $\bm{q}$
and the mediator mass, $m_{\phi}$ or $m_{G}$.

####  3.2.2 Fermionic dark matter

The calculation of the DM response functions $\mathcal{R}_{l}(q,v)$ for
fermionic DM proceeds as in the case of spin-0 DM, which we discussed in
detail in the previous subsection. We provide the details of such a
calculation in App. A and App. B. Here, we limit ourselves to listing the
final results, placing the emphasis on models that predict material responses
different from the standard ionisation and crystal form factors. For fermionic
DM with a scalar mediator, we find

$\displaystyle\mathcal{R}_{1}(q,v)=c_{1}^{2}+\frac{1}{16}c_{6}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{4}c_{10}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{1}{4}c_{11}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\,,$
(3.46)

while
$\mathcal{R}_{2}(q,v)=\mathcal{R}_{3}(q,v)=\mathcal{R}_{4}(q,v)=\mathcal{R}_{5}(q,v)=0$.
As in the case of spin-0 DM with a scalar mediator (see Eq. (3.41)), we find
that in a scenario where DM is made of spin-1/2 particles interacting with
electrons through the exchange of a scalar mediator particle the standard
ionisation and crystal form factors can properly describe the response of
materials to an external DM perturbation. In contrast, for fermionic DM with a
vector mediator, we find

$\displaystyle\mathcal{R}_{1}(q,v)$
$\displaystyle=c_{1}^{2}+\frac{3}{16}c_{4}^{2}+\frac{1}{4}c_{7}^{2}|\bm{v}^{\perp}_{\rm
el}|^{2}+\frac{1}{4}c_{8}^{2}|\bm{v}^{\perp}_{\rm
el}|^{2}+\frac{1}{8}c_{9}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}$
$\displaystyle\mathcal{R}_{2}(q,v)$
$\displaystyle=-\left(\frac{\bm{q}}{m_{e}}\cdot\bm{v}^{\perp}_{\rm
el}\right)\frac{m_{e}^{2}}{|\bm{q}|^{2}}\left(\frac{1}{2}c_{7}^{2}+\frac{1}{2}c_{8}^{2}\right)$
$\displaystyle\mathcal{R}_{3}(q,v)$
$\displaystyle=\frac{1}{4}c_{7}^{2}+\frac{1}{4}c_{8}^{2}\,,$ (3.47)

$\mathcal{R}_{4}(q,v)=\mathcal{R}_{5}(q,v)=0$. We therefore conclude that it
is crucial to evaluate the material response functions
$\overline{\mathcal{W}}_{2}(q,\Delta E)$ and
$\overline{\mathcal{W}}_{3}(q,\Delta E)$, in addition to the standard
ionisation and crystal form factors, to consistently compare theoretical
predictions with experimental data in a model with fermionic DM and a vector
mediator. For example, neglecting the $\overline{\mathcal{W}}_{2}(q,\Delta E)$
and $\overline{\mathcal{W}}_{3}(q,\Delta E)$ response functions when computing
an exclusion limit from a null result, or when interpreting a signal at a DM
direct detection experiment, would lead to strongly biased results if nature
is indeed described by fermionic DM with a vector mediator.

Scalar Mediator, | Vector Mediator,
---|---
$u=|\bm{q}|/m_{\phi}$, $t=m_{\chi}m_{e}/m^{2}_{\phi}$ |  $u=|\bm{q}|/m_{G}$, $t=m_{\chi}m_{e}/m^{2}_{G}$ xxxxxxxxxxxxxxx
$\displaystyle c_{1}=\frac{\lambda_{1}h_{1}}{\Gamma(u,t)}$ |  |  | $\displaystyle c_{6}=\frac{\lambda_{2}h_{2}}{\Gamma(u,t)}\frac{m_{e}}{m_{\chi}}$ | $\displaystyle c_{1}=-\frac{\lambda_{3}h_{3}}{\Gamma(u,t)}$ |  |  |  | $\displaystyle c_{4}\ =\frac{4\lambda_{4}h_{4}}{\Gamma(u,t)}$
$\displaystyle c_{10}=-\frac{\lambda_{1}h_{2}}{\Gamma(u,t)}$ |  |  | $\displaystyle c_{11}=\frac{\lambda_{2}h_{1}}{\Gamma(u,t)}\frac{m_{e}}{m_{\chi}}$ | $\displaystyle c_{7}=\frac{2\lambda_{3}h_{4}}{\Gamma(u,t)}$ |  |  |  | $\displaystyle c_{8}=-\frac{2\lambda_{4}h_{3}}{\Gamma(u,t)}$
|  |  |  | $\displaystyle c_{9}=\frac{2\lambda_{3}h_{4}}{\Gamma(u,t)}\frac{m_{e}}{m_{\chi}}+\frac{2\lambda_{4}h_{3}}{\Gamma(u,t)}$ |  |  |
Table 3: Relation between the coupling constants in the Lagrangians of Eq.
(2.3) and Eq. (2.4) and the non-relativistic coupling constants in the DM
response functions for, respectively, fermionic DM with a scalar mediator, Eq.
(3.46), and fermionic DM with a vector mediator, Eq. (3.47). Here, we
introduced the dimensionless function $\Gamma(u,t)^{-1}\equiv
4t(1+u^{2})^{-1}$ to account for the dependence of the non-relativistic
coupling constants on $\bm{q}$ and the mediator mass.

####  3.2.3 Vector dark matter

In the case of spin-1 DM with a scalar mediator we find that just a single
response function contributes to the observed rate of electron transitions
induced by DM scattering in a material, namely

$\displaystyle\mathcal{R}_{1}(q,v)=c_{1}^{2}+\frac{1}{4}c_{10}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\,,$
(3.48)

while
$\mathcal{R}_{2}(q,v)=\mathcal{R}_{3}(q,v)=\mathcal{R}_{4}(q,v)=\mathcal{R}_{5}(q,v)=0$.
We therefore find that in all models where the interaction between DM and the
electrons bound to a material is mediated by the exchange of a scalar mediator
the rate of DM-induced electron transitions can be computed by using the
standard ionisation and crystal form factors. Turning this argument around, we
find that experimental evidence for contributions to the theoretical rate
$\mathscr{R}_{\rm theory}$ from $\mathcal{W}_{l}$, with $l\neq 1$, would point
towards DM-electron interactions mediated by particles with spin different
from zero.

Table 4: The same as in Tab. 2 and Tab. 3, now for vector DM. Scalar Mediator,
$u=|\bm{q}|/m_{\phi}$, $t=m_{X}m_{e}/m^{2}_{\phi}$
---
$\displaystyle c_{1}=\frac{1}{2}\frac{b_{1}h_{1}}{\Gamma(u,t)}$ | $\displaystyle c_{10}=-\frac{1}{2}\frac{b_{1}h_{2}}{\Gamma(u,t)}$
Vector Mediator, $u=|\bm{q}|/m_{G}$, $t=m_{X}m_{e}/m^{2}_{G}$
$\displaystyle c_{1}=-\frac{b_{5}h_{3}}{\Gamma(u,t)}$ |
$\displaystyle c_{4}=\frac{\mathrm{Im}(b_{6})h_{3}}{\Gamma(u,t)}\frac{m_{e}}{2m_{X}}\frac{|\bm{q}\,|^{2}}{m_{e}^{2}}-\frac{2\mathrm{Re}(b_{7})h_{4}}{\Gamma(u,t)}\phantom{xxxxxxx}$ | $\displaystyle c_{11}=\frac{\mathrm{Im}(b_{7})h_{3}}{\Gamma(u,t)}\frac{m_{e}}{2m_{X}}$
$\displaystyle c_{5}=\frac{\mathrm{Im}(b_{6})h_{3}}{\Gamma(u,t)}\frac{m_{e}}{2m_{X}}$ | $\displaystyle c_{14}=-\frac{\mathrm{Im}(b_{7})h_{4}}{\Gamma(u,t)}\frac{m_{e}}{m_{X}}$
$\displaystyle c_{6}=-\frac{\mathrm{Im}(b_{6})h_{3}}{\Gamma(u,t)}\frac{m_{e}}{2m_{X}}$ | $\displaystyle c_{17}=-\frac{\mathrm{Re}(b_{6})h_{3}}{\Gamma(u,t)}\frac{m_{e}}{m_{X}}$
$\displaystyle c_{7}=\frac{2b_{5}h_{4}}{\Gamma(u,t)}$ | $\displaystyle c_{18}=\frac{2\mathrm{Re}(b_{6})h_{4}}{\Gamma(u,t)}\frac{m_{e}}{m_{X}}$
$\displaystyle c_{8}=\frac{\mathrm{Re}(b_{7})h_{3}}{\Gamma(u,t)}$ | $\displaystyle c_{19}=\frac{\mathrm{Im}(b_{6})h_{3}}{\Gamma(u,t)}\frac{m_{e}^{2}}{2m_{X}^{2}}$
$\displaystyle c_{9}=\frac{\mathrm{Im}(b_{6})h_{4}}{\Gamma(u,t)}\frac{m_{e}}{m_{X}}-\frac{\mathrm{Re}(b_{7})h_{3}}{\Gamma(u,t)}$ | $\displaystyle c_{20}=-\frac{\mathrm{Re}(b_{6})h_{3}}{\Gamma(u,t)}\frac{m_{e}}{m_{X}}-\frac{\mathrm{Im}(b_{7})h_{4}}{\Gamma(u,t)}\frac{m_{e}^{2}}{m_{X}^{2}}$

The case of spin-1 DM with a vector mediator exhibits a richer phenomenology,
with several material response functions possibly contributing to the
theoretical rate $\mathscr{R}_{\rm theory}$. Specifically,

$\displaystyle\mathcal{R}_{1}(q,v)$
$\displaystyle=c_{1}^{2}+\frac{1}{2}c_{4}^{2}+\frac{2}{3}c_{5}^{2}\left|\frac{\bm{q}}{m_{e}}\times\bm{v}_{\rm
el}^{\perp}\right|^{2}+\frac{1}{6}c_{6}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{4}c_{7}^{2}|\bm{v}_{\rm
el}^{\perp}|^{2}+\frac{2}{3}c_{8}^{2}|\bm{v}_{\rm el}^{\perp}|^{2}$
$\displaystyle\qquad\,+\frac{1}{3}c_{9}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{2}{3}c_{11}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{1}{6}c_{14}^{2}|\bm{v}_{\rm
el}^{\perp}|^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}$
$\displaystyle\qquad\,+\frac{1}{6}c_{17}^{2}\left(\left|\frac{\bm{q}}{m_{e}}\cdot\bm{v}_{\rm
el}^{\perp}\right|^{2}+\frac{|\bm{q}|^{2}}{m_{e}^{2}}|\bm{v}_{\rm
el}^{\perp}|^{2}\right)+\frac{1}{6}c_{18}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{1}{3}c_{19}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}$
$\displaystyle\qquad\,+\frac{1}{12}c_{20}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{3}c_{4}c_{6}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{2}{3}c_{1}c_{19}\frac{|\bm{q}|^{2}}{m_{e}^{2}}$
$\displaystyle\mathrm{Re}[\mathcal{R}_{2}(v,q)]$
$\displaystyle=-\left(\frac{\bm{q}}{m_{e}}\cdot\bm{v}_{\rm
el}^{\perp}\right)\left[\frac{m_{e}^{2}}{|\bm{q}|^{2}}\left(\frac{1}{2}c_{7}^{2}+\frac{4}{3}c_{8}^{2}\right)+\frac{1}{3}c_{14}^{2}+\frac{2}{3}c_{17}^{2}\right]$
$\displaystyle\mathrm{Im}[\mathcal{R}_{2}(v,q)]$
$\displaystyle=\frac{1}{2}c_{7}c_{10}-\frac{1}{3}c_{4}c_{14}+\frac{4}{3}c_{8}c_{11}-\frac{1}{3}c_{6}c_{14}\frac{|\bm{q}|^{2}}{m_{e}^{2}}$
$\displaystyle\mathcal{R}_{3}(v,q)$
$\displaystyle=\frac{1}{4}c_{7}^{2}+\frac{2}{3}c_{8}^{2}+\frac{|\bm{q}|^{2}}{m_{e}^{2}}\left(\frac{2}{3}c_{5}^{2}+\frac{1}{6}c_{14}^{2}+\frac{1}{6}c_{17}^{2}\right)$
$\displaystyle\mathcal{R}_{4}(v,q)$
$\displaystyle=-\frac{2}{3}c_{5}^{2}+\frac{1}{6}c_{17}^{2}\,,$ (3.49)

with $R_{5}(q,v)=0$. We refer to App. A and App. B for all details underling
the derivation of Eq. (3.49). Interestingly, none of the models considered
here can generate the response function $\overline{\mathcal{W}}_{5}$.

DM00 (70,80) t1,t2 b1,b2 fermiont1,v1 fermionv2,t2 dots,label=$\phi\ $,
label.side=right, tension=.5v1,v3 dots,label=$\ \phi$,tension=.5v3,v2
fermion,label=$e^{-}$,tension=.5v1,v2 dashesb1,v3,b2 $e^{-}$t1 $e^{-}$t2 $S$b1
$S$b2

(a)

DM01 (70,80) t1,t2 b1,b2 fermiont1,v1 fermionv2,t2 dbl_curly,label=$G\ \ $,
label.side=right, tension=.5v1,v3 dbl_curly,label=$\ \ G$,tension=.5v3,v2
fermion,label=$e^{-}$,tension=.5v1,v2 dashesb1,v3,b2 $e^{-}$t1 $e^{-}$t2 $S$b1
$S$b2

(b)

DM10 (70,80) t1,t2 b1,b2 fermiont1,v1 fermionv2,t2 dots,label=$\phi\ $,
label.side=right, tension=.5v1,v3 dots,label=$\ \phi$,tension=.5v3,v2
fermion,label=$e^{-}$,tension=.5v1,v2 curlyb1,v3,b2 $e^{-}$t1 $e^{-}$t2 $X$b1
$X$b2

(c)

DM11 (70,80) t1,t2 b1,b2 fermiont1,v1 fermionv2,t2 dbl_curly,label=$G\ \
$,label.side=right,tension=.5v1,v3 dbl_curly,label=$\ \ G$,tension=.5v3,v2
fermion,label=$e^{-}$,tension=.5v1,v2 curlyb1,v3,b2 $e^{-}$t1 $e^{-}$t2 $X$b1
$X$b2

(d)

Figure 3: One-loop Feynman diagrams for DM-electron scattering assuming: 1)
spin-0 DM with a scalar mediator (top left), 2) spin-0 DM with a vector
mediator (top right), 3) spin-1 DM with a scalar mediator (bottom left) and,
finally, 4) spin-1 DM with a vector mediator (bottom right). The particles in
the processes are denoted by the letters used for the associated fields.

###  3.3 One-loop application to specific models

So far, we focused on DM-electron interactions arising at tree level. These
required a cubic DM-DM-mediator interaction vertex in the Lagragian. It is
interesting to note that even when such cubic vertices are absent at tree
level, DM-electron interactions can be generated at the one-loop level via a
quartic DM-DM-mediator-mediator vertex, as shown in Fig. 3 in four specific
cases: 1) scalar DM with a scalar mediator, 2) scalar DM with a vector
mediator, 3) vector DM with a scalar mediator, and, finally 4) vector DM with
a vector mediator. In App. C, we list the free electron amplitude for the four
diagrams in Fig. 3. They are all proportional to
$\langle\hat{\mathcal{O}}_{1}\rangle$, and therefore contribute to the
observable rate of electron transitions induced by DM scattering in a detector
material, $\mathscr{R}_{\rm theory}$, via the response function
$\overline{\mathcal{W}}_{1}(q,\Delta E)$ only. From the point of view of
material physics and for the purposes of our work, this family of DM-electron
interaction models is therefore equivalent to the dark photon model. For this
reason, here we do not investigate the effects of loop-induced DM-electron
interactions any further.

##  4 Comparison with experimental data

In this section, we compare the rate of electron transitions predicted from
the models of Sec. 2 with experimental data released by the XENON10 [43],
XENON1T[15], EDELWEISS [23] and SENSEI [24] collaborations. As we will see
below, for xenon experiments we adopt the data analysis strategies we used in
[1]. For the analysis of results from crystal detectors, we rely on [2].

XENON10xxxxxx |  |  | XENON1Txxxxxx |
---|---|---|---|---
bin [S2] | observed events |  | bin [S2] | observed events
[14,41) | 126 |  |  [150,200) | 8
[41,68) | 60 |  |  [200,250) | 7
[68,95) | 12 |  |  [250,300) | 2
[95,122) | 3 |  |  [300,350) | 1
[122,149) | 2 |  | - | -
[149,176) | 0 |  | - | -
[176,203) | 2 |  | - | -
Table 5: Events in selected S2 bins observed by XENON10 [43] (left) and
XENON1T [15] (right).

###  4.1 Direct detection data

The XENON10 and XENON1T experiments are dual phase time projection chambers
(TPCs) filled with liquid and gaseous xenon. They have been searching for sub-
GeV DM particles in events where xenon atoms are ionised by DM-electron
scattering in the liquid xenon detector component. In each such event, primary
electrons from the ionised atoms are expected to propagate in the liquid
xenon, producing a number $n_{e}$ of secondary electrons. These secondary
electrons are then extracted from the liquid target by means of an external
electric field, and drifted into the gaseous xenon detector component, where
they eventually generate a scintillation signal. This signal is finally
recorded by photomultiplier tubes through the observation of ${\rm S2}$
associated photoelectrons [43, 15]. Since ${\rm S2}$ is the actually observed
quantity in xenon TPCs experiments, we need to rewrite our theoretical rate,
Eq. (3.9), as a function of ${\rm S2}$.

Let us denote by $\mathscr{P}(n_{e}|E_{e})$ the probability density to produce
$n_{e}$ secondary electrons in a DM-electron scattering event when the energy
of the primary electron is $E_{e}=k^{\prime 2}/(2m_{e})$. Furthermore, let
$\mathscr{P}({\rm S2}|n_{e})$ be the probability density to observe ${\rm S2}$
photoelectrons when $n_{e}$ electrons drifted from the liquid xenon target
into the gaseous phase. Then, the rate of DM-electron scattering events with a
number of observed photoelectrons between ${\rm S2}$ and ${\rm S2}+{\rm d}{\rm
S2}$ is [1]

$\displaystyle{\rm d}\mathscr{R}_{\rm obs}$ $\displaystyle={\rm d}{\rm
S2}\,\varepsilon({\rm S2})\frac{n_{\rm DM}}{128\pi m_{\rm
DM}^{2}m_{e}^{2}}\int\mathrm{d}(\ln\Delta
E)\int\mathrm{d}q\,q\,\widehat{\eta}\left(q,\Delta
E\right)\sum_{n_{e}=1}^{\infty}\mathscr{P}({\rm S2}|n_{e})$
$\displaystyle\times\sum_{l=1}^{r}\mathrm{Re}\left[\mathcal{R}_{l}^{*}(q,v)\overline{\mathcal{W}}_{l}(q,\Delta
E;n_{e})\right]\,,$ (4.1)

where $\varepsilon({\rm S2})$ is the experiment-specific detection efficiency,
and we introduced the modified response

$\displaystyle\overline{\mathcal{W}}_{l}(q,\Delta E;n_{e})\equiv$
$\displaystyle\int\mathrm{d}\Omega_{q}\,\mathcal{W}_{l}(\bm{q},\Delta
E;n_{e})\,,$ (4.2)

with

$\displaystyle\mathcal{W}_{l}(\bm{q},\Delta E;n_{e})$
$\displaystyle=\frac{2}{\pi}\Delta
E\sum_{n,\ell,m}\sum_{\ell^{\prime},m^{\prime}}\int\frac{Vk^{\prime
2}\,\mathrm{d}k^{\prime}}{(2\pi)^{3}}\,\mathscr{B}_{l}\,\mathscr{P}(n_{e}|E_{e})\,\delta(\Delta
E-E_{n\ell m}+E_{k^{\prime}\ell^{\prime}m^{\prime}})\,.$ (4.3)

It is important to note that Eq. (4.1) gives the rate of observable events per
atom in the detector, not the total rate. The total rate of observable events
per unit detector mass is found by multiplying Eq. (4.1) by the number of
xenon atoms in the detector, $n_{\rm Xe}$, and then dividing by the detector
mass, $M_{\rm det}=n_{\rm Xe}m_{\rm Xe}$, where $m_{\rm Xe}$ is the mass of a
xenon atom. One finds,

$\displaystyle{\rm d}\widetilde{\mathscr{R}}_{\rm atoms}\equiv{\rm
d}\mathscr{R}_{\rm obs}/m_{\rm Xe}\,.$ (4.4)

We refer to [1] for further details on $\mathscr{P}(n_{e}|E_{e})$,
$\mathscr{P}({\rm S2}|n_{e})$ and $\varepsilon({\rm S2})$. By integrating Eq.
(4.4) over the ${\rm S2}$ bins in Tab. 5, and multiplying the result by the
experimental exposure (15 kg days for XENON10 and 80755 kg days for XENON1T),
we obtain the expected number of DM signal events in the given ${\rm S2}$ bin.
We can then use the statistical methods reviewed below to compare our
predictions with the XENON10 [43] and XENON1T [15] results we report in Tab.
5.

EDELWEISS |  |  |  | SENSEI |  |
---|---|---|---|---|---|---
$[Q]$ | obs. events | exp. [g-day] |  | $[Q]$ | obs. events | exp. [g-day]
1 | 5814 | 3.23 |  | 1 | 758 | 1.38
2 | 44706 | 17.76 |  | 2 | 5 | 2.09
3 | 2718 | 80.72 |  | 3 | 0 | 9.03
4 | 227 | 80.72 |  | 4 | 0 | 9.10
Table 6: Events corresponding to $Q=1,2,3,4$ observed by the EDELWEISS [23]
and SENSEI [24] experiments, together with the associated $Q$ dependent
effective exposure expressed in [g-day].

A second class of experiments of interest for this work includes EDELWEISS
[23] and SENSEI [24]. These experiments operate, respectively, germanium and
silicon crystal detectors, and search for sub-GeV DM particles in events
associated with the production of electron-hole pairs resulting from
electronic transitions induced by DM-electron scattering in the detector
crystals. In order to compare our theoretical rate, Eq. (3.9), with EDELWEISS
[23] and SENSEI [24] results, we first calculate the total event rate per unit
detector mass by dividing Eq. (3.9) by the detector mass, $m_{\rm det}=N_{\rm
cell}m_{\rm cell}$, where $m_{\rm cell}$ is the mass of a unit cell,

$\displaystyle{\rm d}\widetilde{\mathscr{R}}_{\rm crystals}\equiv{\rm
d}\mathscr{R}_{\rm theory}/m_{\rm det}\,,$ (4.5)

and then convert the deposited energy $\Delta E$ in Eq. (4.5) into a number of
electron-hole pairs produced in a DM-electron scattering event, which we
denote by $Q$. We convert $\Delta E$ into $Q$ using [18]

$\displaystyle Q(\Delta E)=1+\lfloor(\Delta E-E_{\rm
gap})/\varepsilon\rfloor\,,$ (4.6)

where $\lfloor\cdot\rfloor$ is the floor function. The observed band-gap and
mean energy per electron-hole pair are $E_{\rm gap}=0.67$ eV and
$\varepsilon=3.0$ eV for germanium [23], and $E_{\rm gap}=1.20$ eV and
$\varepsilon=3.8$ eV for silicon [24]. We can then calculate the expected
number of DM signal events corresponding to $Q=n$ electron-hole pairs, where
$n\in\mathbb{N}$, by performing the $\Delta E$ integral in Eq. (3.9) over the
interval [$E_{\rm gap}+(n-1)\varepsilon,E_{\rm gap}+n\varepsilon$), and
multiplying the result by the $Q$ dependent effective exposure (that is, the
exposure times the detection efficiency). We report the latter, together with
the events actually observed by the EDELWEISS [23] and SENSEI [24] experiments
in Tab. 6.

Let us now denote by $\mathcal{S}_{i}$ the expected number of DM signal events
in the $i$-th S2 bin, for XENON10 and XENON1T, or in the $i$-th $Q$ bin, for
EDELWEISS and SENSEI. Furthermore, let $\mathcal{D}_{i}$ be the number of
actually observed events in the same (S2 or $Q$) bin. For each of the four
experiments considered here, we calculate 90% confidence level (C.L.)
exclusion limits on the strength of a given DM-electron interaction by
requiring that $S_{i}\leq S_{i}^{90\%}$ for all bins, where

$\displaystyle\sum_{k=\mathcal{D}_{i}+1}^{\infty}\mathcal{P}(k\,|\,\mathcal{S}_{i}^{90\%})=0.9$
(4.7)

defines $S^{90\%}_{i}$ and the number of counts in each bin is assumed to obey
Poisson statistics. Here, $\mathcal{P}(k\,|\,\mathcal{S}_{i}^{90\%})$ is the
Poisson probability of observing $k$ events when an average of
$\mathcal{S}_{i}^{90\%}$ is expected.

By computing exclusion limits using Eq. (4.7), we assume that the observed S2
or $Q$ events do not allow one to reject the “null hypothesis” in favour of a
DM discovery. For a detailed description of possible background contributions
to $\mathcal{D}_{i}$, we refer to the experimental works [23, 24, 43, 15]. For
example, environmental backgrounds and dark current events from thermal
excitations are the main contributions to $\mathcal{D}_{i}$ in the case of
SENSEI. The exclusion limits that we obtain from Eq. (4.7) are conservative,
as they do not rely on any background subtraction.

Simplified models and underlying Lorentz structures |
---|---
Lagrangian | Lorentz structure
$\mathcal{L}_{g_{4}h_{4}}^{(0)}=-ig_{4}(S^{\dagger}\partial_{\mu}S-\partial_{\mu}S^{\dagger}S)G^{\mu}-h_{4}(\overline{e}\gamma_{\mu}\gamma^{5}e)G^{\mu}$ | $\partial_{-}\otimes A$
$\mathcal{L}_{\lambda_{3}h_{4}}^{(1/2)}=-\lambda_{3}\overline{\chi}\gamma^{\mu}\chi G_{\mu}-h_{4}\overline{e}\gamma_{\mu}\gamma^{5}eG^{\mu}$ | $V\otimes A$
$\mathcal{L}_{\lambda_{4}h_{3}}^{(1/2)}=-\lambda_{4}\overline{\chi}\gamma^{\mu}\gamma^{5}\chi G_{\mu}-h_{3}\overline{e}\gamma_{\mu}eG^{\mu}$ | $A\otimes V$
$\mathcal{L}_{b_{5}h_{4}}^{(1)}=-ib_{5}\left[X_{\nu}^{\dagger}\partial_{\mu}X^{\nu}-(\partial_{\mu}X^{\dagger\nu})X_{\nu}\right]G^{\mu}-h_{4}G_{\mu}\overline{e}\gamma^{\mu}\gamma^{5}e$ | $\partial_{-}\otimes A$
$\mathcal{L}_{\mathrm{Re}(b_{6})h_{3}}^{(1)}=-\mathrm{Re}(b_{6})\partial_{\nu}\left(X^{\dagger\nu}X_{\mu}+X^{\dagger}_{\mu}X^{\nu}\right)G^{\mu}-h_{3}G_{\mu}\overline{e}\gamma^{\mu}e$ | $\partial_{+}\otimes V$
$\mathcal{L}_{\mathrm{Im}(b_{6})h_{3}}^{(1)}=-i\mathrm{Im}(b_{6})\partial_{\nu}\left(X^{\dagger\nu}X_{\mu}-X^{\dagger}_{\mu}X^{\nu}\right)G^{\mu}-h_{3}G_{\mu}\overline{e}\gamma^{\mu}e$ | $\partial_{-}\otimes V$
$\mathcal{L}_{\mathrm{Re}(b_{7})h_{3}}^{(1)}=-\mathrm{Re}(b_{7})\varepsilon_{\mu\nu\rho\sigma}\left(X^{\dagger\mu}\partial^{\nu}X^{\rho}+X^{\mu}\partial^{\nu}X^{\dagger\rho}\right)G^{\sigma}-h_{3}G_{\mu}\overline{e}\gamma^{\mu}e$ | $\varepsilon\partial_{+}\otimes V$
$\mathcal{L}_{\mathrm{Im}(b_{7})h_{4}}^{(1)}=-i\mathrm{Im}(b_{7})\varepsilon_{\mu\nu\rho\sigma}\left(X^{\dagger\mu}\partial^{\nu}X^{\rho}-X^{\mu}\partial^{\nu}X^{\dagger\rho}\right)G^{\sigma}-h_{4}G_{\mu}\overline{e}\gamma^{\mu}\gamma^{5}e$ | $\varepsilon\partial_{-}\otimes A$
Table 7: Lagrangians for a selection of simplified models with non-standard
material response functions.

###  4.2 Constraints

We are now ready to set constraints on the models of Sec. 2 from the null
result of XENON10, XENON1T, EDELWEISS and SENSEI. While we have so far
computed the DM response functions one would obtain when all coupling
constants in a given Lagrangian from Sec. 2 are different from zero (Eqs.
(3.41), (3.45), (3.46), (3.47), (3.48) and (3.49)), here we restrict ourselves
to a subset of “simplified models” each characterised by four parameters,
namely: the DM and mediator mass, a coupling constant for the DM-DM-mediator
vertex, and a coupling constant for the electron-electron-mediator interaction
vertex. This choice will simplify the discussion while still allowing us to
highlight the importance of the non-standard response functions,
$\overline{\mathcal{W}}_{l}(q,\Delta E)$, $l\neq 1$, in the non-relativistic
modelling of DM-electron scattering in materials.

Specifically, we focus on the simplified models defined in Tab. 7. The
corresponding Lagrangians follow by integration by parts from the Lagrangians
in Sec. 2 when only one pair of coupling constants at the time is assumed to
be different from zero. We focus on these eight simplified models as they
generate at least one non-standard DM response function and, therefore, can
only be investigated within the response function formalism we developed in
[1, 2].

Figure 4: 90% C. L. exclusion limits on the reference cross section
$\sigma_{\rm ref}$ for the simplified models with couplings constants $g_{4}$
and $h_{4}$ (top left panel), $\lambda_{3}$ and $h_{4}$ (top right panel),
$\lambda_{4}$ and $h_{3}$ (bottom left panel) and $b_{5}$ and $h_{4}$ (bottom
right panel). See Tab. 7 for the corresponding Lorentz structures. Solid,
dotted and dashed lines refer to long-range, “intermediate-range” with $m_{\rm
med}=\alpha m_{e}$ and long-range interactions, while blue, orange, green and
maroon lines are associated with the EDELWEISS, SENSEI, XENON10 and XENON1T
experiments, respectively.

Figure 5: Same as Fig. 4, for the simplified models with couplings constants
$\mathrm{Re}(b_{6})$ and $h_{3}$ (top left), $\mathrm{Im}(b_{6})$ and $h_{3}$
(top right), $\mathrm{Re}(b_{7})$ and $h_{3}$ (bottom left) and
$\mathrm{Im}(b_{7})$ and $h_{4}$ (bottom right).

Figure 6: Same as Fig. 4, for the simplified models defined by the couplings
constants $\lambda_{3}$ and $h_{3}$ (top left), $\lambda_{4}$ and $h_{4}$ (top
right), $\mathrm{Re}(b_{6})$ and $h_{4}$ (bottom left) and $\lambda_{2}$ and
$h_{2}$ (bottom right).

The Lorentz structure of the Lagrangians in Tab. 7 is variegated. For example,
$\mathcal{L}_{g_{4}h_{4}}^{(0)}$ is characterised by a derivative coupling
between scalar DM and a vector mediator, and an axial coupling between the
vector mediator and the electron. Symbolically, we denote this Lorentz
structure by $\partial_{-}\otimes A$, where the minus sign in $\partial_{-}$
refers to the minus sign in parenthesis in the first term of
$\mathcal{L}_{g_{4}h_{4}}^{(0)}$. Analogously,
$\mathcal{L}_{\mathrm{Re}(b_{7})h_{3}}^{(1)}$ is characterised by a derivative
coupling between vector DM and a vector mediator, and a vector coupling
between the mediator and the electron. We denote this Lorentz structure by
$\varepsilon\partial_{+}\otimes V$ to emphasise that the derivative coupling
in question includes a Levi-Civita tensor. By adopting the notation outlined
in these examples, Tab. 7 shows the relation between our simplified models and
the underlying Lorentz structures.

For the models in Tab. 7, we calculate 90% C. L. exclusion limits on the
reference cross section,

$\displaystyle\sigma_{\rm ref}=\frac{\mu^{2}_{\rm
DMe}}{\pi}\frac{\mathcal{P}^{2}_{1}\mathcal{P}^{2}_{2}}{\Lambda^{4}}\,,$ (4.8)

from the null result reported by the XENON10 [43], XENON1T[15], EDELWEISS [23]
and SENSEI [24] collaborations. In Eq. (4.8), $\mu_{\rm DMe}$ is the DM-
electron reduced mass, $\mathcal{P}_{1}$ and $\mathcal{P}_{2}$ are the
coupling constants that define the given simplified model, while $\Lambda$ is
a reference mass scale which depends on the assumed value for the mediator
mass, $m_{\rm med}$. For the simplified models in Tab. 7, we consider three
cases for the mediator mass, each associated with a specific choice for
$\Lambda$: 1) $m_{\rm med}\gg|\bm{q}|$, with $\Lambda=m_{\rm med}$, 2) $m_{\rm
med}=\alpha m_{e}$, with $\Lambda=m_{\rm med}$, and, finally, 3) $m_{\rm
med}=0$, with $\Lambda=\alpha m_{e}$. The first (last) case corresponds to
short-range (long-range) interactions, while the second one corresponds to
$m_{\rm med}$ values that are comparable with a “reference” momentum transfer,
which we set to $q_{\rm ref}=\alpha m_{e}$, where $\alpha$ is the fine
structure constant.

Fig. 4 shows the 90% C. L. exclusion limits on the reference cross section
$\sigma_{\rm ref}$ that we obtain for the simplified models defined by the
couplings constants $g_{4}$ and $h_{4}$ (top left panel), $\lambda_{3}$ and
$h_{4}$ (top right panel), $\lambda_{4}$ and $h_{3}$ (bottom left panel) and
$b_{5}$ and $h_{4}$ (bottom right panel). In all panels, exclusion limits are
presented as a function of the DM particle mass, as well as for different
mediator masses and experiments. More specifically, the solid, dotted and
dashed lines in the figure refer to long-range, “intermediate-range” with
$m_{\rm med}=\alpha m_{e}$ and long-range interactions, while blue, orange,
green and maroon lines are associated with the EDELWEISS, SENSEI, XENON10 and
XENON1T experiments, respectively, as indicated in the legends.

There are no appreciable differences in the exclusion limits presented in the
four panels of Fig. 4 because the underlying simplified models either generate
the operator $\mathcal{O}_{7}$ only in the non-relativistic limit, or linear
combinations of the interaction operators $\mathcal{O}_{7}$ and
$\mathcal{O}_{8}$, and $\mathcal{O}_{7}$ and $\mathcal{O}_{9}$. However, the
operator $\mathcal{O}_{7}$ alone generates three non-standard material
response functions.

As expected from Fig. 1 and Fig. 2, and from the fact that xenon ionisation
energies are larger than germanium and silicon band gaps by an order 10
factor, for $m_{\rm DM}$ below about 5 MeV crystal detectors such as EDELWEISS
and SENSEI place stronger constraints on the cross section $\sigma_{\rm ref}$
than the xenon detectors XENON10 and XENON1T in all panels of Fig. 4. In
contrast, XENON10 and XENON1T set stronger bounds on $\sigma_{\rm ref}$ for
values of the DM particle mass that are larger than about 5 MeV because of
their larger experimental exposure.

In the case of XENON10 and XENON1T, exclusion limits for short-range
interactions are always stronger than for long-range interactions, except near
mass threshold (i.e. for $m_{\rm DM}$ close to the smallest testable DM mass),
where they tend to overlap. Furthermore, exclusion limits for short-range and
intermediate-range interactions essentially coincide. Both results are a
consequence of the differential ionisation rate being proportional to
$\sigma_{\rm ref}$, $\sigma_{\rm ref}(1+|\bm{q}|/q_{\rm ref})^{-4}$ and
$\sigma_{\rm ref}q_{\rm ref}^{4}/|\bm{q}|^{4}$, for long-, intermediate- and
short-range interactions, and $|\bm{q}|^{4}/q_{\rm ref}^{4}$ being typically
larger than one for DM particle masses in the range probed by XENON10 and
XENON1T.

In the case of EDELWEISS and SENSEI, the exclusion limits that we find in Fig.
4 for $m_{\rm DM}$ larger than about 2 MeV and for different mediator masses
exhibit the same hierarchy described above for XENON10 and XENON1T, the only
difference being a non negligible separation between the exclusion limits
obtained for short-range and intermediate-range interactions. In contrast, for
DM particle masses below about 2 MeV, long-range interactions are associated
with significantly stronger exclusion limits than short- or intermediate-range
interactions. This result can be explained by recalling the aforementioned
scaling with $|\bm{q}|/q_{\rm ref}$ of the differential ionisation or
excitation rate, and by noticing that for crystal detectors $|\bm{q}|/q_{\rm
ref}$ is smaller than one for sufficiently small DM particle masses, as one
can infer from Figs. 1 and 2.

Fig. 5 shows the 90% C. L. exclusion limits on the reference cross section
$\sigma_{\rm ref}$ that we find for the simplified models with couplings
constants $\mathrm{Re}(b_{6})$ and $h_{3}$ (top left panel),
$\mathrm{Im}(b_{6})$ and $h_{3}$ (top right panel), $\mathrm{Re}(b_{7})$ and
$h_{3}$ (bottom left panel) and $\mathrm{Im}(b_{7})$ and $h_{4}$ (bottom right
panel). The exclusion limits in Fig. 5 exhibit the same hierarchies and
patterns discussed in details for Fig. 4. Importantly, none of the simplified
models underlying Figs. 4 and 5 could have been investigated accurately
without having first computed the material response functions
$\overline{\mathcal{W}}_{2}(q,\Delta E)$ and
$\overline{\mathcal{W}}_{3}(q,\Delta E)$ for atoms and crystals. Furthermore,
we find that for vector DM models with a vector mediator also the material
response function $\overline{\mathcal{W}}_{4}(q,\Delta E)$ is crucial.

We conclude this analysis by validating our results through a comparison with
the exclusion limits associated with models that generate the material
response function $\overline{\mathcal{W}}_{1}(q,\Delta E)$ only. Fig. 6 shows
our 90% C. L. exclusion limits on the reference cross section $\sigma_{\rm
ref}$ for the simplified models that are identified by the couplings constants
$\lambda_{3}$ and $h_{3}$ (top left panel), $\lambda_{4}$ and $h_{4}$ (top
right panel), $\mathrm{Re}(b_{6})$ and $h_{4}$ (bottom left panel) and
$\lambda_{2}$ and $h_{2}$ (bottom right panel). The top panels in this figure
correspond to the familiar spin-independent and spin-dependent interactions,
arising from the $\mathcal{O}_{1}$ and $\mathcal{O}_{4}$ operators,
respectively. The associated exclusion limits presented here agree with
previous results in the literature, see e.g. [1, 2]. Notice that the
($\lambda_{3},h_{3}$) case coincides with a dark photon model with Dirac DM.
The bottom panels correspond to models that generate the
$\overline{\mathcal{W}}_{1}(q,\Delta E)$ material response function only, but
with a DM response function scaling with $|\bm{q}|^{2}$ (left panel) and
$|\bm{q}|^{4}$ (right panel), respectively.

##  5 Summary and conclusion

A still open question in the field of sub-GeV DM direct detection is whether
experimental results can be interpreted within a theoretical framework where
the response of detector materials is described in terms of a single
ionisation or crystal form factors only, as in the case of the dark photon
model. In this work, we have addressed this question by computing the rate of
electronic transitions induced by DM-electron scattering in isolated xenon
atoms, as well as in silicon and germanium crystals, for a variety of models
where the DM particle can have spin 0, spin 1/2 or spin 1, while the particle
that mediates the DM-electron interaction can have either spin 0 or spin 1. We
have found several examples for which an accurate description of the non-
relativistic scattering of DM particles by the electrons bound in detector
materials requires material response functions that go beyond the standard
ionisation and crystal form factors. For simplicity, we have illustrated this
conclusion by restricting ourselves to the case of “simplified model”, where
just two coupling constants at the time are different from zero. However, the
fact that non-standard material response functions can arise in “minimal
extensions” of the dark photon model is a general conclusion of our study,
which does not rely on this restriction. This result corroborates our previous
findings [1, 2], and shows the importance of a response function formalism in
the interpretation of future DM direct detection data.

For the eight models that generate a non-standard material response, we have
calculated 90% C. L. exclusion limits on the reference cross section for DM-
electron electron scattering, Eq. (4.8), from the null result reported by the
XENON10 [43], XENON1T[15], EDELWEISS [23] and SENSEI [24] experiments. We have
performed this calculation as a function of the DM particle mass and for
different values of the mediator mass. We have not only considered small and
large mediator masses, corresponding to contact and long range interactions,
but also scenarios in which the mediator mass is comparable with the momentum
transfer in the scattering. Importantly, none of these exclusion limits could
have accurately been extracted from data without having first computed the
non-standard material response functions of Sec. 3.

A byproduct of our analysis is Eq. (3.28) which, expanding on our previous
work [1, 2], has enabled us to express the response of xenon and crystal
detectors to an external DM probe in the same, compact form (which was not the
case in [1, 2]). An interesting aspect of this equation is that it allows for
a simple comparison of the performance of different materials in the context
of sub-GeV DM direct detection. We have applied this result to compare the
response of xenon, silicon and germanium detectors finding that these
materials perform best in complementary regions of the plane spanned by the
momentum transfer and the deposited energy.

To conclude, we give a positive answer to the question of whether new material
response functions beyond the standard crystal and ionisation form factors are
needed for the interpretation of present and future direct detection
experiments searching for DM-induced electronic transitions.

## Acknowledgments

RC and TE acknowledge that this work was produced as part of the Knut and
Alice Wallenberg project Light Dark Matter (Dnr. KAW 2019.0080). RC also
acknowledges support from an individual research grant from the Swedish
Research Council (Dnr. 2018-05029). TE thanks the Theoretical Subatomic
Physics group at Chalmers University of Technology for its hospitality. NAS
and MM were supported by the ETH Zurich, and by the European Research Council
(ERC) under the European Union’s Horizon 2020 research and innovation
programme project HERO Grant Agreement No. 810451.

##  Appendix A Squared transition amplitudes

Here, we provide the details underling the derivation of Eqs. (3.46), (3.47),
(3.48) and (3.49). These four equations correspond to the DM response
functions for, respectively, fermionic dark matter with a scalar mediator,
fermionic dark matter with a vector mediator, vector dark matter with a scalar
mediator and, finally, vector dark matter with a vector mediator.

###  A.1 Fermionic dark matter with a scalar mediator

The free amplitude for DM-electron scattering in the case of spin-1/2 DM with
a scalar mediator is

$\displaystyle i\mathcal{M}$
$\displaystyle=\bar{u}^{s^{\prime}}_{\chi}(p^{\prime})i(\lambda_{1}+i\lambda_{2}\gamma^{5})u_{\chi}^{s}(p)\frac{i}{(p-p^{\prime})^{2}-m_{\phi}^{2}}\bar{u}_{e}^{r^{\prime}}(k^{\prime})i(h_{1}+ih_{2}\gamma^{5})u_{e}^{r}(k)\,,$
(A.1)

where $u_{e}^{r}$ and $u^{s}_{\chi}$ are electron and DM free spinors, while
$k$ and $k^{\prime}$ ($p$ and $p^{\prime}$) are the initial and final electron
(DM) four-momenta. By using the spinor bilinear expansions in App. B, we find
that in the non-relativistic limit the amplitude $i\mathcal{M}$ in Eq. (A.1)
can be written as follows

$\displaystyle i\mathcal{M}$ $\displaystyle\simeq
i4m_{\chi}m_{e}\left(\frac{\lambda_{1}h_{1}}{|\bm{q}\,|^{2}+m_{\phi}^{2}}\langle\mathcal{O}_{1}\rangle+\frac{\lambda_{2}h_{2}}{|\bm{q}\,|^{2}+m_{\phi}^{2}}\frac{m_{e}}{m_{\chi}}\langle\mathcal{O}_{6}\rangle-\frac{\lambda_{1}h_{2}}{|\bm{q}\,|^{2}+m_{\phi}^{2}}\langle\mathcal{O}_{10}\rangle\right.$
$\displaystyle+\left.\frac{\lambda_{2}h_{1}}{|\bm{q}\,|^{2}+m_{\phi}^{2}}\frac{m_{e}}{m_{\chi}}\langle\mathcal{O}_{11}\rangle\right)\,.$
(A.2)

The operators appearing in angle brackets in this expression are listed in
Tab. 1. From Eq. (A.2), we find

$\displaystyle\overline{|\mathcal{M}|^{2}}|f_{1\rightarrow 2}|^{2}$
$\displaystyle=\left(c_{1}^{2}+\frac{1}{16}c_{6}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{4}c_{10}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{1}{4}c_{11}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\right)|f_{1\rightarrow
2}|^{2}$ $\displaystyle
2m_{e}\mathrm{Re}\left\\{\overline{\mathcal{M}f_{1\rightarrow
2}[\nabla_{\bm{\ell}}\mathcal{M}^{*}]\cdot\bm{f}^{*}_{1\rightarrow
2}}\right\\}$ $\displaystyle=0$ $\displaystyle
m_{e}^{2}\overline{|\nabla_{\bm{\ell}}\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})]\cdot\bm{f}_{1\rightarrow 2}(\bm{q})|^{2}}$ $\displaystyle=0\,,$ (A.3)

which implies Eq. (3.46), once the azimuthal average in Eq. (3.33) is computed
by using Eq. (3.35). We give the relation between the non-relativistic
coupling constants in the equations above, i.e. $c_{1}$, $c_{6}$, $c_{10}$ and
$c_{11}$, and the particle masses and coupling constants in Eq. (A.2) in Tab.
3. Here and below, we implicitly assume that $\mathcal{M}$ bilinears are
evaluated at $\bm{\ell}=0$ and $\cos\theta=\xi$.

###  A.2 Fermionic dark matter with a vector mediator

In the case of spin-1/2 DM with a vector mediator, the free amplitude for DM-
electron scattering is

$\displaystyle i\mathcal{M}$
$\displaystyle=\left[\bar{u}^{s^{\prime}}(p^{\prime})_{\chi}i(\lambda_{3}\gamma_{\mu}+\lambda_{4}\gamma_{\mu}\gamma^{5})u_{\chi}^{s}(p)\right]\frac{-ig^{\mu\nu}}{(p-p^{\prime})^{2}-m_{G}^{2}}\left[\bar{u}_{e}^{r^{\prime}}(k^{\prime})i(h_{3}\gamma_{\nu}+h_{4}\gamma_{\nu}\gamma^{5})u_{e}^{r}(k)\right]\,,$
(A.4)

which, reduces to

$\displaystyle i\mathcal{M}$ $\displaystyle\simeq
i4m_{\chi}m_{e}\Bigg{[}-\frac{\lambda_{3}h_{3}}{|\bm{q}|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{1}\rangle+\frac{4\lambda_{4}h_{4}}{|\bm{q}|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{4}\rangle+\frac{2\lambda_{3}h_{4}}{|\bm{q}|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{7}\rangle-\frac{2\lambda_{4}h_{3}}{|\bm{q}|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{8}\rangle$
$\displaystyle+\left(\frac{2\lambda_{3}h_{4}}{|\bm{q}|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{\chi}}+\frac{2\lambda_{4}h_{3}}{|\bm{q}|^{2}+m_{G}^{2}}\right)\langle\hat{\mathcal{O}}_{9}\rangle\Bigg{]}\,,$
(A.5)

in the non-relativistic limit. We give explicit expressions for the operators
and coupling constants in this equation in Tab. 1 and Tab. 3, respectively.
Finally, by applying Eq. (A.5), we obtain

$\displaystyle\overline{|\mathcal{M}|^{2}}|f_{1\rightarrow 2}|^{2}$
$\displaystyle=\left(c_{1}^{2}+\frac{3}{16}c_{4}^{2}+\frac{1}{4}c_{7}^{2}|\bm{v}_{\rm
el}^{\perp}|^{2}+\frac{1}{4}c_{8}^{2}|\bm{v}^{\perp}_{\rm
el}|^{2}+\frac{1}{8}c_{9}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\right)|f_{1\rightarrow
2}|^{2}$ $\displaystyle
2m_{e}\mathrm{Re}\left\\{\overline{\mathcal{M}f_{1\rightarrow
2}[\nabla_{\bm{\ell}}\mathcal{M}^{*}]\cdot\bm{f}^{*}_{1\rightarrow
2}}\right\\}$
$\displaystyle=-\left(\frac{1}{2}c_{7}^{2}+\frac{1}{2}c_{8}^{2}\right)\bm{v}_{\rm
el}^{\perp}\cdot\mathrm{Re}\left(\bm{f}^{*}_{1\rightarrow 2}f_{1\rightarrow
2}\right)$ $\displaystyle
m_{e}^{2}\overline{|\nabla_{\bm{\ell}}\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})]\cdot\bm{f}_{1\rightarrow 2}(\bm{q})|^{2}}$
$\displaystyle=\left(\frac{1}{4}c_{7}^{2}+\frac{1}{4}c_{8}^{2}\right)|\bm{f}_{1\rightarrow
2}|^{2}\,.$ (A.6)

By taking the azimuthal average of the equations above, we obtain the DM
response functions in Eq. (3.47).

###  A.3 Vector dark matter with a scalar mediator

The free amplitude for DM-electron scattering in the case of spin-1 DM with a
scalar mediator is

$\displaystyle i\mathcal{M}$
$\displaystyle=\left[i(b_{1}m_{X}g^{\mu\nu})\varepsilon^{*s^{\prime}}_{\mu}(p^{\prime})\varepsilon^{s}_{\nu}(p)\right]\frac{i}{(p-p^{\prime})^{2}-m_{\phi}^{2}}\left[\bar{u}_{e}^{r^{\prime}}(k^{\prime})i(h_{1}+ih_{2}\gamma^{5})u_{e}^{r}(k)\right]\,,$
(A.7)

and, in the non-relativistic limit,

$\displaystyle i\mathcal{M}$ $\displaystyle\simeq
i4m_{X}m_{e}\left(\frac{1}{2}\frac{b_{1}h_{1}}{|\mathbf{q}|^{2}+m_{\phi}^{2}}\langle\hat{\mathcal{O}}_{1}\rangle-\frac{1}{2}\frac{b_{1}h_{2}}{|\mathbf{q}|^{2}+m_{\phi}^{2}}\langle\hat{\mathcal{O}}_{10}\rangle\right)\,,$
(A.8)

from which

$\displaystyle\overline{|\mathcal{M}|^{2}}|f_{1\rightarrow 2}|^{2}$
$\displaystyle=\left(c_{1}^{2}+\frac{1}{4}c_{10}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\right)|f_{1\rightarrow
2}|^{2}$ $\displaystyle
2m_{e}\mathrm{Re}\left\\{\overline{\mathcal{M}f_{1\rightarrow
2}[\nabla_{\bm{\ell}}\mathcal{M}^{*}]\cdot\bm{f}^{*}_{1\rightarrow
2}}\right\\}$ $\displaystyle=0$ $\displaystyle
m_{e}^{2}\overline{|\nabla_{\bm{\ell}}\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})]\cdot\bm{f}_{1\rightarrow 2}(\bm{q})|^{2}}$ $\displaystyle=0\,.$ (A.9)

The operators $\mathcal{O}_{1}$ and $\mathcal{O}_{10}$ are defined in Tab. 1,
whereas an explicit relation between $c_{1}$ and $c_{10}$ and the masses and
coupling constants in Eq. (A.8) can be found in Tab. 4. As a result, Eq.
(3.48) follows from Eq. (A.9) by a straightforward application of the
definition (3.33) and Eq. (3.35).

###  A.4 Vector dark matter with a vector mediator

Finally, in the case of spin-1 DM with a vector mediator the free amplitude
for DM-electron scattering is

$\displaystyle i\mathcal{M}=$
$\displaystyle\bigg{[}-ib_{5}\left(p^{\prime\mu}+p^{\mu}\right)\varepsilon^{s^{\prime}\nu*}(p\,^{\prime})\varepsilon^{s}_{\nu}(p)$
$\displaystyle+\mathrm{Re}(b_{6})\left(p_{\nu}^{\prime}-p_{\nu}\right)\left(\varepsilon^{s^{\prime}\nu*}(p\,^{\prime})\varepsilon^{s\mu}(p)+\varepsilon^{s^{\prime}\mu*}(p\,^{\prime})\varepsilon^{s\nu}(p)\right)$
$\displaystyle+i\mathrm{Im}(b_{6})\left(p_{\nu}^{\prime}-p_{\nu}\right)\left(\varepsilon^{s^{\prime}\nu*}(p\,^{\prime})\varepsilon^{s\mu}(p)-\varepsilon^{s^{\prime}\mu*}(p\,^{\prime})\varepsilon^{s\nu}(p)\right)$
$\displaystyle-\mathrm{Re}(b_{7})\varepsilon_{\mu\nu\rho\sigma}\left(p^{\prime\nu}+p^{\nu}\right)\varepsilon^{s^{\prime}\mu*}(p\,^{\prime})\varepsilon^{s\rho}(p)+i\mathrm{Im}(b_{7})\varepsilon_{\mu\nu\rho\sigma}\left(p^{\prime\nu}-p^{\nu}\right)\varepsilon^{s^{\prime}\mu*}(p\,^{\prime})\varepsilon^{s\rho}(p)\bigg{]}$
$\displaystyle\times\left[\frac{-1}{(p-p^{\prime})^{2}-m_{G}^{2}}\right]\overline{u}_{e}^{r^{\prime}}(k^{\prime})\gamma_{\mu}(h_{3}+h_{4}\gamma_{5})u_{e}^{r}(k)\,.$
(A.10)

By using the results listed in App. B, we find that in the non-relativistic
limit Eq. (A.10) reduces to

$\displaystyle i\mathcal{M}\simeq i4m_{X}m_{e}$
$\displaystyle\Bigg{[}-\frac{b_{5}h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{1}\rangle+\left(\frac{\mathrm{Im}(b_{6})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{2m_{X}}\frac{|\bm{q}\,|^{2}}{m_{e}^{2}}-\frac{2\mathrm{Re}(b_{7})h_{4}}{|\bm{q}\,|^{2}+m_{G}^{2}}\right)\langle\hat{\mathcal{O}}_{4}\rangle$
$\displaystyle+\frac{\mathrm{Im}(b_{6})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{2m_{X}}\langle\hat{\mathcal{O}}_{5}\rangle-\frac{\mathrm{Im}(b_{6})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{2m_{X}}\langle\hat{\mathcal{O}}_{6}\rangle+\frac{2b_{5}h_{4}}{|\bm{q}\,|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{7}\rangle$
$\displaystyle+\frac{\mathrm{Re}(b_{7})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\langle\hat{\mathcal{O}}_{8}\rangle+\left(\frac{\mathrm{Im}(b_{6})h_{4}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{X}}-\frac{\mathrm{Re}(b_{7})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\right)\langle\hat{\mathcal{O}}_{9}\rangle$
$\displaystyle+\frac{1}{2}\frac{\mathrm{Im}(b_{7})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{X}}\langle\hat{\mathcal{O}}_{11}\rangle-\frac{\mathrm{Im}(b_{7})h_{4}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{X}}\langle\hat{\mathcal{O}}_{14}\rangle$
$\displaystyle-\frac{\mathrm{Re}(b_{6})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{X}}\langle\hat{\mathcal{O}}_{17}\rangle+\frac{2\mathrm{Re}(b_{6})h_{4}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{X}}\langle\hat{\mathcal{O}}_{18}\rangle+\frac{\mathrm{Im}(b_{6})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}^{2}}{2m_{X}^{2}}\langle\hat{\mathcal{O}}_{19}\rangle$
$\displaystyle-\left(\frac{\mathrm{Re}(b_{6})h_{3}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}}{m_{X}}+\frac{\mathrm{Im}(b_{7})h_{4}}{|\bm{q}\,|^{2}+m_{G}^{2}}\frac{m_{e}^{2}}{m_{X}^{2}}\right)\langle\hat{\mathcal{O}}_{20}\rangle\Bigg{]}\,.$
(A.11)

As above, we define the operators in Eq. (A.11) in Tab. 1. From Eq. (A.11), we
obtain the expressions

$\displaystyle\overline{|\mathcal{M}|^{2}}|f_{1\rightarrow 2}|^{2}$
$\displaystyle=\Bigg{(}c_{1}^{2}+\frac{1}{2}c_{4}^{2}+\frac{2}{3}c_{5}^{2}\big{|}\frac{\bm{q}}{m_{e}}\times\bm{v}_{\rm
el}^{\perp}\big{|}^{2}+\frac{1}{6}c_{6}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{4}c_{7}^{2}|\bm{v}_{\rm
el}^{\perp}|^{2}$ $\displaystyle+\frac{2}{3}c_{8}^{2}|\bm{v}_{\rm
el}^{\perp}|^{2}+\frac{1}{3}c_{9}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{2}{3}c_{11}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{1}{6}c_{14}^{2}|\bm{v}_{\rm
el}^{\perp}|^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}$
$\displaystyle+\frac{1}{6}c_{17}^{2}\Big{(}\big{|}\frac{\bm{q}}{m_{e}}\cdot\bm{v}_{\rm
el}^{\perp}\big{|}^{2}+\frac{|\bm{q}|^{2}}{m_{e}^{2}}|\bm{v}_{\rm
el}^{\perp}|^{2}\Big{)}+\frac{1}{6}c_{18}^{2}\frac{|\bm{q}|^{2}}{m_{e}^{2}}$
$\displaystyle+\frac{1}{3}c_{19}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{12}c_{20}^{2}\frac{|\bm{q}|^{4}}{m_{e}^{4}}+\frac{1}{3}c_{4}c_{6}\frac{|\bm{q}|^{2}}{m_{e}^{2}}+\frac{2}{3}c_{1}c_{19}\frac{|\bm{q}|^{2}}{m_{e}^{2}}\Bigg{)}|f_{1\rightarrow
2}|^{2}$ $\displaystyle
2m_{e}\mathrm{Re}\left\\{\overline{\mathcal{M}f_{1\rightarrow
2}[\nabla_{\bm{\ell}}\mathcal{M}^{*}]\cdot\bm{f}^{*}_{1\rightarrow
2}}\right\\}$ $\displaystyle=\left[\frac{\bm{q}}{m_{e}}\cdot\bm{v}_{\rm
el}^{\perp}\left(\frac{4}{3}c_{5}^{2}-\frac{1}{3}c_{17}^{2}\right)\right]\frac{\bm{q}}{m_{e}}\cdot\mathrm{Re}\left(\bm{f}^{*}_{1\rightarrow
2}f_{1\rightarrow 2}\right)$
$\displaystyle-\left[\frac{1}{2}c_{7}^{2}+\frac{4}{3}c_{8}^{2}+\frac{|\bm{q}|^{2}}{m_{e}^{2}}\left(\frac{4}{3}c_{5}^{2}+\frac{1}{3}c_{14}^{2}+\frac{1}{3}c_{17}^{2}\right)\right]$
$\displaystyle\times\bm{v}_{\rm
el}^{\perp}\cdot\mathrm{Re}\left(\bm{f}^{*}_{1\rightarrow 2}f_{1\rightarrow
2}\right)$ $\displaystyle
m_{e}^{2}\overline{|\nabla_{\bm{\ell}}\mathcal{M}(\bm{q},\bm{v}^{\perp}_{\rm
el})]\cdot\bm{f}_{1\rightarrow 2}(\bm{q})|^{2}}$
$\displaystyle=\left[\frac{1}{4}c_{7}^{2}+\frac{2}{3}c_{8}^{2}+\frac{|\bm{q}|^{2}}{m_{e}^{2}}\left(\frac{2}{3}c_{5}^{2}+\frac{1}{6}c_{14}^{2}+\frac{1}{6}c_{17}^{2}\right)\right]|\bm{f}_{1\rightarrow
2}|^{2}$
$\displaystyle+\left(-\frac{2}{3}c_{5}^{2}+\frac{1}{6}c_{17}^{2}\right)\left|\frac{\bm{q}}{m_{e}}\cdot\bm{f}_{1\rightarrow
2}\right|^{2}\,,$ (A.12)

from which Eq. (3.49) follows by computing an azimuthal average, as shown
explicitly in Eq. (3.35). As in the previous examples, the coupling constants
in Eq. (A.12) are defined in Tab. 4.

##  Appendix B Useful identities

Below, we list expressions that we use in the non-relativistic reduction of
free scattering amplitudes, and in the calculation of the squared modulus of
the transition amplitude for the models in Sec. 2. We start from the non-
relativistic reduction of free spinors and spinor bilinears,

$\displaystyle
u_{\chi}^{r}(p)\simeq\frac{1}{\sqrt{4m_{\chi}}}\begin{pmatrix}(2m_{\chi}-\bm{p}\cdot\bm{\sigma})\xi^{r}\\\
(2m_{\chi}+\bm{p}\cdot\bm{\sigma})\xi^{r}\\\ \end{pmatrix}\,,$

and

$\displaystyle\bar{u}_{\chi}^{s^{\prime}}(p^{\prime})u_{\chi}^{s}(p)$
$\displaystyle\simeq 2m_{\chi}\,\xi^{\dagger s^{\prime}}\xi^{s}$ (B.1)
$\displaystyle\bar{u}_{\chi}^{s^{\prime}}(p^{\prime})i\gamma^{5}u_{\chi}^{s}(p)$
$\displaystyle\simeq 2i\bm{q}\cdot\bm{s}^{s^{\prime}s}$ (B.2)
$\displaystyle\bar{u}_{\chi}^{s^{\prime}}(p^{\prime})\gamma^{\mu}u_{\chi}^{s}(p)$
$\displaystyle\simeq\begin{pmatrix}2m_{\chi}\,\xi^{\dagger
s^{\prime}}\xi^{s}\\\ \bm{P}\,\xi^{\dagger
s^{\prime}}\xi^{s}+2i\bm{q}\times\bm{s}^{s^{\prime}s}\end{pmatrix}$ (B.3)
$\displaystyle\bar{u}_{\chi}^{s^{\prime}}(p^{\prime})\gamma^{\mu}\gamma^{5}u_{\chi}^{s}(p)$
$\displaystyle\simeq\begin{pmatrix}2\bm{P}\cdot\bm{s}^{s^{\prime}s}\\\
4m_{\chi}\,\bm{s}^{s^{\prime}s}\end{pmatrix}$ (B.4)
$\displaystyle\bar{u}_{\chi}^{s^{\prime}}(p^{\prime})\sigma^{\mu\nu}u_{\chi}^{s}(p)$
$\displaystyle\simeq\begin{pmatrix}0&i\bm{q}\,\xi^{\dagger
s^{\prime}}\xi^{s}-2\bm{P}\times\bm{s}^{s^{\prime}s}\\\ -i\bm{q}\,\xi^{\dagger
s^{\prime}}\xi^{s}+2\bm{P}\times\bm{s}^{s^{\prime}s}&4m_{\chi}\,\varepsilon_{ijk}\,s^{s^{\prime}sk}\end{pmatrix}\,,$
(B.5)

respectively, which apply to the case of fermionic DM. Here,
$\bm{s}^{s^{\prime}s}\equiv\xi^{\dagger s^{\prime}}\bm{S}_{\rm DM}\xi^{s}$,
$\bm{P}=\bm{p}+\bm{p}^{\prime}$, $\bm{q}=\bm{p}-\bm{p}^{\prime}$, $\xi^{s}$,
$s=1,2$ and $\xi^{\dagger s^{\prime}}$, $s^{\prime}=1,2$ are two component
spinors, while $\bm{\sigma}=(\sigma_{1},\sigma_{2},\sigma_{3})$ is a vector
the component of which are the three Pauli matrices. Similar expressions apply
to the free electron spinor $u_{e}^{r}$. In the non-relativistic reduction of
scattering amplitudes, we also use the polarisation vectors of massive vector
DM or spin-1 mediator particles. They can be expanded as

$\displaystyle\varepsilon^{s\mu}_{\mu}(\bm{p})\simeq\begin{pmatrix}\frac{1}{2m_{X}}(\bm{P}+\bm{q}\,)\cdot\bm{e}_{s}\\\
\bm{e}_{s}\end{pmatrix}\quad\quad\varepsilon^{s^{\prime}\mu*}(\bm{p}^{\prime})\simeq\begin{pmatrix}\frac{1}{2m_{X}}(\bm{P}-\bm{q}\,)\cdot\bm{e}^{\prime}_{s^{\prime}}\\\
\bm{e}^{\prime}_{s^{\prime}}\end{pmatrix}\,,$ (B.6)

where $s,s^{\prime}=1,2,3$, $m_{X}$ can be either the DM or the mediator mass,
$\bm{e}_{3}=\bm{p}/|\bm{p}|$,
$\bm{e}^{\prime}_{3^{\prime}}=\bm{p}^{\prime}/|\bm{p}^{\prime}|$, and
$\bm{e}_{s}\cdot\bm{e}^{\prime}_{s^{\prime}}\simeq\delta_{ss^{\prime}}$. In
addition, we use the non-relativistic expansion of the momentum transfer four-
vector,

$\displaystyle
q\equiv(E_{\bm{p}}-E_{\bm{p}^{\prime}},\bm{p}-\bm{p}^{\prime})\simeq\left(\frac{\bm{P}\cdot\bm{q}}{2m_{\rm
DM}}\,,\,\bm{q}\right)\,,$ (B.7)

which allows us to write the non-relativistic reduction of scattering
amplitudes in terms of $\bm{q}$ and

$\displaystyle\bm{v}_{\rm el}^{\perp}$
$\displaystyle=\frac{\bm{p}+\bm{p}\,^{\prime}}{2m_{\rm
DM}}-\frac{\bm{k}+\bm{k}\,^{\prime}}{2m_{e}}=\bm{v}-\frac{\bm{k}}{m_{e}}-\frac{\bm{q}}{2\mu_{\rm
DMe}}\,,$ (B.8)

where $\bm{k}$ ($\bm{k}^{\prime}$) is the initial (final) free electron
momentum, and $\bm{v}$ the DM-electron relative velocity. We conclude this
appendix by listing useful spin sums for electrons or spin-1/2 DM particles,

$\displaystyle\frac{1}{2}\sum_{r^{\prime}r}|\bm{s}^{r^{\prime}r}|^{2}=\frac{3}{4}$
$\displaystyle\frac{1}{2}\sum_{r^{\prime}r}(\bm{s}^{r^{\prime}r}\times\bm{A})\cdot\bm{s}^{r^{\prime}r}=0$
(B.9)
$\displaystyle\frac{1}{2}\sum_{r^{\prime}r}(\bm{s}^{r^{\prime}r}\cdot\bm{A})(\bm{s}^{r^{\prime}r}\cdot\bm{B})=\frac{1}{4}\bm{A}\cdot\bm{B}$
$\displaystyle\frac{1}{2}\sum_{r^{\prime}r}(\bm{s}^{r^{\prime}r}\cdot\bm{A})\delta^{r^{\prime}r}=0$
(B.10)
$\displaystyle\frac{1}{2}\sum_{r^{\prime}r}(\bm{s}^{r^{\prime}r}\times\bm{A})(\bm{s}^{r^{\prime}r}\times\bm{B})=\frac{1}{2}\bm{A}\cdot\bm{B}$
$\displaystyle\frac{1}{2}\sum_{r^{\prime}r}(\bm{s}^{r^{\prime}r}\times\bm{A})\delta^{r^{\prime}r}=0$
(B.11)

and for spin-1 DM particles,

$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}\bm{A}\cdot\bm{s}_{V}^{s^{\prime}s}=0$
(B.12)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}\bm{A}\cdot\mathcal{S}^{s^{\prime}s}\cdot\bm{s}_{V}^{s^{\prime}s}=0$
(B.13)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}(\bm{A}\cdot\bm{s}_{V}^{s^{\prime}s})(\bm{B}\cdot\bm{s}_{V}^{s^{\prime}s})=\frac{2}{3}\bm{A}\cdot\bm{B}$
(B.14)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}(\bm{A}\times\bm{s}_{V}^{s^{\prime}s})\cdot(\bm{B}\times\bm{s}_{V}^{s^{\prime}s})=\frac{4}{3}\bm{A}\cdot\bm{B}$
(B.15)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}(\bm{A}\cdot\mathcal{S}^{s^{\prime}s}\cdot\bm{B})\delta^{s^{\prime}s}=\frac{1}{3}\bm{A}\cdot\bm{B}$
(B.16)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}(\bm{A}\cdot\mathcal{S}^{s^{\prime}s}\cdot\bm{B})(\bm{s}_{V}^{s^{\prime}s}\cdot\bm{C})=0$
(B.17)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}(\bm{A}\cdot\mathcal{S}^{s^{\prime}s})\cdot(\bm{B}\cdot\mathcal{S}^{s^{\prime}s})=\frac{2}{3}\bm{A}\cdot\bm{B}$
(B.18)
$\displaystyle\frac{1}{3}\sum_{s^{\prime}s}|\bm{A}\cdot\mathcal{S}^{s^{\prime}s}\cdot\bm{B}|^{2}=\frac{1}{6}|\bm{A}\cdot\bm{B}|^{2}+\frac{1}{6}|\bm{A}|^{2}|\bm{B}|^{2}\,,$
(B.19)

where $\bm{A}$ and $\bm{B}$ are arbitrary three-dimensional vectors, $r=1,2$,
$r^{\prime}=1,2$, $s=1,2,3$, $s^{\prime}=1,2,3$, and

$\displaystyle\bm{s}_{V}^{\ s^{\prime}s}$
$\displaystyle\equiv-i\bm{e}^{\prime}_{s^{\prime}}\times\bm{e}_{s}$ (B.20)
$\displaystyle\mathcal{S}^{s^{\prime}s}_{ij}$
$\displaystyle\equiv\frac{1}{2}(e_{si}e^{\prime}_{s^{\prime}j}+e_{sj}e^{\prime}_{s^{\prime}i})\,.$
(B.21)

The expressions above are useful in the evaluation of the squared modulus of
transition amplitudes.

##  Appendix C One-loop amplitudes

Below, we list the free electron amplitudes for the four diagrams in Fig. 1 in
the non-relativistic limit.

1). Scalar DM with a scalar mediator. For the top left diagram in Fig. 1, we
obtain the amplitude

$\displaystyle i\mathcal{M}\simeq-i4m_{S}m_{e}$
$\displaystyle\bigg{[}\frac{g_{2}h_{1}^{2}}{2}\frac{m_{e}}{m_{S}}\frac{\delta_{\lceil
h_{2}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y+1}{\Delta}$
$\displaystyle+\frac{g_{2}h_{2}^{2}}{2}\frac{m_{e}}{m_{S}}\frac{\delta_{\lceil
h_{1}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y-1}{\Delta}\bigg{]}\langle\hat{\mathcal{O}}_{1}\rangle\,,$
(C.1)

where we assume that the only coupling constants different from zero are
either $g_{2}$ and $h_{1}$, or $g_{2}$ and $h_{1}$, while

$\displaystyle\Delta=(x^{2}+y^{2})m_{e}^{2}+(x+y)(m_{M}^{2}-m_{e}^{2})+xy|\bm{q}\,|^{2}+(1-x-y+2xy)m_{e}^{2}\,,$
(C.2)

with $m_{M}=m_{\phi}$.

2) Scalar DM with a vector mediator. The amplitude for the top right diagram
in Fig. 1 is equal to

$\displaystyle i\mathcal{M}\simeq i4m_{S}m_{e}$
$\displaystyle\bigg{[}g_{3}h_{3}^{2}\frac{m_{e}}{m_{S}}\frac{\delta_{\lceil
h_{4}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y-2}{\Delta}$
$\displaystyle+g_{3}h_{4}^{2}\frac{m_{e}}{m_{S}}\frac{\delta_{\lceil
h_{3}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y+2}{\Delta}\bigg{]}\langle\hat{\mathcal{O}}_{1}\rangle$
(C.3)

where $\Delta$ is given in Eq. (C.2) with $m_{M}=m_{G}$. In this latter
equation, we assume that the coupling constants $g_{3}$ and $h_{3}$, or
$g_{3}$ and $h_{4}$, are the only ones being different from zero in the
Lagrangian.
3) Vector DM with a scalar mediator. For the bottom left diagram in Fig. 1, we
find the amplitude

$\displaystyle i\mathcal{M}\simeq-i4m_{X}m_{e}$
$\displaystyle\bigg{[}\frac{b_{2}h_{1}^{2}}{2}\frac{m_{e}}{m_{X}}\frac{\delta_{\lceil
h_{2}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y+1}{\Delta}$
$\displaystyle+\frac{b_{2}h_{2}^{2}}{2}\frac{m_{e}}{m_{X}}\frac{\delta_{\lceil
h_{1}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y-1}{\Delta}\bigg{]}\langle\hat{\mathcal{O}}_{1}\rangle\,,$
(C.4)

where $m_{M}=m_{\phi}$ in $\Delta$, and either $b_{2}$ and $h_{1}$, or $b_{2}$
and $h_{1}$ are the only non-zero coupling constants.

4) Vector DM with a vector mediator. Finally, the amplitude for the bottom
right diagram in Fig. 1 is

$\displaystyle i\mathcal{M}\equiv i4m_{X}m_{e}$
$\displaystyle\bigg{[}b_{3}h_{3}^{2}\frac{m_{e}}{m_{X}}\frac{\delta_{\lceil
h_{4}\rceil 0}\delta_{\lceil b_{4}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y-2}{\Delta}$
$\displaystyle+b_{3}h_{4}^{2}\frac{m_{e}}{m_{X}}\frac{\delta_{\lceil
h_{3}\rceil 0}\delta_{\lceil b_{4}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y+2}{\Delta}$
$\displaystyle+\frac{b_{4}h_{3}^{2}}{2}\frac{m_{e}}{m_{X}}\frac{\delta_{\lceil
h_{4}\rceil 0}\delta_{\lceil b_{3}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y-1}{\Delta}$
$\displaystyle+\frac{b_{4}h_{4}^{2}}{2}\frac{m_{e}}{m_{X}}\frac{\delta_{\lceil
h_{3}\rceil 0}\delta_{\lceil b_{3}\rceil
0}}{(4\pi)^{2}}\int_{0}^{1}\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y+1}{\Delta}\bigg{]}\langle\hat{\mathcal{O}}_{1}\rangle\,,$
(C.5)

where only $b_{3}$ and $h_{3}$, $b_{3}$ and $h_{4}$, $b_{4}$ and $h_{3}$, or
$b_{4}$ and $h_{4}$ are simultaneously different from zero, and $m_{M}=m_{G}$.
Inspection of the above equations shows that all non-relativistic amplitudes
listed here can be written in terms of a single integral function of the
integer number $a$, namely

$\displaystyle
I(a)=\int_{0}^{1}\\!\\!\\!\\!\,\mathrm{d}x\int_{0}^{1-x}\\!\\!\\!\\!\\!\\!\\!\\!\,\mathrm{d}y\frac{x+y+a}{\Delta}\,.$
(C.6)

## References

  * [1] R. Catena, T. Emken, N. A. Spaldin and W. Tarantino, _Atomic responses to general dark matter-electron interactions_ , _Phys. Rev. Res._ 2 (2020) 033195 [1912.08204].
  * [2] R. Catena, T. Emken, M. Matas, N. A. Spaldin and E. Urdshals, _Crystal responses to general dark matter-electron interactions_ , _Phys. Rev. Res._ 3 (2021) 033149 [2105.02233].
  * [3] “Astroparticle Physics European Consortium (APPEC).” http://www.appec.org/roadmap.
  * [4] A. Drukier and L. Stodolsky, _Principles and Applications of a Neutral Current Detector for Neutrino Physics and Astronomy_ , _Phys. Rev._ D30 (1984) 2295.
  * [5] M. W. Goodman and E. Witten, _Detectability of Certain Dark Matter Candidates_ , _Phys. Rev._ D31 (1985) 3059.
  * [6] T. Marrodán Undagoitia and L. Rauch, _Dark matter direct-detection experiments_ , _J. Phys._ G43 (2016) 013001 [1509.08767].
  * [7] G. Arcadi, M. Dutra, P. Ghosh, M. Lindner, Y. Mambrini, M. Pierre et al., _The waning of the WIMP? A review of models, searches, and constraints_ , _Eur. Phys. J._ C78 (2018) 203 [1703.07364].
  * [8] L. Roszkowski, E. M. Sessolo and S. Trojanowski, _WIMP dark matter candidates and searches - current status and future prospects_ , _Rept. Prog. Phys._ 81 (2018) 066201 [1707.06277].
  * [9] A. Mitridate, T. Trickle, Z. Zhang and K. M. Zurek, _Snowmass White Paper: Light Dark Matter Direct Detection at the Interface With Condensed Matter Physics_ , in _2022 Snowmass Summer Study_ , 3, 2022, 2203.07492.
  * [10] M. Battaglieri et al., _US Cosmic Visions: New Ideas in Dark Matter 2017; Community Report_ , _FERMILAB-CONF-17-282-AE-PPD-T_ (2017) [1707.04591].
  * [11] R. Essig, J. Mardon and T. Volansky, _Direct Detection of Sub-GeV Dark Matter_ , _Phys. Rev._ D85 (2012) 076007 [1108.5383].
  * [12] Y. Kahn and T. Lin, _Searches for light dark matter using condensed matter systems_ , _Rept. Prog. Phys._ 85 (2022) 066901 [2108.03239].
  * [13] R. Essig, T. Volansky and T.-T. Yu, _New Constraints and Prospects for sub-GeV Dark Matter Scattering off Electrons in Xenon_ , _Phys. Rev._ D96 (2017) 043017 [1703.00910].
  * [14] DarkSide collaboration, P. Agnes et al., _Constraints on Sub-GeV Dark-Matter–Electron Scattering from the DarkSide-50 Experiment_ , _Phys. Rev. Lett._ 121 (2018) 111303 [1802.06998].
  * [15] XENON collaboration, E. Aprile et al., _Light Dark Matter Search with Ionization Signals in XENON1T_ , _Phys. Rev. Lett._ 123 (2019) 251801 [1907.11485].
  * [16] XENON collaboration, E. Aprile et al., _Excess electronic recoil events in XENON1T_ , _Phys. Rev. D_ 102 (2020) 072004 [2006.09721].
  * [17] P. W. Graham, D. E. Kaplan, S. Rajendran and M. T. Walters, _Semiconductor Probes of Light Dark Matter_ , _Phys. Dark Univ._ 1 (2012) 32 [1203.2531].
  * [18] R. Essig, M. Fernandez-Serra, J. Mardon, A. Soto, T. Volansky and T.-T. Yu, _Direct Detection of sub-GeV Dark Matter with Semiconductor Targets_ , _JHEP_ 05 (2016) 046 [1509.01598].
  * [19] S. Derenzo, R. Essig, A. Massari, A. Soto and T.-T. Yu, _Direct Detection of sub-GeV Dark Matter with Scintillating Targets_ , _Phys. Rev._ D96 (2017) 016026 [1607.01009].
  * [20] SuperCDMS collaboration, R. Agnese et al., _First Dark Matter Constraints from a SuperCDMS Single-Charge Sensitive Detector_ , _Phys. Rev. Lett._ 121 (2018) 051301 [1804.10697].
  * [21] N. A. Kurinsky, T. C. Yu, Y. Hochberg and B. Cabrera, _Diamond Detectors for Direct Detection of Sub-GeV Dark Matter_ , _Phys. Rev. D_ 99 (2019) 123005 [1901.07569].
  * [22] DAMIC collaboration, A. Aguilar-Arevalo et al., _Constraints on Light Dark Matter Particles Interacting with Electrons from DAMIC at SNOLAB_ , _Phys. Rev. Lett._ 123 (2019) 181802 [1907.12628].
  * [23] EDELWEISS collaboration, Q. Arnaud et al., _First germanium-based constraints on sub-MeV Dark Matter with the EDELWEISS experiment_ , _Phys. Rev. Lett._ 125 (2020) 141301 [2003.01046].
  * [24] SENSEI collaboration, L. Barak et al., _SENSEI: Direct-Detection Results on sub-GeV Dark Matter from a New Skipper-CCD_ , _Phys. Rev. Lett._ 125 (2020) 171802 [2004.11378].
  * [25] S. M. Griffin, Y. Hochberg, K. Inzani, N. Kurinsky, T. Lin and T. Chin, _Silicon carbide detectors for sub-GeV dark matter_ , _Phys. Rev. D_ 103 (2021) 075002 [2008.08560].
  * [26] S. M. Griffin, K. Inzani, T. Trickle, Z. Zhang and K. M. Zurek, _Extended calculation of dark matter-electron scattering in crystal targets_ , _Phys. Rev. D_ 104 (2021) 095015 [2105.05253].
  * [27] S. Knapen, J. Kozaczuk and T. Lin, _Dark matter-electron scattering in dielectrics_ , _Phys. Rev. D_ 104 (2021) 015031 [2101.08275].
  * [28] Y. Hochberg, Y. Kahn, N. Kurinsky, B. V. Lehmann, T. C. Yu and K. K. Berggren, _Determining Dark-Matter–Electron Scattering Rates from the Dielectric Function_ , _Phys. Rev. Lett._ 127 (2021) 151802 [2101.08263].
  * [29] R. Lasenby and A. Prabhu, _Dark matter–electron scattering in materials: Sum rules and heterostructures_ , _Phys. Rev. D_ 105 (2022) 095009 [2110.01587].
  * [30] H.-Y. Chen, A. Mitridate, T. Trickle, Z. Zhang, M. Bernardi and K. M. Zurek, _Dark matter direct detection in materials with spin-orbit coupling_ , _Phys. Rev. D_ 106 (2022) 015024 [2202.11716].
  * [31] Y. Hochberg, Y. Zhao and K. M. Zurek, _Superconducting Detectors for Superlight Dark Matter_ , _Phys. Rev. Lett._ 116 (2016) 011301 [1504.07237].
  * [32] Y. Hochberg, B. V. Lehmann, I. Charaev, J. Chiles, M. Colangelo, S. W. Nam et al., _New Constraints on Dark Matter from Superconducting Nanowires_ , 2110.01586.
  * [33] Y. Hochberg, Y. Kahn, M. Lisanti, K. M. Zurek, A. G. Grushin, R. Ilan et al., _Detection of sub-MeV Dark Matter with Three-Dimensional Dirac Materials_ , _Phys. Rev._ D97 (2018) 015004 [1708.08929].
  * [34] R. M. Geilhufe, F. Kahlhoefer and M. W. Winkler, _Dirac Materials for Sub-MeV Dark Matter Detection: New Targets and Improved Formalism_ , _Phys. Rev. D_ 101 (2020) 055005 [1910.02091].
  * [35] A. Coskuner, A. Mitridate, A. Olivares and K. M. Zurek, _Directional dark matter detection in anisotropic Dirac materials_ , _Phys. Rev. D_ 103 (2021) 016006.
  * [36] Y. Hochberg, Y. Kahn, M. Lisanti, C. G. Tully and K. M. Zurek, _Directional detection of dark matter with two-dimensional targets_ , _Phys. Lett._ B772 (2017) 239 [1606.08849].
  * [37] G. Cavoto, M. G. Betti, C. Mariani, F. Pandolfi, A. D. Polosa, I. Rago et al., _Carbon nanotubes as anisotropic target for dark matter_ , _J. Phys. Conf. Ser._ 1468 (2020) 012232 [1911.01122].
  * [38] S. Knapen, T. Lin, M. Pyle and K. M. Zurek, _Detection of Light Dark Matter With Optical Phonons in Polar Materials_ , _Phys. Lett._ B785 (2018) 386 [1712.06598].
  * [39] T. Trickle, Z. Zhang, K. M. Zurek, K. Inzani and S. M. Griffin, _Multi-Channel Direct Detection of Light Dark Matter: Theoretical Framework_ , _JHEP_ 03 (2020) 036 [1910.08092].
  * [40] T. Trickle, Z. Zhang and K. M. Zurek, _Detecting Light Dark Matter with Magnons_ , _Phys. Rev. Lett._ 124 (2020) 201801 [1905.13744].
  * [41] B. Holdom, _Two U(1)’s and Epsilon Charge Shifts_ , _Phys. Lett._ 166B (1986) 196.
  * [42] (XENON Collaboration)§, XENON collaboration, E. Aprile et al., _Emission of single and few electrons in XENON1T and limits on light dark matter_ , _Phys. Rev. D_ 106 (2022) 022001 [2112.12116].
  * [43] XENON10 collaboration, J. Angle et al., _A search for light dark matter in XENON10 data_ , _Phys. Rev. Lett._ 107 (2011) 051301 [1104.3088].
  * [44] R. Catena, J. Conrad and M. B. Krauss, _Compatibility of a dark matter discovery at XENONnT or LZ with the WIMP thermal production mechanism_ , _Phys. Rev._ D97 (2018) 103002 [1712.07969].
  * [45] B. A. Dobrescu and I. Mocioiu, _Spin-dependent macroscopic forces from new particle exchange_ , _JHEP_ 11 (2006) 005 [hep-ph/0605342].
  * [46] J. Fan, M. Reece and L.-T. Wang, _Non-relativistic effective theory of dark matter direct detection_ , _JCAP_ 1011 (2010) 042 [1008.1591].
  * [47] A. L. Fitzpatrick, W. Haxton, E. Katz, N. Lubbers and Y. Xu, _The Effective Field Theory of Dark Matter Direct Detection_ , _JCAP_ 1302 (2013) 004 [1203.3542].
  * [48] R. Catena and P. Ullio, _A novel determination of the local dark matter density_ , _JCAP_ 1008 (2010) 004 [0907.0018].
  * [49] D. Baxter et al., _Recommended conventions for reporting results from direct dark matter searches_ , _Eur. Phys. J. C_ 81 (2021) 907 [2105.00599].
  * [50] T. Emken, “Dark Atomic Response Tabulator (DarkART)[Code, v0.1.0].” The code can be found under https://github.com/temken/darkart., 2021.  DOI:10.5281/zenodo.6046225.
  * [51] E. Urdshals and M. Matas, _QEdark-EFT_ , _Zenodo_ (2021) .