The Majorana nature of neutrinos, as favored by many theoretical models, is a key for understanding of tiny neutrino masses observed in neutrino oscillation experiments. A golden process for answering this open question of particle physics is the neutrinoless double beta decay (0νββ-decay) [1–3],

in which an atomic nucleus with Z protons decays to another one with two more protons and the same mass number A, by emitting two electrons and nothing else. The observation of this process, which violates total lepton number conservation and is forbidden in the Standard Model, guaranties that neutrinos are Majorana particles, i.e., their own antiparticles [4].

The searches for the 0νββ-decay have not yielded any evidence for Majorana neutrinos yet. This could be because neutrinos are Dirac particles, i.e., not their own antiparticles. In this case we will never observe the decay. However, it is assumed that the reason for it is not sufficient sensitivity of previous and current 0νββ-decay experiments to the occurrence of this rare process.

Due to the evidence for neutrino oscillations and therefore for 3 neutrino mixing and masses the 0νββ-decay mechanism of primary interest is the exchange of 3 light Majorana neutrinos interacting through the left-handed V-A weak currents (m_ββ mechanism). In this case, the inverse 0νββ-decay half-life is given by Vergados et al. [1], DellOro et al. [2] and Vergados et al. [3]

where G₀₁, g_A and M_ν represent an exactly calculable phase space factor, the axial-vector coupling constant and the nuclear matrix element (whose calculation represents a severe challenge for nuclear theorists), respectively. m_e is the mass of an electron. The effective neutrino mass,

is a linear combination of the three neutrino masses m_i, weighted with the square of the elements U_ei of the first row of the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) neutrino mixing matrix. The measured value of m_ββ would be a source of important information about the neutrino mass spectrum (normal or inverted spectrum), absolute neutrino mass scale and the CP violation in the neutrino sector. However, that is not the only possibility.

There are several different theoretical frameworks that provide various 0νββ-decay mechanisms, which generate masses of light Majorana neutrinos and violate the total lepton number conservation. One of those theories is the left-right symmetric model (LRSM) [5–9], in which corresponding to the left-handed neutrino, there is a parity symmetric right-handed neutrino. The parity between left and right is restored at high energies and neutrinos acquire mass through the see-saw mechanism, what requires presence of additional heavy neutrinos. In general one cannot predict the scale where the left-right symmetry is realized, which might be as low as a few TeV—accessible at Large Hadron Collider, or as large as GUT scale of 10¹⁵ GeV.

The LRSM, one of the most elegant theories beyond the Standard Model, offers a number of new physics contributions to 0νββ-decay, either from right-handed neutrinos or Higgs triplets. The main question is whether these additional 0νββ-mechanisms can compete with the m_ββ mechanism and affect the 0νββ-decay rate significantly. This issue is a subject of intense theoretical investigation within the TeV-scale left-right symmetry theories [10–14]. In analysis of heavy neutrino mass mechanisms of the 0νββ-decay an important role plays a study of related lepton number and lepton flavor violation processes in experiments at Large Hadron Collider [2, 15–19].

The goal of this article is to discuss in details the W_L − W_R exchange mechanism of the 0νββ-decay mediated by light neutrinos (λ mechanism) and its coexistence with the standard m_ββ mechanism. For that purpose the corresponding nuclear matrix elements (NMEs) will be calculated within the quasiparticle random phase approximation with a partial restoration of the isospin symmetry [20] by taking the advantage of improved formalism for this mechanism of the 0νββ-decay of Štefánik [21]. A possibility to distinguish m_ββ and λ mechanisms in the case of observation of the 0νββ-decay on several isotopes will be analyzed. Further, the dominance of any of these two mechanisms in the 0νββ-decay rate will be studied within seesaw model with right-handed gauge boson at TeV scale. We note that a similar analysis was performed by exploiting a simplified 0νββ-decay rate formula and different viable particle physics scenarios in Tello et al. [10], Barry and Reodejohann [11], Bhupal Dev et al. [12], Deppisch et al. [13] and Borah et al. [14].

Recently, the 0νββ-decay with the inclusion of right-handed leptonic and hadronic currents has been revisited by considering exact Dirac wave function with finite nuclear size and electron screening of emitted electrons and the induced pseudoscalar term of hadron current, resulting in additional nuclear matrix elements [21]. In this section we present the main elements of the revisited formalism of the λ mechanism of the 0νββ-decay briefly. Unlike in Štefánik et al. [21] the effect the weak-magnetism term of the hadron current on leading NMEs is taken into account.

If the mixing between left and right vector bosons is neglected, for the effective weak interaction hamiltonian density generated within the LRSM we obtain

Here, G_β = G_F cos θ_C, where G_F and θ_C are Fermi constant and Cabbibo angle, respectively. The coupling constant λ is defined as

