Constraints on light neutrino parameters derived from the study of neutrinoless double beta decay

The study of the neutrinoless double beta ($0 \beta\beta$) decay mode can provide us with important information on the neutrino properties, particularly on the electron neutrino absolute mass. In this work we revise the present constraints on the neutrino mass parameters derived from the $0 \beta\beta$ decay analysis of the experimentally interesting nuclei. We use the latest results for the phase space factors (PSFs) and nuclear matrix elements (NMEs), as well as for the experimental lifetimes limits. For the PSFs we use values computed with an improved method reported very recently. For the NMEs we use values chosen from literature on a case-by-case basis, taking advantage of the consensus reached by the community on several nuclear ingredients used in their calculation. Thus, we try to restrict the range of spread of the NME values calculated with different methods and, hence, to reduce the uncertainty in deriving limits for the Majorana neutrino mass parameter. Our results may be useful to have an up-date image on the present neutrino mass sensitivities associated with $0 \beta\beta$ measurements for different isotopes and to better estimate the range of values of the neutrino masses that can be explored in the future double beta decay (DBD) experiments.


Introduction
Neutrinoless double beta decay is a beyond Standard Model(BSM) process by which an even-even nucleus transforms into another even-even nucleus with the emission of two electrons/positrons but no antineutrinos/neutrinos in the final states. Its study is very attractive because it would clarify the question about the lepton number conservation, decide on the neutrinos character (are they Dirac or Majorana particles?) and give a hint on the scale of their absolute masses. Moreover, the study of the oνββ decay has a broader potential to search for other BSM phenomena. The reader can find up-todate information on these studies from several recent reviews [1]- [6], which also contain therein a comprehensive list of references in the domain.
The scale of the absolute mass of neutrinos is a key issue for understanding the neutrino properties. It cannot be derived from neutrino oscillation experiments which can only measure the square of the neutrino mass differences between different flavors [7]- [12]. Analysis of 0νββ decay and cosmological data are at present the most sensitive ways to investigate this issue.
The lifetimes of the 0νββ decay modes can be expressed, in a good approximation, as a product of a phase space factor (depending on the atomic charge and energy released in the decay, Q ββ ), a nuclear matrix element (related to the nuclear structure of the parent and daughter nuclei), and a lepton number violation (LNV) parameter (related to the BSM mechanism considered). Thus, to extract reliable limits for the LNV parameters we need accurate calculations of both PSFs and NMEs, as well as reliable measurements of the lifetimes.
The largest uncertainties in theoretical calculations for DBD are related to the NMEs values. That is why there is a continuous effort in literature to develop improved nuclear structure methods for their computation. At present, the NMEs are computed by several methods which differ conceptually, the most employed being the proton-neutron Quasi Random Phase Approximation (pnQRPA) [13]- [17], Interacting Shell Model(ISM) [18]- [21], Interacting Boson Approximation (IBA) [22]- [24], Projected Hartree Fock Bogoliubov (PHFB) [25] and Energy Density Functional (EDS) method [26]. There are still large differences between the NMEs values computed with different methods and by different groups, and these have been largely discussed in the literature (see for example [5]- [6]). On the other side, there is a consensus in the community on the way that several nuclear effects and nuclear parameters should be used in calculations. In this work, we take advantage of this consensus when we chose the NMEs values, trying to restrict their range of spread, and consequently, to reduce the uncertainty in deriving the neutrino Majorana mass parameters.
Unlike the NMEs, the PSFs have been calculated a long time ago [27]- [33] and were considered to be computed with enough precision. However, recently, they were recalculated within an improved approach, by using exact electron Dirac wave functions (w.f.) taking into account the finite nuclear size and electron screening effects [34]. The authors found differences between their results and those calculated previously with approximate electron w.f., especially for heavier nuclei. We have also independently recalculated the PSFs by developing new routines for computing the relativistic (Dirac) electron w.f. by taking into account the nuclear finite size and screening effects. In addition, we use a Coulomb potential derived from a realistic proton density distribution in the daughter nucleus [35]- [36]. In this work we use new PSFs values obtained by improving the numerical precision of our routines as compared with our previous works. The obtained values are very close to those reported in refs. [34]- [35].
Finally, for the lifetimes limits, we take the most recent results found in literature. The paper is organized as follows. In the next section we shortly remind the general formalism for the derivation of the neutrino mass parameters from 0νββ decay analysis, highlighting the nuclear ingredients involved in calculations. In Section 3 we discuss the way of choosing the NME values and report our results for the light neutrino Majorana mass parameters, while in Section 4 we formulate the conclusions of our work.

