Pushing the frontier of WIMPy inelastic dark matter: journey to the end of the periodic table

We explore the reach of low-background experiments made of small quantities of heavy nuclear isotopes in probing the parameter space of inelastic dark matter that is kinematically inaccessible to classic direct detection experiments. Through inelastic scattering with target nuclei, dark matter can yield a signal either via nuclear recoil or nuclear excitation. We present new results based on this approach, using data from low-energy gamma quanta searches in low-background experiments with Hf and Os metal samples, and measurements with CaWO$_4$ and PbWO$_4$ crystals as scintillating bolometers. We place novel bounds on WIMPy inelastic dark matter up to mass splittings of about 640 keV, and provide forecasts for the reach of future experiments.


I. INTRODUCTION
Dark matter (DM) direct detection experiments are designed to record the small amounts of energy expected to be deposited by galactic DM particles as they pass through the detector and interact with nuclei or electrons in the target material. Most such detectors are optimised to search for WIMP-like particles, with masses in the GeV to TeV range, and weak-scale DM-nucleon elastic scattering cross sections σ n 10 −40 cm 2 . These typically employ target nuclei such as C, O, Si, Ar, Ge or Xe, which are kinematically well-matched to the 10-100 GeV range. Experiments are typically conducted deep underground to mitigate cosmic ray and radiogenic background components.
A simple inelastic dark matter model consists of a lighter state χ 1 and a heavier state χ 2 with a mass splitting δ ≡ M χ2 − M χ1 . Generically, the lighter state will dominate the relic abundance, and two-to-two scattering processes with baryons must lead to an excitation or deexcitation of the DM state. In contrast to elastic scattering, the kinematics of inelastic dark matter exhibits two features. First, the dark matter kinetic energy must be large enough to overcome the mass splitting, imposing a lower limit on the relative velocity required for an interaction to take place. Second, the minimum momentum exchange required for a collision to occur increases with the mass splitting, meaning that IDM with typical halo velocities can only scatter with heavier nuclei and deposit large enough recoil energy. In other words, even with sufficient kinetic energy, there is a momentum transfer threshold ∼ µ χN v below which interactions may not occur. As a consequence, inelastic dark matter may evade the searches in most direct detection experiments, either because the analysis region is limited to low ( 10 − 50 keV) recoil energies, or because of the limited target nucleus mass.
We thus turn to heavier target nuclei (Z > 54), which have received little attention in the search for DMinduced recoils. Many elements between Z ∼ 58 (cerium) and Z ∼ 83 (bismuth) are sensitive to IDM models with 100 keV mass splittings, while remaining sufficiently stable not to overwhelm a detector with alpha or beta decay backgrounds. Such elements, in particular tungsten and lead, have been utilized in scintillating crystals to search for rare alpha decays [60,61]. These experiments feature both light and heat measurements, which pro-FIG. 1. The two processes considered in this work. Nuclear Recoil (NR): dark matter is up-scattered to a heavier state and the target nucleus recoils from the impact, leading to heat, scintillation, and/or ionization signals. Nuclear Excitation (NE): in addition to NR, the nucleus can also be excited and then deexcite to ground state with the emission of one or several gamma quanta. vide a powerful tool to distinguish DM-induced nuclear recoil (NR) events due to high momentum-transfer coherent DM-nucleus interactions (Figure 1 top) from background events.
In addition to nuclear recoil signals from coherent scattering, heavier elements can contribute a second type of signal, thanks to their comparatively loosely bound outer nucleon shells: gamma photons from induced nuclear excitation (NE, Figure 1 bottom) and subsequent decay, or deexcitation of a long-lived excited nuclear state. The goal of this paper is thus to use these two approaches, i.e. nuclear recoil searches and nuclear transition gammas to place novel constraints on inelastic DM models.
Recently, a search for collisional deexcitation of the metastable nuclear isomer 180m Ta induced by DM has led to novel constraints on IDM [62,63]. 180m Ta is fairly unique: it is stable on time scales much larger than the age of the Universe, and its natural abundance is well established. As the isomers of other nuclei are shortlived, this method is not generally applicable. However, the excitation of ground state nuclei to excited states remains possible in the scattering χ 1 + N → χ 2 + N * if the excitation energy is low enough [64][65][66]. We focus on a few isotopes, with possible nuclear transitions that we will explore listed in Table I below. Once scattered into an excited state, the gamma quanta emitted during subsequent decay to the ground state can be measured. By measuring the activity of these transition lines, we may obtain a conservative upper bound on the total inelastic scattering rate from DM interactions, and thus on the DM-nucleon scattering cross section. In comparison with traditional direct detection experiments, this method appears rather insensitive. Indeed, nuclear excitation is more kinematically suppressed, and the scattering does not enjoy the usual A 2 enhancement factor because the DM is interacting with a single nucleon rather than scattering coherently with the full nucleus. However, we identify two advantages: 1) The high mass of transition nuclei better matches heavier DM candidates, and 2) the energy of the outgoing gamma is independent of the recoil energy, since it is fixed by the structure of the target nucleus. This means that, for a given transition, there is no danger of the signal being lost outside of the analysis region of interest. Previously, we have performed such an analysis in the context of a search for rare decays of Hf [67] and set stronger bounds on IDM mass splittings than previously-reported results from direct detection experiments. In this paper, we extend our approach, using previously-published data collected within low-background measurements with Hf and Os metal samples, along with the limits derived from bolometeric measurements using CaWO 4 and PbWO 4 scintillating crystals. Although the NE bounds that we find are not as strong as the NR bounds that we derive from bolometer measurements, we will identify a path to leading constraints on WIMPy inelastic dark matter using the nuclear properties of target materials. This paper is organized as follows. In Sec. II we present the kinematics of inelastic dark matter scattering, including dark matter-induced nuclear transitions and their nuclear response. A number of experiments are listed in Sec. III and we show the resulting bounds on inelastic dark matter in Sec. IV. Finally, we conclude in Sec. V.

