Towards the meV limit of the effective neutrino mass in neutrinoless double-beta decays

Neutrino oscillation experiments have firmly established that neutrinos are massive particles and lepton flavors are significantly mixed [1]. As expected in a class of seesaw models for neutrino mass generation [2], massive neutrinos are Majorana particles and they may provide a natural and elegant explanation for the observed matter-antimatter asymmetry in our Universe [3]. If this is indeed the case, it will be a great challenge to determine two associated Majorana-type CP-violating phases, which are measurable only in the lepton-number-violating processes. Currently, the experimental search for neutrinoless double-beta ($ 0\nu\beta\beta $) decays $ {^{A}_Z}N \to {^A_{Z+2}}N + 2e^- $ of some heavy nuclei $ {^{A}_Z}N $, which possess an even atomic number $ Z $ and an even mass number $ A $, is the most promising way to demonstrate the Majorana nature of massive neutrinos and to prove the existence of lepton number violation in nature [4]. Assuming that three light Majorana neutrinos are responsible for the $ 0\nu\beta\beta $ decays of an even-even nuclear isotope, we can find the half-life [5]

where $ G_{0\nu} $ denotes the relevant phase-space factor, $ {\cal M}_{0\nu} $ is the nuclear matrix element (NME), and $ m_e = 0.511\; {\rm MeV} $ is the electron mass. In the standard parametrization of lepton flavor mixing matrix, the effective neutrino mass $ |m_{\beta \beta}| $ for $ 0\nu\beta\beta $ decays appearing in Eq. (1) reads

where $ m_i $ (for $ i = 1, 2, 3 $) stand for the absolute masses of three ordinary neutrinos. Out of three neutrino mixing angles only two $ \{\theta_{12}, \theta_{13}\} $ are involved in the effective neutrino mass in Eq. (2), where $ \{\rho, \sigma\} $ are two Majorana CP phases. The other neutrino mixing angle $ \theta_{23} $ and the Dirac-type CP-violating phase $ \delta $ are irrelevant for $ 0\nu\beta\beta $ decays.

The main purpose of the present study is to explore the physics potential of pinning down the fundamental parameters in future $ 0\nu\beta\beta $ decay experiments, in particular the absolute neutrino masses and two Majorana phases that are not accessible at all in neutrino oscillation experiments. Our motivation is three-fold:

● Neutrino oscillation experiments have measured with reasonably good precisions the relevant two neutrino mixing angles $ \{\theta_{12}, \theta_{13}\} $, and two independent neutrino mass-squared differences $ \Delta m^2_{21} \equiv m^2_2 - m^2_1 $ and $ |\Delta m^2_{31}| \equiv |m^2_3 - m^2_1| $ [1]. In the near future, the JUNO experiment [6,7] will be able to offer an unambiguous answer to whether neutrino mass ordering is normal $ m_1 < m_2 < m_3 $ (NO) or inverted $ m_3 < m_1 < m_2 $ (IO), and to improve the precisions of three parameters $ \{\sin^2 \theta_{12}, \Delta m^2_{21}, \Delta m^2_{31}\} $ to the level below one percent. In addition, the ultimate precision on $ \sin^2\theta_{13} $ from the Daya Bay experiment will be 3% [8,9]. According to neutrino oscillation data, at least two neutrino masses should be above the meV level, e.g., $ m_3 > m_2 = \sqrt{m^2_1 + \Delta m^2_{21}} \geqslant \sqrt{\Delta m^2_{21}} \approx 8.6\; {\rm {meV}} $ (for NO). Given oscillation parameters, the observation of $ 0\nu\beta\beta $ decays will be extremely important in the determination of the lightest neutrino mass $ m_{\rm 1} $ (for NO) or $ m_3 $ (for IO) and two Majorana CP phases $ \{\rho, \sigma\} $. It should be noticed that the latest global-fit analysis of neutrino oscillation data already hints at a $ 2\sigma $ preference for the NO [10], when the atomspheric neutrino data from Super-Kamiokande are not taken into account. The hints for the NO have also been found by other independent analyses [11,12]. The neutrino mass ordering will be unambiguously clarified by future neutrino oscillation experiments with a much higher statistical significance [13].

