Inelastic Fermion Dark Matter Origin of XENON1T Excess with Muon $(g-2)$ and Light Neutrino Mass

Motivated by the recently reported excess in electron recoil events by the XENON1T collaboration, we propose an inelastic fermion dark matter (DM) scenario within the framework of a gauged $L_{\mu}-L_{\tau}$ extension of the standard model which can also accommodate tiny neutrino masses as well as anomalous muon magnetic moment $(g-2)_{\mu}$. A Dirac fermion DM, naturally stabilised due to its chosen gauge charge, is split into two pseudo-Dirac mass eigenstates due to Majorana mass term induced by singlet scalar which also takes part in generating right handed neutrino masses responsible for type I seesaw origin of light neutrino masses. The inelastic down scattering of heavier DM component can give rise to the XENON1T excess for keV scale mass splitting with lighter DM component. We fit our model with XENON1T data and also find the final parameter space by using bounds from $(g-2)_{\mu}$, DM relic, lifetime of heavier DM, inelastic DM-electron scattering rate, neutrino trident production rate as well as astrophysical, cosmological observations. A tiny parameter space consistent with all these bounds and requirements will face further scrutiny in near future experiments operating at different frontiers.

Introduction: XENON1T collaboration has recently reported an excess of electron recoil events near 1-3 keV energy [1]. While this excess is consistent with the solar axion model at 3.5σ significance and with neutrino magnetic moment signal at 3.2σ significance, both these interpretations are in strong tension with stellar cooling constraints. While XENON1T collaboration can neither confirm or rule out the possible origin of this excess arising due to beta decay occurring in trace amount of tritium present in the xenon container, it has generated a great deal of interest among the particle physics community to look for possible new physics interpretations. Different dark matter (DM) interpretations of this excess have been proposed in several works . For other interpretations and discussions related to this excess, please refer to . In the present work, we adopt the idea of inelastic DM in the context of XENON1T excess within the framework of a well motivated gauge extension of the standard model (SM).
One popular extension of the SM is the implementation of an Abelian gauge symmetry L α − L β where L α is the lepton number of generation α = e, µ, τ . Interestingly, such a gauge extension is anomaly free and can have very interesting phenomenology related to neutrino mass, DM as well as flavour anomalies like the muon anomalous magnetic moment (g − 2) µ [72]. While there can be three different combination for this gauge symmetry, we particularly focus on L µ − L τ gauge symmetry. For earlier works in different contexts, please see [73][74][75][76]  If the mass splitting between these two mass eigenstates is appropriately tuned, the heavier component can be long-lived and can comprise a significant fraction of total DM density in the present universe. Here we show that inelastic fermion DM can give rise to the required electron recoil events observed by XENON1T while at the same time being consistent with relic abundance, muon (g − 2) and light neutrino mass criteria. Similar idea of addressing muon (g − 2) and XENON1T excess within a scalar extension of the SM was recently proposed in [40]. Another recent work, particularly in the context of L µ − L τ gauge symmetry, showed that solar neutrinos with such new gauge interactions can not be responsible for XENON1T excess [80]. Our proposal in this work provides an alternative way to address the excess in gauged L µ − L τ model augmented by inelastic fermion DM.
Gauged L µ − L τ Symmetry: As mentioned before, we consider an L µ − L τ gauge extension of the SM. The SM fermion content with their gauge charges under Lµ−Lτ gauge symmetry are denoted as follows.
Note that the chiral fermion content of the model mentioned above keeps the model free from triangle anomalies. The DM field is represented by a Dirac fermion χ L,R ∼ (1, 1, 0, 1 2 ) where the choice of L µ − L τ charge is made in such a way that stabilise it without requiring any additional symmetries. In order to break the gauge symmetry spontaneously as well as to generate the desired fermion mass spectrum, the scalar fields are chosen as follows.
While the neutral component of the Higgs doublet H breaks the electroweak gauge symmetry, the singlets break L µ − L τ gauge symmetry after acquiring non-zero vacuum expectation values (VEV). Denoting the VEVs of singlets φ 1,2 as v 1,2 , the new gauge boson mass can be found to be M Z = g x (v 2 1 + 4v 2 2 ) with g x being the L µ − L τ gauge coupling. Please note that, in principle, the symmetry of the model allows a kinetic mixing term between U (1) Y of SM and U (1) Lµ−Lτ of the form is the mixing parameter. This kinetic mixing plays a crucial role in giving rise to the XENON1T excess as we discuss later.
The relevant part of the DM Lagrangian is The DM field is identified as χ which has a bare mass as well as coupling to φ 1 . While the bare mass term is of Dirac type, the coupling to φ 1 introduces a Majorana mass term after φ 1 acquires a non-zero VEV. Thus, the Dirac fermion χ is split into two Majorana fermions χ 1 , χ 2 . The DM Lagrangian in this physical basis is As will be discussed below, the mass splitting between χ 1 , χ 2 is chosen to be very small δ = M 2 − M 1 = 2m + ∼ O(keV) in order to give the required fit to XENON1T excess. This ensures M 1 ≈ M 2 ≈ M χ while leaving m − as a free parameter.
In the above Lagrangian for DM, θ is a mixing angle given by tan θ ≈ m − /M χ .
Light Neutrino Masses: In order to account for tiny non-zero neutrino masses for light neutrinos, we extend the minimal gauged U (1) Lµ−Lτ model with additional neutral fermions as where the quantum numbers in the parentheses are the gauge charges under SU Also, the chosen gauge charges of right handed neutrinos do not introduce any new contribution to triangle anomalies. The relevant Yukawa interaction terms are given by Which clearly predict diagonal charged lepton and Dirac neutrino mass matrices M , M D . Thus, the non-trivial neutrino mixing will arise from the structure of right handed neutrino mass matrix M R only which is generated by the chosen scalar singlet fields. The ight handed neutrino, Dirac neutrino and charged lepton mass matri-ces are given by Using type I seesaw approximation (M D M R ), the light neutrino mass matrix can be found from the following seesaw formula Since the charged lepton mass matrix is diagonal in our model, the light neutrino mass matrix can be diagonalised by using the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) mixing matrix as where m i are the light neutrino mass eigenvalues. Since M R has a very general structure, one can fit the model parameters with light neutrino data in several ways, independently of rest of our analysis.
Anomalous Muon Magnetic Moment: The magnetic moment of muon is given by where g µ is the gyromagnetic ratio and its value is 2 for a structureless, spin 1 2 particle of mass m and charge q. Any radiative correction, which couples the muon spin to the virtual fields, contributes to its magnetic moment and is given by The anomalous muon magnetic moment has been measured very precisely while it has also been predicted in the SM to a great accuracy. At present the difference between the predicted and the measured value is given by which shows there is still room for NP beyond the SM (for details see [72]). In a recent article, the status of the SM calculation of muon magnetic moment has been updated [81]. According to this study ∆a µ = (27.9 ± 7.6) × 10 −10 which is a 3.7σ discrepancy. Quoting the errors in 3σ range, we have ∆a µ = (27.9 ± 22.8) × 10 −10 In our model, the additional contribution to muon magnetic moment comes from one loop diagram mediated by Z boson. The contribution is given by [82] where α = g 2 x /(4π).
Relic Abundance of DM: Relic abundance of two component DM in our model χ 1,2 can be found by numerically solving the corresponding Boltzmann equations. Let n 2 = n χ2 and n 1 = n χ1 are the total number densities of two dark matter candidates respectively. The two coupled Boltzmann equations in terms of n 2 and n 1 are given below, (n eq 2 ) 2 (n eq 1 ) 2 n 2 1 − σv χ2χ1→XX (n 1 n 2 − n eq 1 n eq 2 ) , dn 1 dt + 3n 1 H = − σv χ1χ1→XX n 2 1 − (n eq 1 ) 2 + σv χ2χ2→χ1χ1 n 2 2 − (n eq 2 ) 2 (n eq 1 ) 2 n 2 1 − σv χ2χ1→XX (n 1 n 2 − n eq 1 n eq 2 ) , where, n eq i is the equilibrium number density of dark matter species i and H denotes the Hubble parameter. The thermally averaged annihilation and coannihilation processes (χ iχi → XX) are denoted by σv where X denotes all particles to which DM can annihilate into.
Since we consider GeV scale DM, the only annihilations into light SM fermions can occur. We consider all the singlet scalars to be heavier than DM masses. Also, the singlet mixing with SM Higgs are assumed to be tiny so that singlet mediated annihilation channels are neg-ligible and only the annihilations mediated by Z gauge boson dominate. Additionally, the keV scale mass splitting between the two DM candidates lead to efficient coannihilations while keeping their conversions into each other sub-dominant. We have solved these two coupled Boltzmann equations using micrOMEGAs [83]. Due to tiny mass splitting, almost identical annihilation channels and sub-dominant conversion processes, we find almost identical relic abundance of two DM candidates. Thus each of them constitutes approximately half of total DM relic abundance in the universe. We constrain the model parameters by comparing with Planck 2018 limit on total DM abundance Ω DM h 2 = 0.120 ± 0.001 [84]. Here Ω DM is the density parameter of DM and h = Hubble Parameter/(100 km s −1 Mpc −1 ) is a dimensionless parameter of order one.
Since the mass splitting between χ 2 and χ 1 is kept at keV scale δ ∼ O(keV), there can be decay modes like χ 2 → χ 1 νν primarily mediated by Z . If both the DM components are to be there in the present universe, this lifetime has to be more than the age of the universe that is τ χ2 > τ age ≈ 4 × 10 17 s. The decay width of this process is Γ χ2→χ1νν ≈ g 4 x δ 5 /(160π 3 M 4 Z ). Thus, imposing the lifetime constraint on heavier DM component puts additional constraints on the model parameters.
XENON1T Excess: Similar to the proposal in [9], here also we consider the down-scattering of heavier DM component χ 2 e → χ 1 e as the process responsible for XENON1T excess of electron recoil events near 1-3 keV energy [1].
For a fixed DM velocity v,the differential cross section is given by where m e is the electron mass, a 0 = 1 αme is the Bohr radius, α = e 2 4π = 1 137 is the fine structure constant, E r is the recoil energy, q is the transferred momentum, K(E r , q) is the atomic excitation factor, and σ e is the free electron cross section. The atomic excitation factor is taken from [85]. We assume the DM form factor to be unity. The free electron cross-section is given by where α z = g 2 4π , α x = g 2 x 4π and is the kinetic mixing parameter between Z and Z mentioned earlier which we take to be ≤ 10 −4 . Here in the inelastic scattering case, the limits of integration in Eq. (11) are determined depending on the relative values of recoil energy (E r ) and the mass splitting between the DM particles (δ = M 2 − M 1 ). It should be noted that σ e is independent of DM mass as the reduced mass of DM-electron is almost equal to electron mass for GeV scale DM mass we are considering.
And for E r ≤ δ The differential event rate for the inelastic DM scattering with electrons in xenon is given by where n T = 4 × 10 27 Ton −1 is the number density of xenon atoms and n χ2 is the density of the dark matter χ 2 . As mentioned before n χ2 ≈ n χ1 ≈ n DM /2.