II. KINEMATICS OF INELASTIC DARK MATTER SCATTERING
In models of inelastic DM, the dark sector comprises two states χ 1 and χ 2 with a mass splitting δ = M χ2 − M χ1 . At tree level, χ 2 couples to χ 1 via a scalar or a gauge boson field, but the self-couplingχ 1 χ 1 orχ 2 χ 2 is loop-suppressed. In the simplest case, we consider the Lagrangian L ⊃ gφχ 2 χ 1 + h.c.. Since χ 2 is heavier than χ 1 , after freeze-out we expect χ 2 to be depleted either though decay or through the interconversion process χ 2 χ 2 → χ 1 χ 1 when the temperature of the universe falls below the mass splitting, T < δ [44,73]. As a consequence, the dark matter today is dominated by χ 1 .
We consider two types of dark matter interactions in the detector. The first is nuclear recoil (NR) where the target nucleus displaces, and χ 1 transforms to χ 2 . In addition, χ 1 may also interact with the valence nucleons and instigate nuclear excitation (NE) into an unstable state N * . As long as the mediator mass M φ is much larger than the momentum transfer q, the interaction can be described by the effective Lagrangian where the UV scale Λ is set by the mediator mass and the coupling between the dark sector and the SM. We define the quantity ∆, which represents the sum of the DM mass splitting δ and the excitation energy ∆E ≡ E N * − E N : ∆E vanishes in NR-only scattering. The kinematics of both interactions follows where the momentum transfer q = q f − q i . Eq. (3) can be written where is the reduced mass of the DM-nucleus system, and θ is the laboratory frame scattering angle. We also define M χ ≡ M χ1 M χ2 which is valid so long as M χ δ. We will assume the dark matter velocity v follows a Maxwellian distribution with v 0 = 220 km/s truncated at the escape velocity v esc = 600 km/s in the Earth's frame, assuming the Earth's velocity is v e = 240 km/s [74]. We refer the reader to Refs. [74][75][76][77] for detailed discussions of uncertainties on the halo parameter values. The minimum/maximum momentum transfer is obtained at | cos θ| = 1 which gives and the minimum/maximum recoil energy is simply obtained via the relationship E R = q 2 /(2M N ). The minimum dark matter velocity required to scatter off a nucleus with any recoil energy is thus Larger dark matter kinetic energy is therefore required for inelastic scattering as ∆ increases. Consequently, the maximum mass splitting that can be probed in an experiment is which grows for larger reduced mass µ χN and hence heavier target nuclei. For a specific nuclear recoil energy, the minimum dark matter velocity for an interaction to occur can also be obtained from Eq. (4): The kinematics of NR scattering is shown in Fig. 2 for a few nuclei. For each isotope, the region to the left of the corresponding solid curve represents allowed values of E R , given the IDM mass splitting δ. The shape of these curves is determined by the DM velocity distribution. Focusing first on xenon experiments (solid green line), the limitation is twofold: 1) the Xe atomic mass limits the maximum transfer as in Eq. (7). 2) The analyzed recoil energy, shown here for XENON1T as a dashed blue line, is restricted in most xenon experiments. For the recoil energy range of 1 − 40 keV the maximum mass splitting that can be constrained is thus around 250 keV. Recently, XENON1T reported the event spectrum up to the electron recoil energies of 210 keV [43]. A dedicated analysis of nuclear recoils at high energy would be required to constrain mass splitting above 200 keV, but the sensitivity will still be limited to δ < 450 keV due to the limited mass number of Xe isotopes. Conversely, tungsten (solid grey line) is a heavy and stable element, Kinematics of inleastic dark matter with mass splitting δ scattering with a target nucleus with/without nuclear transition. Kinematically allowed recoil energies are shown for 132 Xe, 184 W, 208 Pb and 189 Os with the solid lines respectively. The region for 208 Pb is further shaded. We consider nuclear recoil (NR) for the first three target nuclei, and nuclear excitation (NE) to the 36.17 keV excited state for 189 Os. The maximum recoil energy analyzed in XENON1T 2018 is about 40 keV and that of CRESST-II [79] is about 120 keV, as depicted by the horizontal dashed lines. These limit the experiments' respective sensitivities to higher mass splittings. In this figure, the dark matter mass is fixed at well-suited to exploring high mass splittings. CRESST-II (CaWO 4 , dashed grey line) analyzed recoil data up to 120 keV, which translates to a mass splitting of about 420 keV. Although the CRESST-III [78] data are also available, the detector was optimized to sub-GeV mass dark matter detection with better performance at low recoil energies. Therefore, we still expect CRESST-II to set the leading bound at higher mass splittings.
In addition to NR interactions in traditional direct detection experiments, we are also interested in deexcitation gamma quanta in NE scattering: χ 1 + N → χ 2 + N * → χ 2 + N + γ, whose kinematics is also illustrated in Figure 2 for Os. Depending on the IDM model, χ 2 may also decay to SM or invisible particles. If χ 2 decays to photon(s) promptly inside the detector, the decay photon could be resolved as a signature of inelastic scattering; otherwise χ 2 may leave the detector. We conservatively assume the latter. For both NR and NE, the scattering rate is given by where N Ti is the number of the target nuclei of a given isotope in the detector, ρ χ = 0.4 GeV cm −3 is the local dark matter density, and f (v) is the one-dimensional Maxwellian velocity distribution. The minimum and maximum recoil energies can be found from Eq. (5). The differential cross section is where σ n is the per-nucleon scattering cross section and the nuclear response is encapsulated in S(E R ). For NR, this is just the Helm form factor. In the case of NE, dark matter only interacts with the few valence nucleons on the surface shell of the nucleus, where the nuclear response depends on the momentum transfer without the A 2 enhancement factor in S(E R ). As first derived in [80], the nuclear response function is written as where the sum runs over all allowed even angular momentum states within the range Here, J i and J f are respectively the inital and final angular momenta, j L are spherical Bessel functions and Y LM are spherical harmonics. The calculation of the response function is usually quite involved for heavy nuclei due to the complexity of the valence structure and the deformation of the nucleus. Fortunately, we may use the reduced transition probability in the nuclear E2 transition to approximately eliminate the angular matrix element, which is defined as with Y 2 the spherical tensor operator. Assuming the nuclear transition probability density is peaked near the surface of the nucleus, we find [81] with the nucleus radius R 1.2 × A 1/3 fm. The measured reduced transition probabilities B(E2) are listed in Table I. These are expressed in terms of the Weisskopf unit (W.u.) [82] defined as [83] for a general Eλ transition. Since higher multipole transitions involving j 2 L are strongly suppressed, we only consider the leading contribution to the nuclear response where L = 2.