● The upper bound on the absolute scale of neutrino masses extracted from the tritium beta decays is $ m_\beta \leqslant 2.3\; {\rm {eV}} $ (Mainz [14]) and $ m_\beta \leqslant 2.2\; {\rm {eV}} $ (Troitsk [15]) at the 95% confidence level, where the effective neutrino mass $ m_\beta $ for beta decays is defined as $ m_\beta \equiv \left( m^2_1 |U_{e1}|^2 +\right. $$\left. m^2_2 |U_{e2}|^2 + m^2_3 |U_{e3}|^2\right)^{1/2} $ with the moduli of the matrix elements of lepton flavor mixing matrix being $ |U_{e1}| = $ $ \cos \theta_{13} \cos \theta_{12} $, $ |U_{e2}| = \cos \theta_{13} \sin \theta_{12} $ and $ |U_{e3}| = \sin \theta_{13} $. The first measurement of tritium beta decays by the KATRIN experiment has improved the limit to $ m_\beta \leqslant 1.1\; {\rm {eV}} $ at the 90% confidence level [16,17]. In future, the full operation of the KATRIN [18,19] and Project 8 [20] experiments will be capable of bringing the upper limit down to $ m_\beta \leqslant 200\; {\rm {meV}} $ and $ m_\beta \leqslant 40\; {\rm {meV}} $, respectively. On the other hand, the cosmological observations of cosmic microwave background by the Planck satellite gives the most restrictive bound on the sum of three neutrino masses $ \Sigma \equiv m_1 + m_2 + m_3 < 120\; {\rm {meV}} $ [21,22]. However, there is still a long way to go until the neutrino mass region of a few meV is accessed.

● Future large and ultra-low background liquid scintillator (LS) detectors, such as JUNO, have great potential of searching for $ 0\nu\beta\beta $ decays, by dissolving the $ 0\nu\beta\beta $-decaying isotope ¹³⁰Te or ¹³⁶Xe into LS. This concept has been discussed in Ref. [23], where the Xe-loaded LS is taken as a target. It has been demonstrated that a sensitivity (at the 90% confidence level) to $ T^{1/2}_{0\nu} $ of $ 1.8 \times 10^{28}\; {\rm {yr}} $ is achievable with 50 tons of fiducial ¹³⁶Xe and 5 years of exposure, while the corresponding sensitivity to the effective neutrino mass $ |m_{\beta\beta}| $ could reach $ (5\cdots 12)\; {\rm {meV}} $ depending on the NME value. It has also been pointed out that ¹³⁰Te may be an advantageous candidate due to its high natural abundance. If the nuclear isotope ¹³⁰Te with a maximum fraction of 4% is loaded, a total target mass of $ 400 $ tons could be obtained, leading to an improvement on the sensitivity by a factor of $ (400/50)^{1/4} \approx 1.68 $, namely, $ |m_{\beta\beta}| \approx (2.3\cdots 6.0)\; {\rm {meV}} $ for the same background index.^① Moreover, if the nominal value of the background index in Ref. [23] is further significantly reduced, e.g., by two orders of magnitude, the ultimate sensitivity will hopefully be close to $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $. See, e.g., Ref. [24], for an earlier discussion about the sensitivity of future $ 0\nu\beta\beta $ decay experiments. While it is not clear that $ 1\; {\rm {meV}} $ sensitivity could be practically achieved, it at least points towards a direction for experimental development that could extend sensitivity beyond the currently planned next generation of instruments into the NO region.

In the literature, it has been noticed [25-29] that $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ could be set as a target value and useful information on the absolute neutrino masses and the Majorana CP phases can be obtained. We emphasize that with the sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ in future $ 0\nu\beta\beta $ decay experiments, it is likely to observe the signal events. Even if this is not the case, a stringent constraint on or even a precise determination of the absolute neutrino masses can be achieved. The present study differs from previous works in two aspects. First, we concentrate on the effective neutrino mass $ |m_{\beta\beta}| $ in the NO case and update its value with both the latest global-fit results of all the relevant neutrino oscillation parameters and the future measurements from neutrino oscillation experiments. These inputs help better define how much the effective mass $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ can tell us the information about the lightest neutrino mass $ m_1 $ and the Majorana phases $ \{\rho, \sigma\} $. Second, we further explore the implications for the neutrino mass spectrum, the effective mass $ m_\beta $ for beta decays and the sum of three neutrino masses $ \Sigma $, and stress that the determination of the lightest neutrino mass from the $ 0\nu\beta\beta $ decay experiment with a sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ sets up a challenging goal for future beta-decay experiments and cosmological observations. As it is unclear whether the sensitivity at the level of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ can be ultimately accomplished in a realistic $ 0\nu\beta\beta $ decay experiment, we take this sensitivity as a working scenario and investigate its phenomenological implications. Our main conclusions are summarized as follows.

● With the help of precision measurements of two neutrino mass-squared differences from oscillation experiments, we can completely determine the neutrino mass spectrum, namely, $ m_1 \in [0.7, 8]\; {\rm {meV}} $, $ m_2 \in [8.6, 11.7]\; {\rm {meV}} $ and $ m_3 \in [50.3, 50.9]\; {\rm {meV}} $.

● One of the two Majorana CP-violating phases is limited to $ \rho \in [130^\circ, 230^\circ] $, which cannot be obtained from any other realistic experiments. Unfortunately, the other Majorana CP phase $ \sigma $ is essentially unconstrained, and the observations of different lepton-number-violating processes are necessary to fix both phases in the far future.

● The sum of three neutrino masses is found to be $ \Sigma \in [59.2, 72.6]\; {\rm {meV}} $, while the effective neutrino mass for beta decays turns out to be $ m_\beta \in [8.9, 12.6]\; {\rm {meV}} $. It will be a challenging task for future cosmological observations and beta-decay experiments to achieve these values.