Formalism
We shortly present the general formalism for the derivation of neutrino mass parameters from 0νββ decay analysis. We start with the lifetime formula and then describe the main steps and ingredients used in the theoretical calculation of their components, i.e. PSFs and NMEs.
Assuming that the dominant mechanism of occurrence for the 0νββ decay mode is the exchange of Majorana left-handed light neutrinos between two nucleons from the parent nucleus, the lifetime reads: where G 0ν are the PSFs for this decay mode, depending on the energy decay Q ββ and nuclear charge Z, M 0ν are the corresponding NMEs, depending on the nuclear structure of the parent and daughter nuclei involved in the decay, m e is electron mass and m ν is the light neutrino Majorana mass parameter. This parameter can be expressed as a (coherent) linear combination of the light neutrino masses: where U ek are the elements of the first row in the PMNS (Pontecorvo-Maki-Nakagawa-Sakata) neutrino matrix and m k are the light neutrino masses [37]. From Eq. (1) the expression of m ν reads: For deriving m ν we need accurate calculations of both PSFs and NMEs, for each isotope for which there are experimental lifetime limits. The PSFs have been calculated a long time ago in some approximations [27]- [33] and were considered, until recently, to be computed with enough precision. However, they were recalculated recently in refs. [34]- [36] using more advanced approaches for the numerical evaluation of the Dirac wave functions with the inclusion of nuclear finite size and screening effects. In addition, in ref. [35] the usual Coulomb spherical potential was replaced by another one, derived from a more realistic proton density distribution in the daughter nucleus. These recent PSF calculations led to significant differences in the comparison to the older calculations, especially for the heavier isotopes, that should be taken into account for a precise derivation of the neutrino mass parameters.
The computation of the NMEs is a subject of debate in the literature for long time, because they bring the large uncertainties in the theoretical calculations for DBD. Different groups have developed several conceptually different nuclear structure methods [13]- [26], as we have mentioned in the previous Section. The expression of the NMEs can be written, in general, as a sum of three components: where M 0ν GT , M 0ν F and M 0ν T are the Gamow-Teller (GT ), Fermi(F ) and Tensor(T ) components, respectively. These are defined as follows: where O α mn are transition operators (α = GT, F, T ) and the summation is performed over all the nucleon states. An important part of the NME calculation is the computation of the reduced matrix elements of the two-body transition operators O α . Their calculation can be decomposed into products of reduced matrix elements within the spin and relative coordinate spaces. Their explicit expressions are [3], [21]: The most difficult is the computation of the radial part of the two-body transition operators, which contains the neutrino potentials. These potentials depend weakly on the intermediate states and are defined by integrals of momentum carried by the virtual neutrino exchanged between the two nucleons [16]: where R = r 0 A 1/3 fm (r 0 = 1.2f m), ω = q 2 + m 2 ν is the neutrino energy and j i (qr) is the spherical Bessel function (i = 0, 0 and 2 for GT, F, and T, respectively). Usually, in calculations one uses the closure approximation which consists of a replacement of the excitation energies of the states in the intermediate odd-odd nucleus contributing to the decay, by an average expression E . This approximation works good in the case of 0νββ decay modes and simplifies much the calculations. The expressions of h α (α = GT, F, T ) are: and where m π is the pion mass, m p is the proton mass and with (µ p − µ n ) = 4.71. The expressions (9)-(10) include important nuclear ingredients that should be taken into account for a precise computation of the NMEs, such as the higher order currents in the nuclear interaction (HOC) and finite nucleon size effect (FNS). Inclusion of HOC brings additional terms in the H GT component and leads to the appearance of the H T component in the expressions of the neutrino potentials. FNS effect is taken into account through G V and G A form factors: For the vector and axial coupling constants the majority of the calculations take either the quenched value, g V = 1 or the unquenched one, g A = 1.25, while the values of the vector and axial vectors form factors are Λ V = 850MeV and Λ A = 1086MeV [1], respectively. As one can see, when HOC and FNS corrections are included in the calculations, the dependence of NMEs expression on g A is not trivial and the NMEs values obtained with the quenched or the unquenched value of this parameter can not be obtained from by simply re-scaling.
To compute the radial matrix elements nl|H α |n ′ l ′ an important ingredient is the adequate inclusion of SRCs, induced by the nuclear interaction. The way of introducing the SRC effects has also been subject of debate ( [14], [16]). The SRC effects are included by correcting the single particle w. f. as follows: The correlation function f (r) can be parametrized in several ways. There are three parameterizations which are used, Miller-Spencer (MS), UCOM and CCM (with CD-Bonn and AV18 potentials). The Jastrow prescription [32] for the correlation function is: and includes all these parameterizations, depending of values of the a, b, c parameters. Including HOC and FNS effects, the radial matrix elements of the neutrino potentials become: where ν is the oscillator constant and V α (q) is an expression containing the q dependence of the neutrino potentials. From Eqs. (4)- (15) one can see that a set of approximations and parameters are involved in the NMEs expressions, as the HOC and FNS and SRC effects, and g A , r 0 , (Λ A , Λ B ), < E > parameters. Are there any recommendations on how should they be included in the calculations? At present there is a general consensus in the community in this respect, that will be discussed in the next section.