III. EXPERIMENTS
Apart from a few notable exceptions, there have been no dedicated experiments searching for DM-nucleus interactions in heavy elements Z 54. However, these isotopes have been considered in many calibration experiments as well as searches for rare nuclear decays with the objective of better understanding their nuclear properties. Here, we list a number of promising target compounds, along with published data. We separate these out into two categories: NR, leading to nuclear recoils identifiable in scintillating bolometers, and NE, using gamma ray detectors.
Our strategy is to identify events which could have resulted from a DM-nucleus interaction, and translate the corresponding rate into a limit on the DM interaction rate. This is a conservative approach which, in the absence of a detailed background model for each experiment, can only be used for limit-setting rather than discovery. We will recast these constraints on the rates into constraints on IDM model parameters in Sec. IV.

A. Nuclear recoil measurements
CaWO 4 . Crystals used by the CRESST (II-III) [84,85] experiments contain a sizable fraction of tungsten, and are thus a good target for IDM searches. However, the analysis region for the CRESST DM runs was limited to low recoil energies, so we turn to the results of the crystal characterization runs. The radiopurity of a variety of CaWO 4 crystals was investigated at the Gran Sasso National Laboratory (LNGS). We employ measurements between 2009 and 2011 with the CRESST main setup where the scintillation light and phonon signals were both collected. The radiopurity results for a single crystal named Daisy are shown in Figure 2 of [60], representing 90.10 kg·days of exposure. The two-channel readout (scintillation and phonon) facilitates background event discrimination. The phonon channel was calibrated with α particles, which typically produce a light yield that is around 22% of the β/γ light yield in CaWO 4 crystals. Even less light is expected for nuclear recoils: for example, the light yield for nuclear recoils from a neutron calibration run of CRESST using CaWO 4 was found to be below 11% [86]. Conservatively, we will set bounds by counting all events observed in Daisy below the α−band as dark matter candidates. Although the nuclear recoil energy is not available due to the lack of dedicated neutron calibration, we approximate it by the equivalent α energy up to an uncertainty of about 20%, which does not impact our results in a notable way. With this criterion we identify three dark matter event candidates between 300 keV and 2 MeV, which yield an upper limit of 6.7 events at 90% CL. Since the bolometric detector is sensitive to the total energy deposition in the heat channel, the efficiency is close to 1 in the energy range we consider.
PbWO 4 . Similarly, the α decay of lead isotopes was studied with 11.09 kg·days of PbWO 4 background measurements, by using the crystal as a scintillating bolometer at the LNGS [61]. A lower light yield is expected from crystals of this type, which for α particles can be estimated via [61] LY α (keV/MeV) = (0.28±0.01)+(2.93±0.14)·10 −5 · E α keV , where E α is the α particle energy. Taking the events with a lower light yield than expected from α particles as dark matter candidates, we identify 12 events between 600 keV and 2 MeV, which corresponds to an upper limit of 17.8 events at 90% CL. Likewise we approximate the nuclear recoil energy by the α-equivalent energy since most energy is deposited in the heat channel in both cases. The detection efficiency is also close to 1. The results of CaWO 4 and PbWO 4 experiments are summarized in Table II