Here, M_WL and M_WR are masses of the Standard Model left-handed W_L and right-handed W_R gauge bosons, respectively. The left- and right-handed leptonic currents are given by

The weak eigenstate electron neutrinos ν_eL and ν_eR are superpositions of the light and heavy mass eigenstate Majorana neutrinos ν_j and N_j, respectively. We have

Here, U, S, T, and V are the 3 × 3 block matrices in flavor space, which constitute a generalization of the Pontecorvo-Maki-Nakagawa-Sakata matrix, namely the 6 × 6 unitary neutrino mixing matrix [22]

The nuclear currents are, in the non-relativistic approximation, [23]

Here, m_N is the nucleon mass. qV≡qV(q2), qA≡qA(q2), qM≡qM(q2) and qP≡qP(q2) are, respectively, the vector, axial-vector, weak-magnetism and induced pseudoscalar form-factors. The nucleon recoil terms are given by

where q_n = p_n − pn′ is the momentum transfer between the nucleons. The initial neutron (final proton) possesses energy En′ (E_n) and momentum pn′ (p_n). r_n, τn+ and σ→n, which act on the n-th nucleon, are the position operator, the isospin raising operator and the Pauli matrix, respectively.

By assuming standard approximations [21] for the 0νββ-decay half-life we get

The effective lepton number violating parameters η_ν (W_L − W_L exchange), η_λ (W_L − W_R exchange) and their relative phase Ψ are given by

The coefficients C_I (I = mm, mλ and λλ) are linear combinations of products of nuclear matrix elements and phase-space factors:

The explicit form and calculated values of phase-space factors G_0i (i = 1, 2, 3, 4, 10 and 11) of the 0νββ-decaying nuclei of experimental interest are given in Štefánik et al. [21]. The NMES, which constitute the coefficients C_I in Equation (13), are defined as follows:

The partial nuclear matrix elements M_I, where I = GT, F, T, ωF, ωGT, ωT, qF, qGT, and qT are given by

Here, O_F,GT,T are the Fermi, Gamow-Teller and tensor operators 1,σ→1·σ→2 and 3(σ→1·r^12)(σ→2·r^12). The two-nucleon exchange potentials h_I(r) with I = F, GT, T, ωF, ωGT, ωT, qF, qGT, and qT can be written as

where

and

Here, E_i, E_f and Ē_n are energies of the initial and final nucleus and averaged energy of intermediate nuclear states, respectively. r = (r_r − r_s), r_{r, s} is the coordinate of decaying nucleon and j_i(qr) (i = 1, 2, 3) denote the spherical Bessel functions. pr+pr′≃0, Er - Er′≃0 and pr - pr′≃q, where q is the momentum exchange. The form factors gV(q2), gA(q2), gM(q2) and gP(q2) are defined in Simkovic et al. [24]. We note that factor 4 in definition of the two-nucleon exchange potentials h_I(r) with I = ωF, ωGT, and ωT in Equation (48) of Štefánik et al. [21] needs to be replaced by factor 2.

The nuclear matrix elements are calculated in proton-neutron quasiparticle random phase approximation with partial restoration of the isospin symmetry for ⁴⁸Ca, ⁷⁶Ge, ⁸²Se, ⁹⁶Zr, ¹⁰⁰Mo, ¹¹⁰Pd, ¹¹⁶Cd, ¹²⁴Sn, ¹³⁰Te and ¹³⁶Xe, which are of experimental interest. In the calculation the same set of nuclear structure parameters is used as in Simkovic et al. [20]. The pairing and residual interactions as well as the two-nucleon short-range correlations derived from the realistic nucleon-nucleon Argonne V18 potential are considered [26]. The closure approximation for intermediate nuclear states is assumed with (Ē_n − (E_i + E_f)/2) = 8 MeV. The free nucleon value of axial-vector coupling constant (g_A = 1.25 − 1.27) is considered.

In Table 1 the calculated NMEs are presented. The values of M_F,GT,T and M_ν differ slightly (within 10%) with those given in Simkovic et al. [20], which were obtained without consideration of the closure approximation. By glancing Table 1 we see that M_{Fω, GTω, Tω} ≃ M_F,GT,T and M_νω ≃ M_ν as for the average neutrino momentum q = 100 MeV and used average energy of intermediate nuclear states we have q/(q + Ē_n − (E_i + E_f)/2) ≃ 1. The absolute value of M_{Fq, GTq, Tq} is smaller in comparison with M_F,GT,T by about 50% for Fermi NMEs and by about factor two in the case of Gamow-Teller and tensor NMEs. From Table 1 it follows that there is a significant difference between results of this work and the QRPA NMEs of Muto et al. [25], especially in the case of ¹⁰⁰Mo. This difference can be attributed to the progress achieved in the 0νββ-decay formalism due to inclusion of higher order terms of nucleon currents [21, 24], the way of adjusting the parameters of nuclear Hamiltonian [27], description of short-range correlations [26] and restoration of the isospin symmetry [20].