The remaining part of this paper is organized as follows. In Sec. 2, the three-dimensional description of the effective neutrino mass $ |m_{\beta\beta}| $ as a function of the lightest neutrino mass $ m_1 $ and the Majorana CP phase $ \rho $ is given, where the future precision on neutrino oscillation parameters is implemented and the latest global-fit results from Ref. [10] are also considered for comparison. Furthermore, the implications for the neutrino mass spectrum, the effective neutrino mass $ m_\beta $ in beta decays and the sum of three neutrino masses from cosmological observations are discussed in Sec. 3. Finally, we give some further remarks and summarize our main conclusions in Sec. 4.

The conventional way of graphically showing the possible range of $ |m_{\beta\beta}| $ is to plot it as a function of the lightest neutrino mass ($ m_1 $ for the NO case or $ m_3 $ for the IO case) by varying $ \rho $ and $ \sigma $ in the whole range of $ [0,360^\circ) $, as first suggested in Ref. [30]. In Fig. 1, the allowed range of $ |m_{\beta\beta}| $ in the NO case is shown as the gray region. The boundaries of the allowed range are denoted by the dashed curves, which are obtained by using the best-fit values of $ \{\theta_{12}, \theta_{13}\} $ and $ \{\Delta m^2_{21}, \Delta m^2_{31}\} $ from the latest global-fit analysis of neutrino oscillation data in Ref. [10]. In the left panel, the red bands along the dashed curves are caused by the $ 1\sigma $ uncertainties of the oscillation parameters from the global-fit analysis, while those in the right panel are due to the $ 1\sigma $ uncertainties after the JUNO measurements [7]. In each panel, the horizontal dashed line corresponds to $ |m_{\beta\beta}| = 1\; {\rm {meV}} $ and the allowed range of $ m_1 $ is indicated by two vertical dashed lines. More explicitly, we quote the best-fit values and current uncertainties of the relevant oscillation parameters from Ref. [10] in the NO case as below

Figure 1.  (color online) The effective neutrino mass $ |m_{\beta\beta}| $ is shown as a function of the lightest neutrino mass $ m_1 $ in the NO case, where the gray region is allowed and the dashed curves refer to its boundaries. In the left panel, the boundaries are obtained by using the best-fit values of $ \{\theta_{12}, \theta_{13}\} $ and $ \{\Delta m^2_{21}, \Delta m^2_{31}\} $ from the latest global-fit analysis of neutrino oscillation parameters in Ref. [10] and the red region is caused by the $ 1\sigma $ uncertainties of these parameters. In the right panel, the future precisions on those oscillation parameters after the JUNO measurements [7] are implemented. In each panel, the horizontal dashed line corresponds to $ |m_{\beta\beta}| = 1\; {\rm {meV}} $ and the allowed range of $ m_1 $ is indicated by two vertical dashed lines.

For the future measurements of these parameters, we assume that the best-fit values are the same, but the precision on $ \sin^2 \theta_{12} $, $ \Delta m^2_{21} $ and $ \Delta m^2_{31} $ will be improved after the JUNO experiment to 0.54%, 0.24% and 0.27%, respectively. It is worth mentioning that the ultimate precision of 3% on $ \sin^2\theta_{13} $ from the Daya Bay experiment will be adopted, which is comparable to that in Eq. (3).

In order to clarify the dependence of $ |m_{\beta\beta}| $ on the oscillation parameters, we consider its upper ("U") and lower ("L") boundaries that are derived by varying the Majorana CP phase $ \sigma $, namely,

where the sign "+" (or "-") corresponds to "U" (or "L"), and $ \overline{m}_{12} \equiv m_1 \cos^2 \theta_{13} \cos^2 \theta_{12} {\rm e}^{{\rm i}\rho} + m_2 \cos^2 \theta_{13} \sin^2 \theta_{12} $ is the sum of the first two terms in $ m_{\beta\beta} $ defined in Eq. (2). Notice that the phase $ \sigma $ has been properly chosen to match exactly (or differ by $ \pm \pi $ from) the phase of $ \overline{m}_{12} $ to draw the upper (lower) boundary. Some comments on the upper and lower boundaries of $ |m_{\beta\beta}| $ in Fig. 1 are in order.

In the NO case, given two neutrino mass-squared differences $ \Delta m^2_{21} \!\simeq \! 7.39\!\times \!10^{-5}\; {\rm {eV}}^2 $ and $ \Delta m^2_{31} \!\simeq \!2.523\!\times \!10^{-3}\; {\rm {eV}}^2 $, we have $ m_3 = \sqrt{m^2_1 + \Delta m^2_{31}} $ and $ m_2 = \sqrt{m^2_1 + \Delta m^2_{21}} $. The upper boundary is determined by