Numerical results and discussions
The neutrino mass parameters are derived from Eq. (3). To get < m ν > in the same units as m e we take the NMEs dimensionless and the PSFs (G 0ν ) in units of [yr] −1 .
The PSF values were obtained by recalculating them with our code, developed in [35], but with improved numerical precision. At this point we mention that the improved PSF values come, on the one hand, by the use of a Coulomb potential describing a more realistic proton charge density in the daughter nucleus instead of the (usual) constant charge density one, to solve the Dirac equations for obtaining the electron w.f. On the other hand, we got a better precision of our numerical routines that compute the PSFs by enhancing the number of the interpolation points on a case-to-case basis until the results become stationary. The obtained values are very close to both those reported previously in refs. [34] and [35]. This gives us confidence on their reliability. We mention that these PSFs values differ from older calculations as, for example, those reported in refs. [29], [31], [33] by up to 28%. Such differences are important for precise estimations and justify the re-actualization of the PSFs values in extracting Majorana neutrino mass parameters.
For the experimental lifetimes we took the most recent results reported in literature. Especially, we remark the newest results for 76 Ge from GERDA [38] and for 136 Xe from EXO [39]. The largest uncertainty in the derivation of m ν comes from the values of the NMEs. Fortunately, at present there is a general consensus in the community on the employment of the different nuclear effects (approximations) and parameters which appear in the NMEs expressions (Eqs. (4)-(15)) [40]. Thus, one can restrict the range of spread of the NMEs values for a particular nucleus, if one takes into account some recommendations resulting from the analysis of many NMEs calculations. For example, one recommends the inclusion in calculation of the HOC, FNS and SRC effects (although their effects can partially compensate each other [41]). For SRCs, softer parametrizations like UCOM [14] and CCM [42]- [43] are recommended, while the MS produces a too sever cut of the w.f. for very short inter-nucleon distances, which reflects into smaller NMEs values. Concerning the nuclear parameters, one recommends the use of an unquenched value for the g A axial vector constant, the values specified above for the vector and axial vector form factors (Λ V , Λ A ), and a value of r 0 = 1.2f m for the nuclear radius constant. The value for the average energy ( E ), used in the closure approximation, is a function of atomic mass A, but the results are less sensitive to changes within a few MeV. The use in different ways of these ingredients can result in significant differences between the NMEs values. Hence, a consensus is useful to approach the results obtained by different groups. Having agreement on these nuclear ingredients, the differences in the NMEs values should be searched in the features of the different nuclear structure methods. These methods use different ways of building the wave functions, different specific model spaces and type of nucleon-nucleon correlations and use some specific parameters [5], [20], [41]. Unfortunately, the uncertainties in the NMEs calculation associated with a particular nuclear structure method can not be easily fixed and they are still a subject of debate. As a general feature, ShM calculations underestimate the NMEs values (due to the limitations of the model spaces used), while the other methods overestimate them. There are, however, a few hints on how to understand/bring closer the NMEs results obtained with different methods. One idea would be to analyze the structure of the wave functions used in terms of the seniority scheme [40]. Another one, is to (re)calculate the NME values as to reproduce s.p. occupancies numbers measured recently for 76 Ge and 82 Se nuclei [44]. For example, when the QRPA calculations have modified with the s.p. energies that reproduce the experimental occupancies, the new QRPA NMEs values are much closer to the ShM ones.
In Table 1 we display the NMEs values obtained with different nuclear methods. For uniformity and in agreement with the consensus discussed above, we chose those results that were performed with the inclusion of HOC, FNS and SRC(UCOM and CD-Bonn) effects, and with unquenched g A = 1.25, as nuclear ingredients. We mention that the newest experimental determinations of this parameter report values even larger (1.269, 1.273) [45]. However, the differences between NMEs values obtained with g A = 1.25 -1.273 are not significant [41]. The NMEs values for 76 Ge and 82 Se written in parenthesis represent the adjusted NMEs values obtained with QRPA method by the Tuebingen and Jyvaskyla groups, when the s.p. energies were adjusted to the occupancy numbers reported in ref. [44]. One remarks, the QRPA calculations with s.p. occupancies in accordance with experiment, get significantly close to the ShM results, which is remarkable. In the future, one expects measurements of the occupancy numbers for other nuclei, as well. Also, it would be interesting if other methods, besides QRPA, would try to recalculate the NMEs by adjusting s.p. energies to experimental occupancy numbers.
We also make some remarks about the NMEs values on a case-by-case basis. For 48 Ca we appreciate that ShM calculations give more realistic results than the other methods. In support of this claim we mention, in the case of this isotope, ShM calculations are performed within the full pf shell and using good effective NN interactions, checked experimentally on other spectroscopic quantities [41], [46]- [47]. Also, we remark that NMEs values obtained with ShM for this isotope were used to correctly predict T 2ν 1/2 , before its experimental measurement [48]. For the isotopes with A =96 -130 there is a larger spread of the NMEs values calculated with different methods and, consequently, a larger uncertainty in predicting the m ν parameters. For 136 Xe there are new ShM large scale calculations with inclusion of a larger model space than the older calculations [49]. For this isotope the NMEs values are more grouped. Corroborated with a quite good experimental lifetime, from this isotope one gets the most stringent constraint for the m ν parameter.
In Table 2 we present our results for the Majorana neutrino mass parameters ( m ν ) together with the values of Q ββ , the PSFs (G 0ν ), NMEs (M 0ν ) and experimental lifetimes (T 0ν 1/2 ) for all the isotopes for which data exists. Making a sort of the NMEs values from literature according to the considerations presented, we reduce the interval of their spread to about a factor of 2, even less (with one exception). This results in reducing the uncertainty in deriving the constraints on the light neutrino Majorana mass parameters, while taking into account NME values obtained with all the main nuclear methods existent on the market. One observes that the stringent constraints are obtained from the 136 Xe isotope, followed by the 76 Ge one. This is due to both the experimental sensitivity of the experiments measuring these isotopes and to the reliability of the PSFs and NMEs theoretical calculations. The experiments measuring these isotopes are already exploring the quasi-degenerate scenarios for the neutrino mass hierarchy (which around 0.5 eV). With the ingredients presented in Table 2 (PSFs and NMEs) one can appreciate, as well, the sensitivities, translated into neutrino mass parameters, of the future generation of DBD experiments.