B. Nuclear excitation measurements
Os. Recently, the search for rare alpha decay of osmium isotopes was conducted at the LNGS [89]. 117.96 grams of high purity Os metal were cut into thin slices to reduce gamma absorption inside the metal. The slices were installed on a Broad-Energy Germanium detector with high sensitivity and ultra-low background. The accumulated energy spectrum was reported, with an exposure of 15,851 hours. The relevant transition levels and the corresponding background rates of Os isotopes are summarized in Table I. The deexcitation gamma quanta will appear as a peak in the continuous spectrum. Nonobservation of the peak implies that the 68% C.L. upper limit on the transition rate (marked as "Bkg") can be estimated as where s is the accumulated counts per keV as indicated by the measured energy spectrum. The 90% C.L. limit can be estimated by assuming gaussian distribution, where the transition rate is multiplied by factor of 1.645. The energy resolution or full width at half maximum (FWHM) is approximated by [89] σ keV = 0.57(5) + 0.029(2) × E γ /keV , where E γ is the transition energy. M and t exp are given by the mass of an Os isotope and the exposure time  [87], CRESST-II [79], XENON1T [88] and 180 Ta experiment [63] have been rederived accordingly and are shown with the shaded regions. The green, blue and magenta lines show the cross section limits derived from inelastic dark matter scattering which induces nuclear transitions of Hf [67], Os [89], Hg and the deexcitation gamma rays. Dashed lines are for future projections. The solid purple and orange lines depict the limits on inelastic dark matter from the recoil of heavy nuclei, derived in this work using data from CaWO4 [60] and PbWO4 [61] crystals respectively (see text for details).
respectively. The efficiency of HPGe detector depends on details of the sample including geometry, density and chemical composition, and the energy of gamma quanta. In the absence of the original experimental data, we construct a polynomial fit to the efficiency curve as a function of ln(E γ ) [90] between 46.5 keV and 200 keV, and take the efficiency to be constant below 46.5 keV 1 . A dedicated analysis may yield an O(1) correction to the efficiency adopted here. If intermediate states exist between the excited nucleus and its ground state, deexcitation may take place through the emission of gamma cascades. Picking a specific gamma energy for background analysis entails the branching ratio I γ , which is 100% for all transitions listed in Table I. We therefore search for the deexcitation gamma and set bounds on the cross section by requiring that the dark matter scattering rate in Eq. (9) not exceed the measured transition rate.
Hf. Rare decay of hafnium isotopes were also investigated at the LNGS by measuring the internal background of 55.38 grams of Hf foil using modified HPGe detector [67]. Low-background data has been collected over 310.6 days. The most constraining transitions and the background rates are also listed in Table I. The corresponding limits on inelastic dark matter have been derived in Ref. [67].
Hg. 201 Hg is an ideal target NE gamma search because of its low-lying excited state at 1.6 keV. Since Hg is toxic it is preferable to use mercury-containing compounds as the target for dark matter scattering. One such example is to use mercury cadmium telluride (Hg 1−x Cd x Te). While this would suffer from a high contamination from 113 Cd beta decays, which will dominate over other background below 316 keV, increasing the fraction of Hg in the crystal would somewhat diminish this background component. Conservatively we consider x = 0.2. The 113 Cd decay background at 1.6 keV can be obtained by rescaling the experimental data at 100 keV according to the predicted decay spectrum in Fig. 5 of Ref. [91]. This 1.6 keV gamma can be measured with an external lowthreshold gamma detector such as Ge cryogenic bolometer [92][93][94], which is to be installed in close contact with the Hg 0.8 Cd 0.2 Te target. The energy resolution of the Ge gamma detector is estimated to be 50 eV around 1.6 keV [94]. Assuming 10 kg·yr of Hg 0.8 Cd 0.2 Te exposure and ∼ 50% efficiency for gamma detection, we estimate the background rate to be 0.0004 mBq per kg of Hg 0.8 Cd 0.2 Te (or equivalently 0.0056 mBq per kg of 201 Hg) at 1.6 keV transition energy.
Although the DM analysis in this work is based on the measurement of the 1.6 keV gamma in DM induced NE, the DM limit can be further improved by investigating nuclear recoil signals. Since Hg 0.8 Cd 0.2 Te as a semiconductor is also sensitive to energy deposition, we can improve this setup by measuring the recoil energy of Hg with Hg 0.8 Cd 0.2 Te. The coincident detection of nuclear recoil at Hg 0.8 Cd 0.2 Te and 1.6 keV gamma at Ge bolometer as DM signals provides a powerful way for background rejection. The exploration of this method is left for further work.