implying $ |m_{\beta\beta}|_{\rm U} \to \sin^2 \theta_{12} \cos^2\theta_{13} \sqrt{\Delta m^2_{21}} + \sin^2\theta_{13} \sqrt{\Delta m^2_{31}} $ in the limit of $ m_1 \to 0 $. In this limit, the first term $ \sin^2 \theta_{12} \cos^2 \theta_{13} \sqrt{\Delta m^2_{21}} \simeq 2.6\; {\rm {meV}} $ is about twice larger than the second term $ \sin^2\theta_{13} \sqrt{\Delta m^2_{31}}\simeq 1.1\; {\rm {meV}} $, where $ \sin^2 \theta_{12} \simeq 0.310 $ and $ \sin^2 \theta_{13} \simeq 0.02241 $ have been used. Hence the uncertainties from $ \sin^2 \theta_{12} $ and $ \Delta m^2_{21} $ dominate over those from $ \sin^2 \theta_{13} $ and $ \Delta m^2_{31} $. This can be clearly seen by comparing the uncertainty of the upper boundary in the left panel of Fig. 1 with that in the right panel, as the JUNO measurements improve the precisions on $ \sin^2 \theta_{12} $ and $ \Delta m^2_{21} $. When the lightest neutrino mass $ m_1 $ increases, in particular for $ m_1 \gtrsim \sqrt{\Delta m^2_{21}} \simeq 8.6\; {\rm {meV}} $, we obtain $ m_1 \simeq m_2 $ and the first term on the right-hand side of Eq. (5) is approximately given by $ m_1 \cos^2 \theta_{13} $, indicating that the dependence on the $ \sin^2 \theta_{12} $ and $ \Delta m^2_{21} $ becomes very weak. On the other hand, the lower boundary in the limit of $ m_1 \to 0 $ reads

which is estimated to be $ |m_{\beta\beta}|_{\rm L} \simeq 1.5\; {\rm {meV}} $ for the best-fit values of neutrino mixing angles and mass-squared differences in Eq. (3). Since the first term on the right-hand side of Eq. (6) is about twice larger than the second term, the dominant uncertainty on the lower boundary comes from $ \sin^2 \theta_{12} $ and $ \Delta m^2_{21} $ [31]. For $ m_1 \gtrsim \sqrt{\Delta m^2_{21}} \simeq 8.6\; {\rm {meV}} $ and thus $ m_1\simeq m_2 $, we obtain $|m_{\beta\beta}|_{\rm L} \simeq \cos^2\theta_{13} \cos2\theta_{12} $$ \sqrt{\Delta m^2_{21} + m^2_1} - \sin^2\theta_{13} \sqrt{\Delta m^2_{31} + m^2_1} $, so the uncertainty from $ \sin^2 \theta_{12} $ continues to be dominant as clearly shown in Fig. 1. In the region of $ 2\; {\rm {meV}} \lesssim m_1 \lesssim 7\; {\rm {meV}} $, the destructive cancellation appears in $ |m_{\beta\beta}| $ due to the Majorana CP phases, leading to a "well"-like structure, which will be examined more carefully in the next subsection.

We have observed that the boundaries of the allowed range of $ |m_{\beta\beta}| $ in the NO case depend crucially on the precision of $ \sin^2 \theta_{12} $ and $ \Delta m^2_{21} $. After the JUNO measurements, as illustrated in the right panel of Fig. 1, the uncertainties from neutrino mixing angles and mass-squared differences could be safely ignored. It will be helpful to quantitatively assess the impact of the uncertainties of neutrino oscillation parameters. If the $ 1\sigma $ uncertainties of neutrino mixing parameters from current neutrino oscillation data are taken, the value of $ |m_{\beta\beta}| = 1\; {\rm {meV}} $ will be reached in the range of $ m_1 \in [0.42, 9.6]\; {\rm {meV}} $. With the improvement by the JUNO measurements, we obtain $ m_1 \in [0.62, 9.0]\; {\rm {meV}} $ for $ |m_{\beta\beta}| = 1\; {\rm {meV}} $, where the lower bound of $ m_1 $ has been modified significantly.

As we have mentioned, in the so-called Vissani graph [30] in Fig. 1, there is a “well”-like structure of $ |m_{\beta\beta}| $ in the region of $ 2\; {\rm {meV}} \lesssim m_1 \lesssim 7\; {\rm {meV}} $ for the NO case. Inside the well, $ |m_{\beta\beta}| $ takes tiny values, where a significant cancellation among three components of $ m_{\beta\beta} $ occurs. The bottom of the well signifies the extreme case of $ |m_{\beta\beta}| \to 0 $, where a complete cancellation takes place. In Ref. [25], the three-dimensional graph, where $ |m_{\beta\beta}| $ is plotted against both the lightest neutrino mass $ m_1 $ and the Majorana CP phase $ \rho $, has been suggested and shown to be extremely useful in revealing the fine structure inside the well [26,27].