Conclusions
We report new values of light Majorana neutrino mass parameters from a 0νββ decay analysis extended to all the isotopes for which theoretical and experimental data exists. We used the most recent results for the experimental lifetimes, T 0ν 1/2 as well as for the theoretical quantities G 0ν and M 0ν . For the PSFs we use newly obtained values, recalculated with an approach described in ref. [35] but with improved numerical accuracy. We use exact electron w.f. obtained by solving a Dirac equation when finite nuclear size and screening effects are included and, in addition, a Coulomb potential derived from a realistic proton distribution in the daughter nucleus has been employed. For choosing the NMEs we take advantage on the general consensus in the community on several nuclear ingredients involved in the calculations (HOC, FHS and SRCs effects, values of several nuclear input parameters) and restrict the range of spread of the NMEs values, reported in the literature. This, in turn, reduces the uncertainty in deriving constraints on the light Majorana neutrino mass parameters, while taking into account NME values obtained with all the main nuclear methods. The stringent constraints are obtained from the 136 Xe and 76 Ge isotopes, due to both the experimental sensitivity and to the reliability of the PSFs and NMEs calculations. The experiments measuring these isotopes are already exploring the quasi-degenerate scenarios for the neutrino mass hierarchy which is around 0.5 eV. Our results may be useful for having an up-to-date image on the current neutrino mass sensitivities associated with 0νββ measurements for different isotopes and to better estimate the range of the neutrino masses that can be explored in the future DBD experiments.