IV. RESULTS
The constraints on dark matter-nucleon scattering based on the experiments described above are shown in Figure 3, where we have set the DM mass to M χ = 1 TeV. For dark matter much heavier than the target nuclei, the cross section bound scales linearly with M χ .
As an illustration of IDM-induced nuclear transition we show the bounds from the NE to different Os excited states in Figure 3. 189 Os sets the leading bound for mass splittings below 520 keV, due to the collisional excitation to the 9.54 keV state. Stronger bounds could be obtained from 505 keV to 553 keV by considering the future measurement of 36.17 keV transition gamma . We also set limits on inelastic dark matter by employing Hf rare decay measurements in a similar approach [67], which yields less stringent bounds compared with osmium owing to relatively higher transition energies and less massive nuclei. We also propose to detect inelastic dark matter by investigating the deexcitation gamma yield of a mercurycontaining compound. With 10 kg·yr of HgCdTe exposure as described above, 201 Hg is expected to yield leading limits amongst nuclear excitation IDM searches for mass splittings between 415 keV and 620 keV.
We have also derived the bounds for IDM NR scattering on CaWO 4 and PbWO 4 target crystals acting as scintillating bolometers. Thanks to the simultaneous light and phonon detection, nuclear recoils can be dis-tinguished on a statistical basis from electron recoil and gamma background. These bounds are more stringent than those from Hf and Os isotopes for all mass splittings. In particular, the CaWO 4 experiment dominates the 386 keV to 573 keV range, and the PbWO 4 measurements constrain the mass splitting onward up to 636 keV. Nevertheless, 201 Hg is still expected to lead the 570 keV and 604 keV mass splitting.
We also show the constraints from XENON1T, PICO-60 [87], CRESST-II [79] and 180m Ta [63]. The XENON collaboration reported the results for WIMP searches from about 1 tonne·yr of data collection [88] in 2018. We reproduce the results following the efficiency and events selection in the Xenon1T-2018 code [95]. 2 events were observed with 0.9 tonne reference mass and 278.8 days livetime, and the background is expected to be 1.62±0.28. This translates to an upper limit of 3.7 events at 90% C.L. Discrimination of nuclear recoil, electron recoil and other backgrounds on an event by event basis requires full likelihood analysis, which is beyond the scope of this work. We derive the bound on inelastic dark matter by integrating over the nuclear recoil energy between 0 and 40.9 keV, with the efficiencies taking care of the low-energy cut. The PICO-60 [87] and CRESST-II [79] bounds in [96] are recomputed with the dark matter velocity distribution assumed in this paper 2 and the Ta bound [63] is reproduced with the same assumptions.
We extend the bound to lower and higher dark matter masses in Figure 4. The colors depict the smallest cross sections that are excluded with the combination of experiments listed above. The mass splitting that can be constrained rises as dark matter mass increases. The white region is not kinematically accessible due to the limitation of µ χN and dark matter velocity, as is evident from Eq. (7). For δ 200 keV leading bounds are set by XENON1T. In the intermediate regime CaWO 4 experiment remains the most sensitive search, which is overtaken by PbWO 4 at high mass splittings. Current limits from dark matter induced nuclear excitations remain subdominant.
Apart from the target nucleus mass, the maximum DM velocity is a limiting factor for probing high mass splittings. Although DM is considered to be virialized, many recent works have shown that a subdominant component could be boosted to high velocity, through cosmicray (e.g. [97][98][99][100]) or neutrino (e.g. [101][102][103]) scattering or stellar acceleration (e.g. [104]), annihilation (e.g. [105,106]), semi-annihilation (e.g. [107,108]) or decay of heavy states (e.g. [109,110]). Cosmic-ray upscattered inelastic dark matter with a vector mediator has recently been studied in Ref. [111]. Regardless of boost mechanism, we illustrate the kinematic reach of future heavy nuclei experiments to boosted IDM particles in Figure 4, for DM velocities of 1000 km/s (dot- FIG. 4. Constraints on inelastic dark matter-nucleon scattering for different dark matter masses and mass splittings. The colored regions to the right of the solid lines are excluded by existing experiments, including XENON1T [88], CRESST-II [79], and the bolometric searches with CaWO4 [60] and PbWO4 [61] scintillating crystals. Colors and contour lines mark the minimum cross section at 90% C.L. that can be excluded by the combination of experiments. The white region on the upper left is not kinematically accessible in existing experiments. Each experiment is sensitive to a limited amount of parameter space as enclosed by the solid lines. If a subcomponent of the dark matter is boosted to higher velocities, more parameter space can be excluded. The dash-dotted lines show the projected reach for maximum dark matter velocities of vmax = 1000 km/s, and the dashed lines for vmax = 2000 km/s. The green, orange and red broken lines respectively correspond to future dark matter search through scattering with tungsten (nuclear recoil), lead (nuclear recoil) and mercury (nuclear excitation), regardless of detection threshold.
dashed lines) and 2000 km/s (dashed lines) assuming W, Pb (NR) and Hg (NE) nuclei as targets. The full parameter space in M χ − δ shown in Fig. 4 can largely be covered by if v max = 2000 km/s. We choose to cap the mass splitting at MeV, beyond which the heavier state may decay to e ± in the detector. An exploration of the prompt decay scenario is left for future work.

