-
The absence of a right-handed neutrino component in the Standard Model (SM), coupled with the imperatives of
$ SU(2)_L $ gauge invariance and renormalizability, results in the prediction that neutrinos should be massless. However, experimental evidence from neutrino oscillation studies [1−4] compellingly demonstrates that neutrinos possess a non-zero, albeit minute, mass. This discrepancy hints the existence of physics that extends beyond the current understanding provided by the SM.In 1937, Majorana proposed that neutrinos might be their own antiparticles, a scenario in which the neutrino and the antineutrino are indistinguishable [5]. This possibility raises an unresolved question: are neutrinos Dirac particles, distinct from their antiparticles, or Majorana particles, which are not? The ``see-saw" mechanism [6−8] is a popular model for generating light neutrino masses, requiring a heavy Majorana neutrino with a mass ranging from a few hundred MeV/
$ c^2 $ to a few GeV/$ c^2 $ . If neutrino masses arise from this mechanism, there could be processes that violate lepton number by two units ($ \Delta L=2 $ ), revealing the Majorana nature of neutrinos. Thus, search for lepton number violation (LNV) is crucial to confirm whether neutrinos are Majorana particles. In 1939, based on Majorana's theory, Furry [9] studied the neutrinoless double beta ($ 0\nu2\beta $ ) decay [10], which is an LNV process. Since then, various theories concerning LNV processes have been proposed, and numerous experiments have been carried out to explore these processes.In addition, the observed matter-antimatter asymmetry in the Universe presents a significant challenge to our understanding of nature. The Big Bang theory, which is the prevailing cosmological model for the Universe's evolution, predicts equal amounts of baryons and anti-baryons in the early Universe. However, the observed amount of baryons exceeds that of the anti-baryons by a factor currently estimated at
$ 10^9 - 10^{10} $ [11]. To give a reasonable interpretation of the baryon-antibaryon asymmetry, Sakharov proposed three principles [12], the first of which is that baryon number conservation must be violated. Many theories believe that if there is baryon number violation, there will also be LNV[13]. Consequently, experimental searches for the LNV$ \Delta L = 2 $ processes are of great interest.Apart from the
$ 0\nu2\beta $ decay experiments, which could search for Majorana neutrinos, several collider experiments have explored related processes. For example, the LHCb experiment investigated the$ B^{-} \to \pi^{+} \mu^{-} \mu^{-} $ decays [14], the CMS experiment studied the$ p\bar{p}\to e^+e^+ X $ process [15], the BaBar experiment analyzed the$ D^{+}_{s} \to \pi^{+} l^{+} l^{+} $ decays [16], and the ATLAS experiment searched for right-handed W bosons and heavy right-handed Majorana or Dirac neutrinos, using final state containing a pair of charged leptons and two jets [17]. The Gerda experiment searched for the rare$ 0\nu2\beta $ decay process in$ ^{76}Ge $ [18]. In addition, the CLEO experiment investigated the$ D \to \pi(K) e^{+}e^{+} $ decays [19] and the FOCUS experiment searched for the decays of$ D^{+} $ and$ D^{+}_{s} $ containing the same-sign di-muons [20]. In recent years, the BESIII experiment searched for the LNV decays of charmed mesons, including$ D \to K \pi e^{+} e^{+} $ [21],$ D^{0} \to \bar{p}e^{+}(pe^{-}) $ [22],$ D^{+} \to \bar{\Lambda}(\bar{\Sigma}^{0}e^{+}) $ [23], and$ D^+(D^-)\to ne^+(\bar n e^-) $ [24]. In all above searches, no evidence has been found.In the τ-charm energy region, the huge data set of charmonia (such as
$ J/\psi $ ) can be analyzed to search for various LNV processes, and it is feasible to search for$ J/\psi \to K^+K^+e^-e^- $ (throughout this paper, charge conjugation is implied) and other LNV processes in the decay products of$ \phi, \omega, K^0, \eta/\eta' $ , and so on. Fig. 1 shows Feynman diagram of$ J/\psi \to K^+K^+e^-e^-+c.c $ through possible high dimension operator [25]. Search for this process would provide hints for new physics, as they are expected to be suppressed due to higher-order processes in the SM.In this paper, by analyzing
$ (10087\pm 44)\times10^{6} $ $ J/\psi $ [26] events collected with the BESIII detector [27] operating at the BEPCII storage ring, we search for the LNV decay$ J/\psi \to K^+K^+e^-e^- $ for the first time. Furthermore, we explore the upper limits on the branching fractions of this decay under different Majorana neutrino mass hypotheses. To avoid possible bias, a blind analysis technique is employed, where the data are analyzed only after the analysis procedure has been finalized and validated with MC simulation. -
The BESIII detector is located at the south collision point of the BEPCII [28] storage ring, which operates at a center-of-mass energy range of
$ \sqrt{s} $ = 1.84 to 4.95$ {\rm{GeV}} $ , with a maximum luminosity of$ 1.1 \times 10^{33} \; \rm cm^{-2}s^{-1} $ at 3.773$ {\rm{GeV}} $ .The cylindrical core of the BESIII detector covers 93% of the full solid angle and consists of a helium-based multilayer drift chamber (MDC), a plastic scintillator time-of-flight system (TOF), and a CsI(Tl) electromagnetic calorimeter (EMC). These components are enclosed within a superconducting solenoidal magnet that provides a 1.0 T magnetic field (0.9 T in 2012). The solenoid is supported by an octagonal flux-return yoke, with resistive plate counter muon identifier modules interleaved with steel. The BESIII coordinate system is defined such that the z-axis aligns with the beam direction, the x-axis is horizontal and points outward from the interaction point, and the y-axis is vertical, forming a right-handed coordinate system.
At a momentum of 1
$ {\rm{GeV}}/c $ , the charged particle momentum resolution is 0.5%, and the specific ionization energy loss (dE/dx) resolution is 6% for electrons from Bhabha scattering. The EMC has a photon energy resolution of 2.5% (5%) at 1$ {\rm{GeV}} $ in the barrel (end-cap) region. The time resolution in the TOF barrel region is 68 ps, while in the end-cap region it is 110 ps. In 2015, the end-cap TOF system was upgraded using multigap resistive plate chamber technology, improving the time resolution to 60 ps [29]. Approximately 87% of the data used in this analysis benefit from this upgrade.Simulated samples produced with the GEANT4-based [30] Monte Carlo (MC) program which includes the geometric description [31, 32] of the BESIII detector and the detector response, are used to determine the detection efficiency and estimate the backgrounds. The simulation includes the beam energy spread and initial state radiation in the
$ e^+e^- $ annihilations modeled with the generator KKMC [33]. The inclusive MC sample includes both the production of the$ J/\psi $ resonance and the continuum processes incorporated in KKMC [33]. All particle decays are modelled with EVTGEN [34] using branching fractions either taken from the Particle Data Group [35], when available, or otherwise estimated with LUNDCHARM [36]. Final state radiation (FSR) from charged final state particles is incorporated using the PHOTOS package [37]. In this work, a sample of 100,000 phase space (PHSP) MC events are generated to study the detection efficiencies for the LNV decay, with equal contributions from each conjugate channel. For background analysis, 10 billion inclusive MC events are used. The continuum background is estimated using off-resonance data samples at$ \sqrt{s} = 3.650 $ ,$ 3.682 $ and$ 3.080 $ $ {\rm{GeV}} $ [38], with integrated luminosities of$ 410\ {\rm{pb^{-1}}} $ ,$ 404\ {\rm{pb^{-1}}} $ and$ 167\ {\rm{pb^{-1}}} $ , respectively. -
All charged tracks are required to satisfy a geometrical acceptance of
$ |{\rm{cos}}\theta| < 0.93 $ , where θ is the polar angle of the charged track, defined with respect to the z-axis, which is the symmetry axis of the MDC. Each track must originate from the interaction region, defined as$ R_{xy} <1.0 $ cm and$ |R_{z}| < 10.0 $ cm, where$ R_{xy} $ and$ R_{z} $ are the distances of the closest approach to the interaction point of the track in the transverse plane and z-axis, respectively. Events with at least four charged tracks with zero net charge are retained for further analysis.Particle identification (PID) for charged tracks combines measurements of the energy deposited in the MDC (dE/dx), the flight time in the TOF, and the energy and shape of clusters in the EMC to calculate the confidence level (CL) for electron, pion, kaon hypotheses (CL
$ _e $ , CL$ _{\pi} $ , CL$ _K $ ). Electron candidates must have CL$ _e $ greater than 0.001, and the ratio of CL$ _e $ to the sum of CLs for electrons, kaons, and pions (CL$ _e $ + CL$ _K $ + CL$ _\pi $ ) should be greater than 0.8. In addition, the ratio of energy deposited in the EMC(E) to the momentum measured in the MDC(p) must lie between 0.8c and 1.2c. For kaon candidates, the criteria include$ {\rm{CL}}_K>0 $ and$ {\rm{CL}}_K> {\rm{CL}}_{\pi} $ . To further reduce the misidentification(mis-PID) between electrons and kaons, an additional requirement of$ E/p $ less than 0.8c is applied [39]. The charges of the two selected kaons must be identical and opposite to the charges of the two selected electrons, which are also required to be the same. The minimum number of kaons and electrons are two. If there are multiple candidates in an event, the combination with the smallest$ \chi^2_{4\rm C} $ is retained for further analysis, where the$ \chi^2_{4\rm C} $ is from a kinematic fit [40] enforcing four-momentum conservation (4C).To suppress potential contamination from mis-PID in processes with the same charged multiplicity, we perform an additional 4C kinematic fit with the hypotheses of
$ e^+e^-\to $ $ K^{+}K^{-}\pi^{+}\pi^{-} $ ,$ K^{+}K^{-}K^{+}K^{-} $ and$ \pi^{+}\pi^{-}\pi^{+}\pi^{-} $ . By requiring that the$ \chi^2_{4\rm C} $ of the$ K^+K^+e^-e^- $ mass assignment is the smallest one among all possible assignments, such backgrounds can be effectively suppressed and only events with$ \chi^2_{K^+K^+e^-e^-,4\rm C}<20 $ are kept for further analysis. This requirement has been optimized with the Punzi significance [41] using the formula$ \frac{\varepsilon}{1.5+\sqrt{b}} $ , where ε is the detection efficiency from the signal MC simulation and b is the number of background events from the inclusive MC sample. Furthermore, the selected$ K^{+}K^{-}e^{+}e^{-} $ and$ \pi^{+}\pi^{-}e^{+}e^{-} $ combinations will be rejected if the$ e^+e^- $ pair comes from gamma conversions in the beam pipe.The
$ J/\psi $ signal region is defined as$ [3.07,\; 3.12] $ GeV/$ c^2 $ , which corresponds to a range of three times of the mass resolution around the$ J/\psi $ mass [35]. The detection efficiency is determined to be 15.68% based on the simulated$ J/\psi \to K^+K^+e^-e^- $ events.With all above selections applied, only one event from data is left in the signal region as shown in Fig. 2, while no event survives from the inclusive MC sample. To investigate the potential contamination from continuum processes, we analyze the off-resonance data samples taken at
$ \sqrt{s}= 3.650 $ GeV,$ 3.682 $ GeV and$ 3.080 $ GeV [38]. After considering the differences in integrated luminosities, cross sections, particle momenta, and center-of-mass energies, no event survives in the signal region,$ N^{\rm{continuum}}=0 $ . Thus, the total number of background events is calculated as$ N^{\rm{bkg}}=N^{\rm{inc}}+N^{\rm{continuum}}=0 $ . -
Systematic uncertainties in the measurement arise from the total number of
$ J/\psi $ events, tracking, PID, kinematic fit, MC modeling, and signal region.● The uncertainty in the total number of
$ J/\psi $ events, determined using inclusive hadronic events, is 0.5% [26].● The uncertainty in the MDC tracking is estimated by comparing efficiency between data and MC control samples of
$ J/\psi \to e^+ e^- (\gamma_{\mathrm{FSR}}) $ and$ J/\psi \to \pi^0 K^+ K^- $ . The resulting uncertainties are 0.7% per electron and 0.3% per kaon. Consequently, the total systematic uncertainty in the MDC tracking is 2.0%.● The uncertainty in the PID is evaluated using the same high-purity control samples of electrons and kaons as mentioned above. The uncertainties are found to be 0.5% per electron and 0.6% per kaon, resulting in a total PID uncertainty of 2.2%.
● The systematic uncertainty from the 4C kinematic fit and the
$ \chi^2_{4\rm C} < 20 $ requirement is studied using a control sample of$ J/\psi \to p \bar{p} \pi^+ \pi^- $ . The relative difference in efficiencies between data and MC simulation is found to be 6.2%, which is taken as the systematic uncertainty from this source.● The systematic uncertainty associated with the signal region is also calculated using the control sample of
$ J/\psi \to p \bar{p} \pi^+ \pi^- $ . The relative difference in efficiencies, 0.4%, is taken as the uncertainty.● The uncertainty due to the signal MC sample size is 0.9%.
● To obtain the systematic uncertainty from the MC modeling, we compare the difference of the detection efficiency between the PHSP generator and two alternative MC models, i.e.,
$ J/\psi\to\nu_M K^+e^-\to (K^+e^-)K^+e^- $ and$ J/\psi\to\nu_M\nu_M\to(K^+e^-)(K^+e^-) $ , where the Majorana neutrino$ \nu_M $ decays into$ K^+e^- $ with different Majorana masses and corresponding detection efficiencies. Possible Feynman diagrams are shown in Fig. 3. The Majorana neutrino mass is varied from the$ Ke $ mass threshold to the maximum phase space allowed by the$ J/\psi $ decay, approximately from 0.575 GeV/$ c^2 $ to 1.525 GeV/$ c^2 $ , with an interval of 0.025 GeV/$ c^2 $ . The largest efficiency difference, 1.0%, is taken as the systematic uncertainty.Figure 3. A possible Feynman diagram for the decay process of
$ J/\psi\to\nu_M\nu_M\to(K^+e^-)(K^+e^-) $ (left), and for the decay process of$ J/\psi\to\nu_M K^+e^-\to (K^+e^-)K^+e^- $ (right).All the uncertainties are summarized in Table 1, and the total systematic uncertainty (
$ \Delta^{\rm{sys}} $ ) is obtained by adding the individual effects in quadrature.Source Uncertainty (%) $ N_{J/\psi} $ 0.5 MDC tracking 2.0 PID 2.2 4C kinematic fit and $ \chi^2 $ cut6.2 signal region 0.4 MC statistics 0.9 MC modeling 1.0 Total $ \Delta^{\rm{sys}} $ 7.0 Table 1. Summary of the systematic uncertainties.
-
One candidate event (
$ N^{\mathrm{obs}} $ ) is observed in the signal region, while zero background event ($ N^{\mathrm{bkg}} $ ) is expected. The efficiency-corrected upper limit (UL) on the signal yield$ N^{\rm{up}}_{K^+K^+e^-e^-} $ is estimated to be$ 23.2 $ at the 90% CL using a frequentist method [42]. This estimation employs an unbounded profile likelihood approach to account for systematic uncertainties ($ \Delta^{sys}=7.0 $ %), where both observed events ($ N^{\mathrm{obs}} = 1 $ ) and background events ($ N^{bkg}=0 $ ) modeled by a Poisson distribution, the detection efficiency ($ \varepsilon_{K^+K^+e^-e^-}=15.83 $ %) following a Gaussian distribution and the systematic uncertainty is taken as its standard deviation. The UL of the branching fraction (BF) of$ {J/\psi}\to K^{+}K^{+}e^{-}e^{-} $ at the 90% CL is determined by$ \begin{aligned}[b] \mathcal{B}(J/\psi \to K^+K^+e^-e^-) &< \mathcal{B}^{\rm{UL}}(J/\psi \to K^+K^+e^-e^-) \\ &= \frac{N^{\mathrm{up}}_{K^+K^+e^-e^-}}{N^{\mathrm{tot}}_{J/\psi}} = 2.1 \times 10^{-9}, \end{aligned} $
(1) where
$ N_{J/\psi}^{{\rm{tot}}} = (10087\pm 44)\times10^{6} $ is the total number of$ J/\psi $ events in data. -
In addition to the overall UL, assuming the LNV signal from Majorana neutrino, we also investigate the ULs on the branching fractions of
$ J/\psi \to 2\nu_M\to (K^+e^-)(K^+e^-) $ with different Majorana masses. Full kinematic ranges of the$ Ke $ invariant mass ($ m_{Ke} $ ) from the threshold to the largest phase space allowed in the$ J/\psi $ decay,$ m_K+m_e \leq m_{\nu_M} \leq \frac{m_{J/\psi}}{2}, $
(2) is considered. We take the Majorana neutrino mass
$ m_{\nu_M} $ from 0.575 GeV/$ c^2 $ to 1.500 GeV/$ c^2 $ and divide this mass range into uniform intervals, with a width of 0.025 GeV/$ c^2 $ . In each interval, the Majorana mass is taken as the center value. Due to extremely small phase space, the first three (0.500 GeV/$ c^2 $ , 0.525GeV/$ c^2 $ , and 0.550 GeV/$ c^2 $ ) are not selected as mass possibilities.The values of the mass dependent detection efficiencies, the signal event counts, the background event counts, the total systematic uncertainties, and the ULs on the BFs at the 90% CL are reported in Table 2, for each
$ m_{Ke} $ bin. The mass dependent ULs are calculated to be from [2.04, 0.92]$ \times10^{-9} $ with the same method, and the results are shown in Fig. 4.$ m_{\nu_M} $ (GeV/$ c^2 $ )$ \varepsilon_{2\nu_M} $ (%)$ N^{\rm{sig}}_{i} $ $ N^{\rm{bkg}}_{i} $ $ \Delta^{\rm{sys}}_{i} $ (%)$ \mathcal{B}^\mathrm{UL} $ ($ 10^{-9} $ )0.575 8.97 0 0 26.8 $ 2.04 $ 0.600 12.34 0 0 25.8 $ 1.48 $ 0.625 13.69 0 0 26.4 $ 1.34 $ 0.650 13.60 0 0 22.9 $ 1.35 $ 0.675 13.79 0 0 19.3 $ 1.33 $ 0.700 14.15 0 0 16.9 $ 1.29 $ 0.725 14.56 0 0 15.5 $ 1.26 $ 0.750 15.32 0 0 14.1 $ 1.19 $ 0.775 15.87 0 0 12.1 $ 1.15 $ 0.800 16.14 0 0 11.3 $ 1.13 $ 0.825 16.18 0 0 8.4 $ 1.13 $ 0.850 16.67 0 0 11.2 $ 1.10 $ 0.875 16.57 0 0 9.0 $ 1.11 $ 0.900 16.26 0 0 7.0 $ 1.11 $ 0.925 16.55 0 0 7.0 $ 1.11 $ 0.950 16.55 0 0 7.0 $ 1.09 $ 0.975 16.72 0 0 7.0 $ 1.09 $ 1.000 16.74 0 0 7.0 $ 1.08 $ 1.025 16.94 0 0 7.1 $ 1.08 $ 1.050 16.96 0 0 7.0 $ 1.07 $ 1.075 17.16 0 0 7.2 $ 1.08 $ 1.100 16.96 0 0 7.5 $ 1.08 $ 1.125 17.01 0 0 7.8 $ 1.06 $ 1.150 17.28 0 0 7.0 $ 1.05 $ 1.175 17.36 0 0 7.4 $ 1.05 $ 1.200 17.47 1 0 7.2 $ 1.93 $ 1.225 17.94 0 0 7.0 $ 1.02 $ 1.250 17.16 0 0 7.0 $ 1.01 $ 1.275 18.30 0 0 7.5 $ 1.00 $ 1.300 18.55 0 0 7.3 $ 0.99 $ 1.325 18.51 0 0 7.6 $ 0.99 $ 1.350 18.75 0 0 8.5 $ 0.98 $ 1.375 18.80 0 0 8.9 $ 0.97 $ 1.400 19.12 0 0 9.4 $ 0.96 $ 1.425 18.82 0 0 8.7 $ 0.97 $ 1.450 19.31 0 0 11.5 $ 0.95 $ 1.475 19.38 0 0 11.1 $ 0.95 $ 1.500 20.00 0 0 12.8 $ 0.92 $ 1.525 20.28 0 0 16.2 $ 0.90 $ Table 2. The mass dependent detection efficiencies(
$ \epsilon_{i} $ ), the signal event counts($ N^{\rm{sig}}_{i} $ ), the background event counts($ N^{\rm{bkg}}_{i} $ ), the total systematic uncertainties ($ \Delta^\mathrm{sys}_\mathrm{i} $ ), and the ULs on the BFs at the 90% CL($ \mathcal{B}^\mathrm{UL} $ ) for different Majorana mass hypotheses from 0.575 GeV/$ c^2 $ to 1.525 GeV/$ c^2 $ with an interval of 0.025 GeV/$ c^2 $ . -
In summary, by using
$ (10087\pm 44)\times10^{6} $ $ J/\psi $ events collected by the BESIII detector, we search for the LNV decay$ J/\psi\to K^+K^+e^-e^- +c.c. $ for the first time. One event survives in the signal region and the UL on the decay branching fraction of$ J/\psi\to K^+K^+e^-e^- +c.c. $ is determined to be$ 2.1 \times 10^{-9} $ at the 90% CL. The Majorana mass dependent ULs on the branching fractions at the 90% CL are estimated, varying within the range of [0.92, 2.04]$ \times10^{-9} $ . By utilizing these results, constraints on the mixed matrix element$ |V_{e\nu_{M}}|^2 $ between the Majorana neutrino and electron can be derived within the theoretical models. -
The BESIII Collaboration thanks the staff of BEPCII (https://cstr.cn/31109.02.BEPC) and the IHEP computing center for their strong support.
Search for the lepton number violating process ${{\boldsymbol J}/{\boldsymbol\psi }{\bf\to}{\boldsymbol K}^{\bf +}{\boldsymbol K}^{\bf +}{\boldsymbol e}^{\bf -}{\boldsymbol e}^{\bf -}{\bf +}{\boldsymbol {c.c.}}}$
- Received Date: 2025-07-09
- Available Online: 2026-02-01
Abstract: Based on