Nuclear matrix elements M2-, M1+ (λ mechanism) and M_ν (m_ββ mechanism) for 10 nuclei under consideration are given in Table 1 and displayed in Figure 1. We note a rather good agreement between M2- and M_ν for all calculated nuclear systems. It is because the contribution of M₁₊ to M₂₋ is suppressed by factor 9 and as a result M₂₋ is governed by the M_νω contribution (see Equation 14). Values of M1+ exhibit similar systematic behavior in respect to considered nuclei as values of M_ν and M2-, but they are suppressed by about factor 2–3 (with exception of ⁴⁸Ca).

The importance of the m_ββ and λ mechanisms depends, respectively, not only on values of η_ν and η_λ parameters, which are unknown, but also on values of coefficients C_I (I = mm, mλ, λλ), which are listed for all studied nuclei in Table 1. They have been obtained by using improved values of phase-space factors G_0k (k = 1, 2, 10 and 11) from Štefánik et al. [21]. We note that the squared value of M_GT and fourth power of axial-vector coupling constant g_A are included in the definition of coefficient C_I unlike in Štefánik et al. [21]. We see that C_λλ is always larger when compared with C_mm. The absolute value of C_mλ is significantly smaller than C_mm and C_λλ. This fact points out on less important contribution to the 0νββ-decay rate from the interference of m_ββ and λ mechanisms.

For 10 nuclei of experimental interest the decomposition of coefficient C_λλ (see Equation 11) on partial contributions CI0k associated with phase-space factors G_0k (k = 2, 10, and 11) is shown in Figure 2. By glancing the plotted ratio CI0k/CI we see that C_λλ is dominated by a single contribution associated with the phase-space factor G₀₂. From this and above analysis it follows that 0νββ-decay half-life to a good accuracy can be written as

with

For a given isotope the factor f_λm reflects relative sensitivity to the m_ββ and λ mechanisms and fλmG is its approximation, which does not depend on NMEs. The values f_λm and fλmG are tabulated in Table 1 and plotted as function of Q_ββ in Figure 3. We see that f_λm depends only weakly on involved nuclear matrix elements (apart for the case of ⁴⁸Ca) what follows from a comparison of f_λm with fλmG. The value of f_λm is mainly determined by the Q-value of double beta decay process. From 10 analyzed nuclei the largest value of f_mλ is found for ⁴⁸Ca and the smallest value for ⁷⁶Ge. A larger value of f_λm means increased sensitivity to m_ββ mechanism in comparison to λ mechanism and vice versa.

Upper bounds on the effective neutrino mass m_ββ and right-handed current coupling strength η_λ are deduced from experimental half-lives of the 0νββ-decay by using the coefficients C_mm, C_mλ and C_λλ of Table 1. The maximum and the value on axis (m_ββ = 0 or η_λ = 0) are listed in Table 2. The decays of ¹³⁶Xe and ⁷⁶Ge set the sharpest limit m_ββ ≤ 0.13 eV and 0.18 eV, and ηλ≤1.7 10-7 and 3.1 10⁻⁷, respectively. These are more stringent than those deduced from other experimental sources.

It is well known that by measuring different characteristics, namely energy and angular distributions of two emitted electrons, it is possible to identify which of m_ββ and λ mechanisms is responsible for 0νββ-decay [21, 23]. It might be achieved only by some of future 0νββ-decay experiments, e.g., the SuperNEMO [37] or NEXT [38]. A relevant question is whether the underlying m_ββ or λ mechanism can be revealed by observation of the 0νββ-decay in a series of different isotopes. In Figure 4 this issue is addressed by an illustrative case of observation of the 0νββ-decay of ¹³⁶Xe with half-life T1/20ν=6.86 1026 years, which can be associated with m_ββ = 50 meV or η_λ = 9.8 10⁻⁸. The 0νββ-decay half-life predictions associated with a dominance of m_ββ and λ mechanisms exhibit significant difference for some nuclear systems. We see that by observing, e.g., the 0νββ-decay of ¹⁰⁰Ge and ¹⁰⁰Mo with sufficient accuracy and having calculated relevant NMEs with uncertainty below 30%, it might be possible to conclude, whether the 0νββ-decay is due to m_ββ or λ mechanism.