V. CONCLUSIONS
As direct detection experiments close in on the window of possible WIMP candidates in the GeV range, inelastic dark matter remains comparatively unconstrained by experiments with light nuclei, due to kinematical suppression, as can be seen from Eq. (7). We have extended the scope of this search by investigating nuclear recoil and nuclear excitation events in experiments with heavy nuclei using two distinct experimental approaches.
In the first scenario, target nuclei recoil due to the momentum transferred by dark matter scattering. In a scintillating bolometer, recoil energies up to about 2 MeV will exhibit distinct features from electron/gamma background due to the low quenching in nuclear recoil. We reinterpret the high phonon energy events with low scintillation yield as dark matter scattering candidate events in experiments with CaWO 4 and PbWO 4 scintillating bolometers and set bounds on mass splitting up to 636 keV for TeV dark matter. This is the largest splitting explored thus far in the literature. Our limits are conservative as we do not have access to the raw data in the experiments. The discrimination of dark matter candidates from background on an event-by-event basis will strengthen the bounds on inelastic dark matter. The constraints can be further improved in the future with low-background scintillating crystals [112]. In particular, RES-NOVA [113] experiment with large exposure and better background modeling using Monte Carlo, may significantly improve the limits on inelastic dark matter if it is repurposed for dark matter search. We have also revisited the bounds from XENON1T, PICO-60, CRESST-II and Ta. A likelihood analysis or the study of the high recoil energy data in XENON1T may notably extend the current Xenon bound below 450 keV.
In the second scenario, target nuclei are excited to higher energy states through inelastic collision with dark matter. In the absence of any background process that could yield such an excitation, the deexcitation gamma photons are unique signals of dark matter. In contrast to the difficult measurements with cryogenic scintillating bolometers calibrated at high recoil energies, the detection of nuclear excitation is simply viable with HPGe detector, at the cost of excitation energy and the loss of coherent enhancement in the interaction cross section. We have shown that a recent Os deexcitation gamma search constrains the mass splitting up to 553 keV, and a future Hg experiment could extend the sensitivity up to 620 keV.
A key feature of WIMPy inelastic dark matter resides in its kinematics tightly coupled to the nuclear masses. As we are marching towards the end of the period table, heavy elements are exploited to tackle the high mass splitting regime, rapidly shrinking this large parameter space. Dedicated high mass, high energy deposition experiments are required to perfect this approach, in combination with more sophisticated analyses of the dark matter flux at different velocities.