In Fig. 2, we reproduce the three-dimensional graph of $ |m_{\beta\beta}| $ in the NO case from Refs. [26,27] and zoom in the region of $ 10^{-2}\; {\rm {meV}} \leqslant |m_{\beta\beta}| \leqslant 1.1\; {\rm {meV}} $, where the best-fit values of neutrino mixing angles and mass-squared differences from Ref. [10] have been used. As we have explained, the upper (or lower) boundary of the allowed range of $ |m_{\beta\beta}| $ can be obtained by properly choosing the Majorana CP phase $ \sigma $ such that it matches (or differs by $ \pm \pi $ from) the phase of $ \overline{m}_{12} $. More explicitly, the values of $ \sigma $ are determined by

Figure 2.  (color online) The zoomed-in three-dimensional graph of $ |m_{\beta\beta}| $ as a function $ m_1 $ and $ \rho $ for the NO case, where the best-fit values $ \Delta m^2_{21} = 7.39 \times 10^{-5}\; {\rm {eV}}^2 $, $ \Delta m^2_{31} = 2.523 \times 10^{-3}\; {\rm {eV}}^2 $, $ \sin^2\theta_{12} = 0.310 $ and $ \sin^2\theta_{13} = 0.02241 $ have been taken [10]. The blue filled volume stands for the allowed range of $ |m_{\beta\beta}| $, while three contour surfaces corresponding to $ |m_{\beta\beta}| = 1.1\; {\rm {meV}} $, $ 0.3\; {\rm {meV}} $ and $ 0.1\; {\rm {meV}} $ are shown and the “bullet”-like region is hollowed out.

With the help of Eq. (4), one finds that the bottom of the well (i.e., $ |m_{\beta\beta}|_{\rm L} = 0 $) would be reached for $ \left| \overline{m}_{12} \right| = m_3 \sin^2\theta_{13} $. Then it is straightforward to verify that $ |m_{\beta\beta}| = 0 $ holds only in the narrow region of $ 2 \; {\rm {meV}} \lesssim m_1 \lesssim 7 \; {\rm {meV}} $ and $ 155^\circ \lesssim \rho \lesssim 205^\circ $. The other Majorana CP phase $ \sigma $ is fixed via Eq. (7) by choosing the minus sign.

As shown in Fig. 2, three contour surfaces corresponding to $ |m_{\beta\beta}| = 1.1\; {\rm {meV}} $, $ 0.3\; {\rm {meV}} $ and $ 0.1\; {\rm {meV}} $ appear above the “bullet”-like structure. The surface of this bullet is described by

whose maximum $ |m_{\beta\beta}|_{*} = m_3 \sin^2 \theta_{13} $ appears at $ \overline{m}_{12} = 0 $ that in turn requires $ \rho = 180^\circ $ and $ m_1/m_2 = \tan^2\theta_{12} $ [26]. Numerically, for the best-fit values of oscillation parameters in Eq. (3), the tip of the bullet is located at $ \left(m_1, \rho, |m_{\beta\beta}|_{*}\right) \simeq \left(4 \; {\rm {meV}}, 180^\circ, 1.1\; {\rm {meV}}\right) $. Note that the region covered by the bullet surface is actually hollowed out, so the allowed parameter space for small values of $ |m_{\beta\beta}| $ decreases rapidly, which becomes clear by comparing among three contour surfaces in Fig. 2. For the following two reasons, $ |m_{\beta\beta}|_{*} \simeq 1.1\; {\rm {meV}} $ can be taken as a threshold value in some sense. On the one hand, as one has already seen, the parameter space for $ |m_{\beta\beta}| \lesssim |m_{\beta\beta}|_{*} $ to hold is very small as compared with the whole parameter space. On the other hand, although it is obviously challenging for the future $ 0\nu\beta\beta $ decay experiments to reach the sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $, this goal is hopefully achievable by further optimizing and improving the experimental setup considered in Ref. [23]. If such a sensitivity is ultimately reached, one can probe the absolute neutrino masses to an unprecedented precision at the meV level and draw a restrictive constraint on the Majorana CP phases.

Now let us assume that the sensitivity of $ |m_{\beta\beta}| = 1\; {\rm {meV}} $ or even below can be reached in the future $ 0\nu\beta\beta $ decay experiments. In this case, the allowed parameter space of $ (m_1, \rho, \sigma) $ is given in Fig. 3. For comparison, we have shown the results for $ |m_{\beta\beta}| = 1\; {\rm {meV}} $, $ 0.3\; {\rm {meV}} $ and $ 0\; {\rm {meV}} $ in the upper panel. In the lower panel, the projections of the parameter space into the $ (m_1, \rho) $-, $ (m_1, \sigma) $- and $ (\rho, \sigma) $-plane have also been presented. Some comments on the numerical results are helpful.