Currently, the uncertainty in calculated 0νββ-decay NMEs can be estimated up to factor of 2 or 3 depending on the considered isotope as it follows from a comparison of results of different nuclear structure approaches [3]. The improvement of the calculation of double beta decay NMEs is a very important and challenging problem. There is a hope that due to a recent progress in nuclear structure theory (e.g., ab initio methods) and increasing computing power the calculation of the 0νββ-decay NMEs with uncertainty of about 30 % might be achieved in future.

The 6 × 6 unitary neutrino mixing matrix U (see Equation 8) can be parametrized with 15 rotational angles and 10 Dirac and 5 Majorana CP violating phases. For the purpose of study different LRSM contributions to the 0νββ-decay the mixing matrix U is usually decomposed as follows [22]

Here, 0 and 1 are the 3 × 3 zero and identity matrices, respectively. The parametrization of matrices A, B, R and S and corresponding orthogonality relations are given in Xing [22].

If A = 1, B = 1, R = 0 and S = 0, there would be a separate mixing of light and heavy neutrinos, which would participate only in left and right-handed currents, respectively. In this case we get η_λ = 0, i.e., the λ mechanism is forbidden.

If masses of heavy neutrinos are above the TeV scale, the mixing angles responsible for mixing of light and heavy neutrinos are small. By neglecting the mixing between different generations of light and heavy neutrinos, the unitary mixing matrix U takes the form

Here, m_D represents energy scale of charged leptons and m_LNV is the total lepton number violating scale, which corresponds to masses of heavy neutrinos. We see that U = U₀ can be identified to a good approximation with the PMNS matrix and V₀ is its analogue for heavy neutrino sector. Due to unitarity condition we find V0=U0†. Within this scenario of neutrino mixing the effective lepton number violating parameters η_ν (m_ββ mechanism) and η_λ (λ mechanism) are given by

with

The importance of m_ββ or λ-mechanism can be judged from the ratio

It is naturally to assume that ζ_m ≈ 1 and to consider the upper bound for the factor ζ_λ, i.e., there is no anomaly cancellation among terms, which constitute these factors. Within this approximation η_λ/η_ν does not depend on scale of the lepton number violation m_LNV and is plotted in Figure 5. The Dirac mass m_D is assumed to be within the range 1 MeV < m_D < 1 GeV. The flavor and CP-violating processes of kaons and B-mesons make it possible to deduce lower bound on the mass of the heavy vector boson M_W2>2.9 TeV [12]. From Figure 5 it follows that within accepted assumptions the λ mechanism is practically excluded as the dominant mechanism of the 0νββ-decay.

In this section the light-heavy neutrino mixing of the strength m_D/m_LNV is considered. However, we note that there are models with heavy neutrinos mixings where strength of the mixing decouples from neutrino masses [39–44]. This subject goes beyond the scope of this paper.

The left-right symmetric model of weak interaction is an attractive extension of the Standard Model, which may manifest itself in the TeV scale. In such case the Large Hadron Collider can determine the right-handed neutrino mixings and heavy neutrino masses of the seesaw model. The LRSM predicts new physics contributions to the 0νββ half-life due to exchange of light and heavy neutrinos, which can be sizable.

In this work the attention was paid to the λ mechanism of the 0νββ-decay, which involves left-right neutrino mixing through mediation of light neutrinos. The recently improved formalism of the 0νββ-decay concerning this mechanism was considered. For 10 nuclei of experimental interest NMEs were calculated within the QRPA with a partial restoration of the isospin symmetry. It was found that matrix elements governing the conventional m_ββ and λ mechanisms are comparable and that the λ contribution to the decay rate can be associated with a single phase-space factor. A simplified formula for the 0νββ-decay half-life is presented (see Equation 19), which neglects the suppressed contribution from the interference of both mechanisms. In this expression the λ contribution to decay rate is weighted by the factor f_λm, which reflects relative sensitivity to the m_ββ and λ mechanisms for a given isotope and depends only weakly on nuclear physics input. It is manifested that measurements of 0νββ-decay half-life on multiple isotopes with largest deviation in the factor f_λm might allow to distinguish both considered mechanisms, if involved NMEs are known with sufficient accuracy.

Further, upper bounds on effective lepton number violating parameters m_ββ (η_ν) and η_λ were deduced from current lower limits on experimental half-lives of the 0νββ-decay. The ratio η_λ/η_ν was studied as function of the mass of heavy vector boson M_WR assuming that there is no mixing among different generations of light and heavy neutrinos. It was found that if the value of Dirac mass m_D is within the range 1 MeV < m_D < 1 GeV, the current constraint on M_WR excludes the dominance of the λ mechanism in the 0νββ-decay rate for the assumed neutrino mixing scenario.

FŠ: calculation of nuclear matrix elements and preparation of manuscript; RD: analysis of lepton number violating parameters and preparation of manuscript; DŠ: derivation of the formalism of neutrinoless double beta decay, and analysis of obtained results.