ACKNOWLEDGMENTS
We thank Joseph Bramante and Peng Wang for helpful discussions, and Harikrishnan Ramani for correspondence on Ta constraints [63]. We also thank the CRESST group at MPI Munich and Alexander Derbin for openly discussing their Tm data and fruitful collaboration, and Stefano Pirro for communications on Ge detector. SN We explore more isotopes in this appendix. Hf, Os and Hg isotopes and their transitions studied in this work are listed in Table III. Bounds on IDM derived accordingly are shown in Figure 5. For constraints from individual Hf isotope, see [67]. We do not include 187 Os 74.4 keV excitation in the analysis since the measured B(E2) is subject to large uncertainty. The projected background rates of transitioning to 201 Hg 26.3 keV and 32.1 keV excited states are derived similarly to 1.6 keV state by rescaling the 113 Cd decay data [91] and by assuming 50% efficiency.
We find the bounds from 187 Os 75.0 keV and 189 Os 95.3 keV transitions are inferior to the 189 Os 69.5 keV constraints. Limits from 201 Hg 26.3 keV and 32.1 keV are also weaker than 1.6 keV due to the kinematical suppression. Other elements that feature low transition energy and high mass include 181 Ta, 183 W, 193 Ir and 197 Au, which could be adopted as nuclear targets for future IDM searches. In Fig. 6 we show the dark matter candidate events in the CaWO 4 scintillating bolometer experiment. The light yields (LYs) of O, Ca and W in nuclear recoil have been studied in Ref. [118], which all approach a constant with LY ∞ being the LY at E R = ∞. The relevant parameters are given in Table IV. The LY of W in nuclear recoil is measured to be 0.0208 ± 0.0024. The LYs are depicted as solid and dashed lines in Fig. 6 at 90% C.L. (again gaussian distribution is assumed to infer the 90% C.L. from 1σ range), which lie below the α-band from Ref. [60]. This reconfirms the low LY of nuclear recoil compared with α's. Conservatively we assume the events below the purple α-band as dark matter candidates, where three such events are found for the phonon energy between 300 keV and 2000 keV. We do not include the events between 1800 keV and 2000 keV, as they are attributable to 144 Nd decays and are far away from the W LY, which dominates the dark matter scattering at high recoil energy due to kinematical suppression.
We show the 1σ range of the corresponding LY with blue lines. To be conservative, we assume events below the 1σ α-band to be dark matter candidates, where 12 such events are identified.
FIG. 6. Event discrimination in the CaWO4 scintillating bolometer experiment. The purple band corresponds to the degraded α-events in Ref. [60]. The blue, and green and magenta lines show the measured light yields of O, Ca and W in nuclear recoil at 90% C.L. [118]. The light yield of W recoil is measured up to ∼240 keV recoil energy, and that of O and Ca recoils are measured up to ∼500 keV, beyond which we depict light yields as dashed lines. Events painted red below the purple alpha-band are considered as dark matter candidates. FIG. 7. Event discrimination in the PbWO4 scintillating bolometer experiment. The green and purple lines correspond to the 3 σ acceptance region for β/γ and α events adapted from Ref. [61] respectively. The blue lines represent the 1 σ region of α events based on the light yield of α particles. Events painted red below the blue lines are considered as dark matter candidates.