Figure 3.  (color online) The allowed parameter space of $ m_1 $, $ \rho $ and $ \sigma $ for $ |m_{\beta\beta}| \leq 1\; {\rm {meV}} $, $ |m_{\beta\beta}| \leq 0.3\; {\rm {meV}} $ and $ |m_{\beta\beta}| = 0 $, where the same best-fit values of $ \Delta m^2_{21} $, $ \Delta m^2_{31} $, $ \sin^2\theta_{12} $ and $ \sin^2\theta_{13} $ as in Fig. 2 have been used. In the lower panel, the projections of the allowed parameter space into the $ (m_1, \rho) $-, $ (m_1, \sigma) $- and $ (\rho, \sigma) $-plane have also been shown.

● In our calculations, the best-fit values of neutrino mixing angles and mass-squared differences in Eq. (3) have been adopted. After taking account of the JUNO measurements, the uncertainties of those parameters will be negligible. In fact, we have also numerically checked that this is indeed true.

● If the sensitivity of $ |m_{\beta\beta}| = 1\; {\rm {meV}} $ is eventually achieved, the allowed ranges of three fundamental parameters $ m_1 $, $ \rho $ and $ \sigma $ can be read off from the red shaded areas in the two-dimensional plots in the lower panel of Fig. 3. More explicitly, we have $ m_1 \in [0.7, 8]\; {\rm {meV}} $, $ \rho \in [130^\circ, 230^\circ] $ and $ \sigma \in [0,360^\circ) $. Note that there is a tiny blank region in the center of the red area in the $ (m_1, \rho) $-plane. This is due to the fact that the tip of the bullet is located at $ |m_{\beta\beta}| = |m_{\beta\beta}|_* \simeq 1.1\; {\rm {meV}} $, which is slightly above the contour surface of $ |m_{\beta\beta}| = 1\; {\rm {meV}} $.

● When the sensitivity is further improved to $ |m_{\beta\beta}| = 0.3\; {\rm {meV}} $, only the green bands in the two-dimensional plots are allowed. As a consequence, the parameter space turns out to be more strictly constrained, namely, $ m_1 \in [1, 7]\; {\rm {meV}} $ and $ \rho \in [150^\circ, 210^\circ] $. However, the whole range of $ \sigma $ is still allowed. In the extreme case of $ |m_{\beta\beta}| = 0\; {\rm {meV}} $, the parameter space is represented by the blue curves. Given the lightest neutrino mass $ m_1 $, one can completely pin down the values of the Majorana CP phases $ \rho $ and $ \sigma $. Even in this case, depending on the exact value of $ m_1 $, any value of $ \sigma $ within $ [0,360^\circ) $ can be taken.

To quantitatively explain the constraining power of the $ 0\nu\beta\beta $ decay experiments, we introduce the probability $ P(|m_{\beta\beta}| < |m_{\beta\beta}|_*) $ as a function of the lightest neutrino mass $ m_1 $. For a given value of $ m_1 $, this probability is calculated as the ratio of the required ranges of $ \rho $ and $ \sigma $ (for $ |m_{\beta\beta}| < |m_{\beta\beta}|_* \simeq 1.1\; {\rm {meV}} $ to be satisfied) to the whole range of $ 360^\circ \times 360^\circ $. In Fig. 4, $ P(|m_{\beta\beta}| < |m_{\beta\beta}|_*) $ has been plotted as the red solid curve. One can observe that this fraction is always below 10% for any viable value of $ m_1 $. Based on the above calculations, we conclude that the fulfillment of $ |m_{\beta\beta}| \lesssim |m_{\beta\beta}|_* \simeq 1.1\; {\rm {meV}} $ requires significant cancellation among three terms in $ m_{\beta\beta} $, and the possibility for such a case to come true is really small [26,27] unless there is some underlying symmetry leading to such a cancellation. In other words, future $ 0\nu\beta\beta $ experiments with a sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ will be able to determine the lightest neutrino mass with a high precision and even to probe the Majorana CP phases.

Figure 4.  (color online) The ratio of the region of $ (\rho, \sigma) $ for $ |m_{\beta\beta}| < |m_{\beta\beta}|_{*} $ to hold to the whole region $ 360^\circ \times 360^\circ $ for any given value of $ m_1 $, where the best-fit values of $ \Delta m^2_{21} $, $ \Delta m^2_{31} $, $ \sin^2\theta_{12} $ and $ \sin^2\theta_{13} $ have been input as in Fig. 2.

Since the lightest neutrino mass $ m_1 $ can be constrained to a narrow range of $ [0.7, 8]\; {\rm {meV}} $ for $ |m_{\beta\beta}| \lesssim 1\; {\rm {meV}} $, it is interesting in the first place to see how well the full neutrino mass spectrum can be determined. In Fig. 5, the absolute neutrino masses $ m_i $ (for $ i = 1, 2, 3 $) have been plotted as three red curves against the lightest neutrino mass $ m_1 $ [32]. The requirement for $ |m_{\beta\beta}| \lesssim 1\; {\rm {meV}} $ leads to