Results and Conclusion:
We first fit our model with XENON1T data using the methodology described above. The result is shown in Fig. 1. The mass splitting is taken to be δ = 2 keV while heavier DM mass is taken in sub-GeV regime consistent with all relevant constraints. DM velocity is taken to be v ≈ 8 × 10 −3 , consistent with its non-relativistic nature. The other relevant parameters used in this fit are g x = 10 −3 , M Z = 0.5 GeV, = 3 × 10 −3 . As we discuss below, this choice of parameters is also consistent with all other relevant bounds. We then calculate the relic abundance of two DM candidates χ 2 , χ 1 using the procedures mentioned above. The left panel of Fig. 2 shows the variation of DM relic abundance with DM mass for a set of fixed benchmark parameters. Clearly, due to tiny mass splitting between two DM candidates and identical gauge interactions, their relic abundances are almost identical. The DM annihilation due to s-channel mediation of Z gauge boson is clearly visible from this figure where correct relic of DM is satisfied near the resonance region M DM ≈ M Z /2.
Final result is summarised in the right panel plot of lines. The black dashed line corresponds to the bound on cross sections for νN → νN µμ measured by CCFR [86], the region above this line is ruled out. This constraint on g x − M Z plane arises purely due to the fact that L µ − L τ gauge boson can contribute to this neutrino trident process. It completely rules out the parameter space satisfying (g − 2) µ at 3σ beyond M Z 1 GeV. The dashed line of dark red colour shows the bound on lifetime of heavier DM χ 2 mentioned earlier so that the region above this line is ruled out. Clearly, this lifetime bound is stronger than the CCFR bound for M Z 1 GeV. Interestingly, the CCFR and lifetime bounds allow only a small triangular region in the (g −2) µ favoured parameter space around M Z ∼ 0.5 GeV (see inset of right panel plot in Fig. 2). The pink solid line corresponds to σ e = 10 −17 GeV −2 required to fit the XENON1T excess.
We also use the strong astrophysical bounds from white dwarf (WD) cooling on such light gauge bosons [87]. This arises as the plasmon inside the WD star can decay into neutrinos through off-shell Z leading to increased cooling efficiency. This leads to a bound in the g x − M Z parameter space as [88]  However, in the region of our interest (triangular region allowed from CCFR and lifetime bounds), the WD cooling constraint remains weaker compared to other relevant bounds, as can be seen from the green dotted line in Fig. 2  (right panel).
We then consider the cosmological bounds on such light DM and corresponding light mediator gauge bo-son Z . A light gauge boson can decay into SM leptons at late epochs (compared to neutrino decoupling temperature T ν dec ∼ O(MeV) increasing the effective relativistic degrees of freedom which is tightly constrained by Planck 2018 data as N eff = 2.99 +0. 34 −0.33 [84]. As pointed out by the authors of [88,89], such constraints can be satisfied if M Z 10 MeV. As can be seen from the right panel plot in Fig. 2, the lifetime requirement of χ 2 already puts a much stronger bound in the region of our interest. Similarly, constraints from cosmic microwave background (CMB) measurements disfavour such light sub-GeV thermal DM production in the early universe through s-channel annihilations into SM fermions [84]. We check that in our model, the DM annihilation processes are p-wave suppressed, relaxing the bounds from CMB, similar to the scenario proposed in [9].
Finally, we perform a random scan for relic abundance of two component DM so that their combined relic satisfy the criteria for observed DM relic abundance. This is shown in terms of scattered points in right panel plot of Fig. 2 where the colour coding is used to denote DM mass. In this random scan, apart from varying g x , M Z we also vary DM mass M DM ∼ M 1 ∼ M 2 in the range (0.05, 3) GeV and the other free parameter m − in the range (0.1, 1) GeV while keeping the tiny mass splitting fixed at δ = 2 keV. Clearly, only a very few points fall in the small triangular region allowed from all constraints and requirements. The density of these points inside the triangular region will increase for a bigger scan size. Since only a tiny region of parameter space is allowed in this model, more precise measurements of (g − 2) µ will be able to confirm or rule out this model as its possible explanation. Future measurements by XENON1T collaboration will also give a clearer picture on the feasibility of this excess. Proposed experiments like NA62 [90] will also be able to rule out or confirm some part of the parameter space discussed in our work. We leave a detailed model building and phenomenological study of such low mass DM scenario in the context of electron recoil signatures to future work.