Figure 5.  (color online) Illustration for the absolute neutrino masses $ m_i $ (for $ i = 1, 2, 3 $) as implied by the sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ for future $ 0\nu\beta\beta $ decay experiments, where the allowed range of $ m_1 \in [0.7, 8]\; {\rm {meV}}$ is lying between two vertical dashed lines.

where the precisions on neutrino oscillation parameters after the JUNO measurements are used. The allowed range of the lightest neutrino mass $ m_1 \in [0.7, 8]\; {\rm {meV}} $ is lying between two vertical dashed lines in Fig. 5, where one can accordingly find out the allowed ranges of $ m_2 $ and $ m_3 $ on the vertical axis by following the intersecting points between the red curves and the dashed lines. From Eq. (9), one immediately observes that $ m_3 $ is already determined with an excellent precision, while the uncertainties on $ m_2 $ and $ m_1 $ are also acceptable. It is worthwhile to stress that the complete determination of neutrino mass spectrum will be very suggestive for the model building of neutrino masses.

Then we explore the implications for the effective neutrino mass $ m_\beta $ in beta decays and the sum of three neutrino masses $ \Sigma $. For this purpose, we have depicted the allowed regions of three observables $ |m_{\beta\beta}| $ from $ 0\nu\beta\beta $ decays, $ m_\beta $ from beta decays, and $ \Sigma $ from cosmological observations in Fig. 6. For three plots in the left column, the global-fit results of neutrino oscillation parameters in Eq. (3) have been input, while for those in right column we have considered the projected precisions after the JUNO measurements. The dashed curves as the boundaries of the gray shaded areas have been obtained by using the best-fit values, and the red bands are caused by the $ 1\sigma $ uncertainties of the input parameters. The improvements after the JUNO measurements can be observed, but the quantitative changes in the upper bounds on $ m_\beta $ and $ \Sigma $ are at the level of $ 0.1\; {\rm {meV}} $. For $ |m_{\beta\beta}| = 1\; {\rm {meV}} $, which has been shown as the horizontal dashed line in the plots in the first two rows of Fig. 6, one can immediately extract the corresponding ranges of $ m_\beta $ and $ \Sigma $. For instance, assuming the $ 1\sigma $ uncertainties of neutrino oscillation parameters after the JUNO experiment, we can obtain

Figure 6.  (color online) The allowed regions of three observables $ |m_{\beta\beta}| $, $ m_\beta $ and $ \Sigma .$, where the global-fit results of current neutrino oscillation data have been used for three plots in the left column while the JUNO measurements are taken into account for those in the right column.

which is respectively below the forecasted sensitivity $ m_\beta \lesssim 40\; {\rm {meV}} $ of the forthcoming beta-decay experiment [20] and that $ \Sigma \lesssim 80\; {\rm {meV}} $ of future observations of cosmic microwave background [33].

As the Majorana CP phases are only involved in the $ 0\nu\beta\beta $ decays, it is impossible to completely pin down all three unknown parameters $ m_1 $, $ \rho $ and $ \sigma $ from such a single type of observations. If the effective neutrino mass $ m_\beta $ can be precisely measured in beta-decay experiments, one will be able to solve the lightest neutrino mass $ m_1 $ with the help of neutrino oscillation data, namely,

and then to predict $ \Sigma = m_1 + \sqrt{m^2_1 + \Delta m^2_{21}} + $$ \sqrt{m^2_1 + \Delta m^2_{31}} $. Moreover, inserting the determined neutrino masses $ m_1 $, $ m_2 = \sqrt{m^2_1 + \Delta m^2_{21}} $ and $ m_3 = \sqrt{m^2_1 + \Delta m^2_{31}} $ into the effective neutrino mass $ |m_{\beta\beta}| $, we can have a good opportunity to probe the Majorana CP phases $ \rho $ and $ \sigma $ by improving the sensitivity of $ 0\nu\beta\beta $ decay experiments to $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $. But this seems to be not the case.

For this reason, we take a different strategy to pin down neutrino masses in the assumption that both $ |m_{\beta\beta}| $ and $ \Sigma $ can be better measured in the foreseeable future. Given the information on the Majorana CP phases, it is straightforward to establish the following relations

where the flavor mixing matrix elements are given by $ U_{e1} = \cos \theta_{13} \cos \theta_{12} e^{{\rm i}\rho/2} $, $ U_{e2} = \cos \theta_{13} \sin \theta_{12} $ and $ U_{e3} = \sin\theta_{13} e^{{\rm i}\sigma/2} $ in accordance with the parametrization adopted in Eq. (2). After directly solving Eq. (12) for neutrino masses, one arrives at

from which we can easily verify that $ m_1 + m_2 + m_3 = \Sigma $ holds as it should do. Notice that the denominator $ {\rm {Im}}\left[(U^2_{e3} - U^2_{e1})^* \cdot (U^2_{e1} - U^2_{e2})\right] $ is supposed to be nonzero in Eq. (13). If it is zero for some specific values of $ \rho $ and $ \sigma $, then we can extract neutrino masses just from the effective mass $ |m_{\beta\beta}| $ and neutrino oscillation parameters.

As an application of Eq. (13), let us consider the special scenario in which $ |m_{\beta\beta}| = 0 $ is realized. In this particular case, both $ \rho $ and $ \sigma $ will be fixed for the complete cancellation in $ |m_{\beta\beta}| $ to happen [26]. Then, we can find three neutrino masses

where one can observe that neutrino masses are evidently proportional to $ \Sigma $. Certainly, the same results can also be obtained by requiring the real and imaginary parts of $ m_{\beta\beta} $ to be vanishing and then deriving two independent neutrino mass ratios [34,35].

In the present work, we have explained why it is important for the future $ 0\nu\beta\beta $ decay experiments to reach the sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $. With such a high sensitivity, it is likely to observe the signal events of $ 0\nu\beta\beta $ decays. Even if there is a null signal, one can place restrictive constraints on the lightest neutrino mass $0.7\; {\rm {meV}} \leqslant m_1 \leqslant $$ 8\; {\rm {meV}} $ and one of the Majorana CP phases $ 130^\circ \leqslant \rho \leqslant 230^\circ $. Interestingly, these constraints further imply that the whole neutrino mass spectrum is almost fixed, namely, $ m_1 \in [0.7, 8]\; {\rm {meV}} $, $ m_2 \in [8.6, 11.7]\; {\rm {meV}} $ and $ m_3 \in [50.3, $$50.9]\; {\rm {meV}} $, the effective neutrino mass for beta decays should be lying in the range $8.9\; {\rm {meV}} \leqslant $$ m_{\beta} \leqslant 12.6\; {\rm {meV}} $ and the sum of three neutrino masses must be $ 59.2\; {\rm {meV}} \leqslant \Sigma \leqslant 72.6\; {\rm {meV}} $. These limits are lying below the forecasted sensitivities from the next-generation beta-decay experiments [19,20] and the future cosmological observations [33]. In this work we have exclusively focused on the Majorana nature of massive neutrinos. If massive neutrinos were the Dirac particles, the lepton-number-violating $ 0\nu\beta\beta $ decays would certainly not take place. In this case we would have no chance to observe any signal of the $ 0\nu\beta\beta $ decays. But it should be kept in mind that a null result of the $ 0\nu\beta\beta $ decays does not mean that massive neutrinos must be of the Dirac nature, simply because our analysis has clearly shown that there is always a small parameter space for $ |m_{\beta\beta}| \to 0 $ even though no other exotic new physics is introduced. That is why one should do everything in one’s power to go beyond the $ 0\nu\beta\beta $ decays and explore some other presumably promising lepton-number-violating processes, although this is no doubt much more challenging if not impossible.

Among all the current $ 0\nu\beta\beta $ decay experiments in operation [36-42], the KamLAND-Zen collaboration has reported the best sensitivity $ |m_{\beta\beta}| < (61\cdots 165)\; {\rm {meV}} $ depending on the NME for the $ 0\nu\beta\beta $ decays of ¹³⁶Xe. The next-generation experiments aim for $ |m_{\beta\beta}| \approx 10\; {\rm {meV}} $, which has been set up as it is the lower boundary of $ |m_{\beta\beta}| $ in the IO case [4,5]. Therefore, an urgent question is whether it is realistic to reach $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ in the near future. The studies in Ref. [23] have demonstrated that if the JUNO-LS detector is upgraded with ¹³⁶Xe-loading in future, a sensitivity of $ |m_{\beta\beta}| \approx 5\; {\rm {meV}} $ is achievable when the most optimistic value of the NME is taken. In order to improve the sensitivity to $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $, one has to remarkably increase the target mass and reduce the backgrounds. Although the radioactive and cosmogenic backgrounds can be rejected by severe radiopurity control and perfect muon veto strategies, the irreducible backgrounds arising from recoiled electrons due to their elastic scattering with solar ⁸B neutrinos and the two-neutrino-emitting double-beta decays could be a serious problem. Nevertheless, inspired by the work in Ref. [23], it is promising to achieve $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ by developing the ¹³⁰Te-loaded LS to reach a sufficiently large target mass and advancing the powerful techniques for background reduction.

Strictly speaking, one should perform a statistical analysis of the experimental sensitivities to neutrino masses and the Majorana CP phases [43-46]. However, we have found that the final results are quite consistent with the simple analysis in the present paper. We believe that our analysis is very suggestive for setting up the future program for $ 0\nu\beta\beta $ decay experiments. If the sensitivity of $ |m_{\beta\beta}| \approx 1\; {\rm {meV}} $ is ultimately realized [47],^② the determination of absolute neutrino masses and the constraints on one of two Majorana CP phases are possible, which cannot be accessible in other types of feasible neutrino experiments. To achieve this goal, it is compulsory to make great efforts in increasing the target mass and reducing the background by two orders of magnitude compared to the present design of next-generation $ 0\nu\beta\beta $ decay experiments.