# Two-neutrino double-beta decay matrix elements based on relativistic nuclear
energy density functional

N. Popara<EMAIL_ADDRESS>Faculty of Veterinary Medicine, University
of Zagreb, Heinzelova c. 102, 10000 Zagreb, Croatia  Department of Physics,
Faculty of Science, University of Zagreb, Bijenička c. 32, 10000 Zagreb,
Croatia  A. Ravlić<EMAIL_ADDRESS>Department of Physics, Faculty of Science,
University of Zagreb, Bijenička c. 32, 10000 Zagreb, Croatia  N. Paar
<EMAIL_ADDRESS>Department of Physics, Faculty of Science, University of Zagreb,
Bijenička c. 32, 10000 Zagreb, Croatia

###### Abstract

Nuclear matrix elements (NMEs) for two-neutrino double-beta decay
($2\nu\beta\beta$) are studied in the framework of relativistic nuclear energy
density functional (REDF). The properties of nuclei involved in the decay are
obtained using the relativistic Hartree-Bardeen-Cooper-Schrieffer theory and
relevant nuclear transitions are described using the relativistic proton-
neutron quasiparticle random phase approximation based on relativistic energy
density functional (REDF-QRPA). Three effective interactions have been
employed, including density-dependent meson-exchange (DD-ME2) and point
coupling interactions (DD-PC1 and DD-PCX), and pairing correlations are
described consistently both in $T=1$ and $T=0$ channels using a separable
pairing interaction. The optimal values of $T=0$ pairing strength parameter
$V_{0pp}$ are constrained by the experimental data on $\beta$-decay half
lives. The $2\nu\beta\beta$ matrix elements and half-lives are calculated for
several nuclides experimentally known to undergo this kind of decay: 48Ca,
76Ge, 82Se, 96Zr, 100Mo, 116Cd, 124Xe, 128Te, 130Te, 136Xe and 150Nd. The
model dependence of the NMEs and their sensitivity on $V_{0pp}$ is
investigated, and the NMEs obtained using optimal values of $V_{0pp}$ are
discussed in comparison to previous studies. The results of the present work
represent an important benchmark for the future applications of the
relativistic framework in studies of neutrinoless double-beta decay.

##  I Introduction

The study of double-beta decays has attracted considerable interest over the
past years. [1, 2]. Most of this interest has been focused on the neutrinoless
mode ($0\nu\beta\beta$), because it offers the possibility of distinguishing
between Dirac and Majorana nature of neutrinos [3]. The two-neutrino decay
mode ($2\nu\beta\beta$), nonetheless, remains of interest, as an important
benchmark of theoretical models for further studies of ($0\nu\beta\beta$)
decays. The $2\nu\beta\beta$ decay is a second-order weak interaction process,
and as such it is allowed by the Standard Model [4], unlike the neutrinoless
mode, which violates lepton conservation and consequently requires physics
beyond the standard model [5]. Furthermore, it is possible to experimentally
observe this mode, as pairing in even-even nuclei makes them more stable than
adjacent odd-odd nuclei, and the transition from an even-even mother to an
odd-odd daughter is forbidden energetically, leaving double-beta decay as the
only allowed decay channel [6]. This allows us to compare the results of our
calculations to experiment directly, and, since the calculations for the
neutrinoless mode use the same ingredients as those for the two-neutrino mode,
their results can be used to constrain some of the model parameters and
benchmark the model for studies of $0\nu\beta\beta$ decay [7].

The contribution of nuclear physics to the calculation of decay rates for both
modes of double-beta decay is contained chiefly in the nuclear matrix element
(NME) [8]. A variety of theoretical approaches has been developed in the past
years for the description of double-beta decay NMEs. More detailed information
is available in the following extensive reviews [9, 10, 6, 11, 12, 13, 14, 15,
16, 2, 17, 18]. The calculations of NMEs have been carried out using various
approximations: the quasiparticle random phase approximation (QRPA) [19, 20,
21], its extension to the renormalized QRPA (RQRPA) [22, 23, 24, 25], the
second QRPA (SQRPA) [26, 27], various energy density functional (EDF)
approaches [28, 29], the Quasiparticle Tamm-Dancoff Approximation [30], the
interacting shell model (ISM) and similar approaches [31, 32, 33], the
interacting boson model (IBM) [34, 35, 36], and others [37, 38, 39]. Among
these approximations, the proton-neutron quasiparticle random-phase
approximation (pn-QRPA) has emerged as one of the approaches that have been
successfully employed in various studies of double-beta decays [40, 1, 41,
42]. Within this approach, the usual isovector pairing in the ground state of
open-shell nuclei is supplemented with an isoscalar proton-neutron pairing in
the residual interaction. This has been shown to be important in nuclei
relevant for double-beta decays (as first established by Ref. [43] and
Ref.[44] and confirmed by subsequent studies [45, 46], including work beyond
the QRPA approximation [47, 48]). In the relativistic framework, various
methods have been employed in studies of $0\nu\beta\beta$ decays. In Ref. [49]
a multireference covariant EDF has been used to determine the wave functions
of the initial and final nuclei, and correlations beyond the mean field have
been described by configuration mixing of both angular momentum and particle
number projected quadrupole deformed mean-field wave functions. Systematic
study of $0\nu\beta\beta$ matrix elements in the relativistic EDF approach has
been reported in Ref. [50], and effects of the the relativity and short-range
correlations have recently been analyzed [51].

The aim of this work is to establish a theory framework for the description of
$2\nu\beta\beta$ decay NMEs based on a modern relativistic nuclear energy
density functional with density dependent point-coupling and meson-exchange
interactions that include density dependence explicitly in the interaction
vertex functions. These interactions have recently been introduced and
successfully implemented in the description of nuclear excitation properties
[52, 53, 54, 55, 56, 57], astrophysically relevant weak interaction processes
[58, 59, 60, 61], and nuclear equation of state [62, 63, 64]. In previous
relativistic EDF studies $2\nu\beta\beta$ decays have not been addressed [49,
50, 51]. We note that the present theory framework includes the density
dependence explicitly in the vertex functions of the meson-exchange and point
coupling interactions.

When modeling various quantities in nuclear physics, it is important to assess
systematic errors. As pointed out in Ref. [65] the systematic error of a
theoretical model may be a consequence of missing physics and/or poor
modeling. Considering that in most cases the perfect model is not available,
systematic errors are rather difficult to estimate. However, some insight
about systematic uncertainties can be obtained from a comparative study of
different theory frameworks and effective nuclear interactions. Therefore, it
is important to address the problem of double beta decays both from non-
relativistic and relativistic frameworks, using various formulations of energy
density functionals and different parameterizations.

Since for the $2\nu\beta\beta$ decay mode experimental data exist for a set of
nuclei, the present study also allows us to benchmark the relativistic model
for the future studies of $0\nu\beta\beta$ decay. The properties of nuclei
involved are described using the relativistic Hartree-Bardeen-Cooper-
Schrieffer (RH-BCS) model, while relevant transitions are obtained using the
proton-neutron relativistic quasiparticle random phase approximation [66],
that has recently been extended by including relativistic density-dependent
point coupling interactions [67, 68]. In the following we denote this method
as REDF-QRPA. Model calculations of the NMEs include various double-beta
emitters: 48Ca, 76Ge, 82Se, 96Zr, 100Mo, 116Cd, 124Xe, 128Te, 130Te, 136Xe,
and 150Nd. In contrast to recent studies based on the pn-QRPA, which tend to
start from a ”realistic” nucleon-nucleon interaction [69, 70], in the present
work in the particle-hole channel we use an interaction derived from a
relativistic nuclear energy density functional. An important aspect of this
study is that we are able to explore the model dependence of the calculated
NMEs, by implementing two different types of relativistic density-dependent
interactions, (i) finite range meson-exchange and (ii) point coupling
interactions. Our model calculations also involve a treatment of the pairing
correlations in open-shell nuclei, both at the level of the nuclear ground
state and in the residual REDF-QRPA interaction. In particular, strength of
the proton-neutron pairing, that is isoscalar $T=0$ pairing in the residual
REDF-QRPA interaction has to be constrained by using experimental data, in a
similar way as already discussed in previous studies in non-relativistic
frameworks [7]. In the present study the experimental data on single $\beta$\-
decays will be employed to constrain $T=0$ pairing strength parameter for
applications in $2\nu\beta\beta$ and in intended future $0\nu\beta\beta$
matrix element calculations.

We note that in this work, that represents the first implementation of the
relativistic density dependent interactions in $2\nu\beta\beta$ decay study,
some effects have not yet been taken into account, from deformation [71, 72]
to renormalization and gauge symmetry [73] and isospin restoration [74, 75].
In forthcoming studies more advanced effects will be taken into account.
Nonetheless we expect that the present work will describe well at least some
subset of the double-beta emitters, and will provide a useful guidance for the
future studies.

This article is organized as follows: in Section II we provide a theoretical
overview concerning the matrix elements involved in double-beta decay and in
Section III we outline the relativistic theory framework. Results are provided
and discussed in Section IV, and the conclusion follows in Section V.

##  II Two-neutrino double-beta decay

Two-neutrino double-beta decay ($2\nu\beta\beta$) is the process whereby two
neutrons in the mother nucleus are converted into protons, accompanied by the
emission of two electrons and two antineutrinos,

$(A,Z)\rightarrow(A,Z+2)+2e^{-}+2\bar{\nu}.$ (1)

The half-life of two-neutrino double-beta decay formally depends on two matrix
elements [76],

$\frac{1}{T_{1/2}^{2\nu}}=G^{2\nu}(Q,Z)g^{4}_{A}\left[\mathcal{M}^{2\nu}_{GT}+\mathcal{M}^{2\nu}_{F}\right]^{2},$
(2)

where $\mathcal{M}^{2\nu}_{GT}$ and $\mathcal{M}^{2\nu}_{F}$ are,
respectively, the Gamow-Teller (GT) and Fermi matrix element, and
$G^{2\nu}(Q,Z)$ is a phase space factor that can be found tabulated, e.g., in
Ref. [77]. The factor $g_{A}$ is the axial-vector coupling constant.

We restrict our consideration to decays from a 0+ state to a 0+ state. Decays
to the states of higher angular momenta, e.g., 2+, are possible, albeit highly
suppressed. Their treatment would require a computational apparatus more
involved than the REDF-QRPA used in this work [78]. The Gamow-Teller matrix
element in this case can be written as [18]:

$\mathcal{M}^{2\nu}_{GT}=m_{e}\sum_{m}\frac{\bra{f}\sum_{a}\sigma_{a}\tau^{-}_{a}\ket{m}\bra{m}\sum_{a}\sigma_{a}\tau^{-}_{a}\ket{i}}{E_{m}+\frac{Q}{2}+m_{e}},$
(3)

where the sum goes over 1+ states in the intermediate nucleus, labelled with
$m$ and having the energy Em as measured from the ground state of the initial
nucleus. The $Q$-value of the double-beta decay reads

$Q=\left(m_{i}-m_{f}\right)c^{2}-2m_{e}.$ (4)

Note that we could have also chosen to group the factor $g_{A}^{2}$ into the
definition of the matrix element in Eq. (3), in order to make the discussion
about the effects of $g_{A}$ and its renormalisation easier. If the two are
grouped together, one obtains an effective NME ${M}^{2\nu}_{GT,eff}$ that
depends on the choice of the axial-vector coupling as well [79]. In this work,
instead of the free nucleon value $g_{A}=1.26$, we employ a quenched value
$g_{A}=1.0$ that is consistent with the previous REDF-QRPA studies of
$\beta$-decay half-lives as well as electron capture rates [80, 81, 82]. Thus
our results for the matrix elements in Eq. (3) can be directly compared with
those of the effective NMEs from Ref. [79]. Likewise, we introduce the factor
me (in MeV) in the definition of the nuclear matrix element to obtain
dimensionless NMEs that can be compared easily to results of recent
theoretical and experimental works. Otherwise the matrix elements would be in
units of MeV-1. In all figures and tables in this work, the REDF-QRPA results
for ${M}^{2\nu}_{GT}$ and ${M}^{2\nu}_{GT,eff}$ are given as dimensionless.

The closure approximation is widely used in calculating $2\nu\beta\beta$
decay. It entails the replacement of the sum over different states in the
intermediary nucleus, as in Eq. (3), with a form that consists of one
denominator containing a suitably-chosen ”average” energy $\langle$E$\rangle$:

$\mathcal{M}^{2\nu}_{GT,closure}=m_{e}\frac{\bra{f}\sum_{a}\sigma_{a}\tau^{-}_{a}\sum_{a}\sigma_{a}\tau^{-}_{a}\ket{i}}{\langle
E\rangle+\frac{Q}{2}+m_{e}}.$ (5)

In this work we calculate the NMEs using the more complete approach given in
Eq. (3), but where relevant for the comparison with other studies we include
the closure approximation as well. It can be shown that the Fermi matrix
element vanishes if the same pairing interaction is consistently used at both
the ground state and excitation levels in the isovector channel [83]. This is
the case for our calculation and as such we will not show the Fermi matrix
element. The REDF-QRPA is not applicable to transitions from the 1+
intermediate nucleus. Therefore, although we conceive of double-beta decay as
a sequence of two $\beta^{+}$, $\beta^{-}$ or electron capture decays [84], in
the present study we calculate $\beta^{-}$ decay from the initial 0+ nucleus
and a $\beta^{+}$ decay from the final 0+ nucleus, resulting in two sets of
intermediate 1+ states. Explicitly, for the GT transitions the matrix element
is given by:

$\mathcal{M}{{}^{2\nu}_{GT}}=m_{e}\sum_{m^{\prime}m}\frac{\bra{m}\sum_{a}\sigma_{a}\tau^{-}_{a}\ket{i}\braket{m}{m^{\prime}}\bra{f}\sum_{a}\sigma_{a}\tau^{-}_{a}\ket{m^{\prime}}}{E_{m}+\frac{Q}{2}+m_{e}}.$
(6)

For the overlap between states belonging to different sets we take the usual
prescription [40]:

$\braket{m}{m^{\prime}}\approx\sum_{pn}\left[X^{m}_{pn}X^{m^{\prime}}_{pn}-Y^{m}_{pn}Y^{m^{\prime}}_{pn}\right],$
(7)

where the quantities X and Y are the REDF-QRPA amplitudes that will be defined
in the next section. This is an approximation in cases where the initial and
final states are not identical, but a reasonable one. For 48Ca, for example,
we calculate the overlap between the lowest states belonging to each of the
different sets of states in the intermediate nucleus to be 0.9855. Higher
states in one set are more likely to have significant overlap with several
states in the other, but the values of the overlaps still go from around 0.14
toward higher ones around 0.89. An alternative prescription from Šimkovic
[40], involves additional factors proportional to occupation numbers:

$\braket{m}{m^{\prime}}\approx\sum_{pn}\left[X^{m}_{pn}X^{m^{\prime}}_{pn}-Y^{m}_{pn}Y^{m^{\prime}}_{pn}\right]\tilde{u}_{p}\tilde{u}_{n},$
(8)

where the quantities $\tilde{u}_{p/n}$ are defined as follows:

$\tilde{u}_{p/n}=u_{p/n}^{m}u_{p/n}^{m^{\prime}}+v_{p/n}^{m}v_{p/n}^{m^{\prime}},$
(9)

and $u_{p/n}$, $v_{p/n}$ are occupation numbers derived, in the present work,
from the relativistic Hartree-BCS model (more details are given in Sec. III).
Several effects are usually neglected in the study of 2$\nu\beta\beta$ decays,
such as higher-order currents and realistic short-range correlations [85],
that are important in the neutrinoless case. This is due to the insensitivity
of the two-neutrino NMEs to the details of the nucleon wavefunctions for low
nucleon separation r12 [85].

##  III Relativistic framework for $2\nu\beta\beta$ decay matrix elements

###  III.1 Relativistic Hartree-BCS for the ground state description

Calculations of nuclear matrix elements for double-beta decay often proceed
starting from a realistic nucleon-nucleon potential [86, 16]. As mentioned in
Sec. I, other theory frameworks based on phenomenological effective
interactions have also been employed, in particular, shell model and nuclear
energy density functionals. In this work we introduce a framework for
2$\nu\beta\beta$ decays based on a relativistic nuclear energy density
functional [87]. In the relativistic nuclear energy density functional (REDF)
framework, the nuclear ground-state density and energy are determined by the
self-consistent solution of relativistic single-nucleon Kohn-Sham equations
[88, 89]. In the present study these equations are implemented through an
interaction Lagrangian density formulated in terms of the relevant degrees of
freedom. Since the REDF has already been extensively used in many previous
studies, here we give only a brief overview of the relevant foundations of
this framework. Two different families of relativistic density-dependent
interactions are used in this work, (i) finite range meson-exchange and (ii)
point coupling interactions. In the former case, pointlike nucleons interact
through the exchange of light mesons, namely $\omega$, $\rho$ and $\sigma$
mesons, in addition to an electromagnetic interaction mediated by photons. The
model is explained in detail in [90, 91, 92, 93, 87]. Meson-nucleon couplings
are established as functions of the vector density, motivated by the
relativistic Brueckner-Hartree-Fock calculations, but introduced through a
phenomenological ansatz with parameters adjusted to the experimental data in
finite nuclei [87]. In this work, the DD-ME2 parameterization of the density-
dependent meson-exchange interaction is used [90], being one of the most
successful parameterizations currently used in the description of a variety of
nuclear properties and astrophysically relevant processes [94, 95, 62].

In the case of point coupling interactions, the effective Lagrangian contains
four fermion contact interaction terms including the isoscalar-scalar,
isoscalar-vector, isovector-vector and isospace-space channels, coupling of
protons to the electromagnetic field, and the derivative term accounting for
the leading effects of finite-range interactions necessary for a quantitative
description of nuclear density distribution and radii (for more details see
Refs. [96, 87]). In our study of 2$\nu\beta\beta$ decays, two
parameterizations of the density dependent point coupling interactions are
used, DD-PC1 [96], and the more recently established DD-PCX [97]. While the
DD-PC1 interaction is adjusted to nuclear binding energies, the DD-PCX
interaction is specifically adjusted both to the nuclear ground state and
excitation properties, in order to constrain not only the nuclear properties
but also the symmetry energy close to the saturation density, and the
incompressibility of nuclear matter by using genuine observables on finite
nuclei in the $\chi^{2}$ minimization [97].

For the description of ground state properties of open-shell nuclei a unified
and self-consistent treatment of mean-field and pairing correlations is
needed. In this work we employ the relativistic Hartree-BCS (RH-BCS) theory
[81, 57] which represents a relativistic extension of the non-relativistic
Hartree-Fock-BCS framework [98, 99]. Spherical symmetry is assumed. The
pairing correlations in the ground state are described using a separable
pairing force, which also includes two parameters for the pairing strength
($G_{p}$ and $G_{n}$) [100]. While in the RH-BCS model with DD-ME2 and DD-PC1
interactions the pairing parameterization from Ref. [100] is used, the DD-PCX
interaction is supplemented with its own parameterization for the separable
pairing force as given in Ref. [97]. Within the RH-BCS model, the isoscalar
pairing is not included, i.e., no proton-neutron mixing is considered.

###  III.2 Proton-neutron relativistic QRPA

Charge-exchange transitions between the states in nuclides involved in
$2\nu\beta\beta$ decay are described in the framework of the proton-neutron
relativistic quasiparticle random phase approximation (REDF-QRPA). An
introduction to the charge-exchange QRPA (or proton-neutron QRPA, pn-QRPA)
calculations can be found in Refs. [101, 102], while in the relativistic
framework it has been introduced in Refs. [103, 66]. The residual REDF-QRPA
interaction is derived from the relativistic formulation of the effective
Lagrangian density, and throughout the calculation of the QRPA matrix
elements, the Dirac wave functions from the RH-BCS model, including both large
and small components are systematically included. In addition, the transition
operators are also extended for their implementation in the relativistic
framework. More details on the relativistic QRPA are given in Refs. [103, 66,
94, 104, 105]. Here we give only a brief overview of the REDF-QRPA adopted for
the purpose of the study of $2\nu\beta\beta$ decay. The states in the
intermediate $(Z+1,N-1)$ nucleus are REDF-QRPA phonons:

$\ket{m}=\sum_{pn}\left[X_{pn}a^{\dagger}_{p}a^{\dagger}_{n}-Y_{pn}a_{n}a_{p}\right]\ket{\text{QRPA}},$
(10)

where the creation operators $a^{\dagger}_{p(n)}$ create a proton(neutron)
state in an orbital labeled $p(n)$, and quantities $X$ and $Y$ are, as noted
earlier, the REDF-QRPA amplitudes. $\ket{\text{QRPA}}$ denotes the QRPA
vacuum, which we take to be the ground state (see Sec.III.1). By linearizing
the time-dependent RH-BCS equations in charge-exchange external field, the
charge-exhange QRPA equations are obtained [78]:

$\begin{pmatrix}A&B\\\ -B^{*}&-A^{*}\end{pmatrix}\begin{pmatrix}X(J)\\\
Y(J)\end{pmatrix}=\omega_{k}\begin{pmatrix}X(J)\\\ Y(J)\end{pmatrix},$ (11)

where the REDF-QRPA matrices $A$ and $B$ are defined as

$\displaystyle\begin{split}&A_{pnp^{\prime}n^{\prime}}(J)=(E_{p}+E_{n})\delta_{pp^{\prime}}\delta_{nn^{\prime}}\\\
&+(u_{p}u_{p^{\prime}}u_{n}u_{n^{\prime}}+v_{p}v_{p^{\prime}}v_{n}v_{n^{\prime}})\bra{pnJ}V\ket{p^{\prime}n^{\prime}J}\\\
&+(u_{p}v_{p^{\prime}}u_{n}v_{n^{\prime}}+v_{p}u_{p^{\prime}}v_{n}u_{n^{\prime}})\bra{pn^{-1}J}V_{res}\ket{p^{\prime}n^{\prime-1}J},\end{split}$
(12)

and

$\displaystyle\begin{split}&B_{pnp^{\prime}n^{\prime}}(J)=(u_{p}u_{p^{\prime}}v_{n}v_{n^{\prime}}+v_{p}v_{p^{\prime}}u_{n}u_{n^{\prime}})\bra{pnJ}V\ket{p^{\prime}n^{\prime}J}\\\
&+(u_{p}v_{p^{\prime}}v_{n}u_{n^{\prime}}+v_{p}u_{p^{\prime}}u_{n}v_{n^{\prime}})\bra{pn^{-1}J}V_{res}\ket{p^{\prime}n^{\prime-1}J},\end{split}$
(13)

where $V_{res}$ is the residual interaction derived from the relativistic
nuclear energy density functional, while $V$ includes the pairing interaction
in the QRPA [66, 67], $u_{p(n)},v_{p(n)}$ being proton(neutron) RH-BCS
amplitudes. The residual interaction $V_{res}$ is derived from the same
effective meson-exchange or point coupling interaction as used in the ground
state calculations (see Sec. III.1 ). In addition, the pseudovector
interaction channel is also included, and its strength parameter has been
previously constrained to the experimental data on Gamow-Teller resonance in
208Pb [66, 94, 67].

The REDF-QRPA equations include both the isovector $(T=1)$ and isoscalar
$(T=0)$ pairing channels, described by the separable interaction [100, 67,
68],

$V^{pn}(r_{1},r_{2})=-fG\delta(\boldsymbol{R}-\boldsymbol{R}^{\prime})P(r)P(r^{\prime})(1-P^{r}P^{\sigma}P^{\tau}),$
(14)

where the projectors $P^{r,\sigma,\tau}$ are defined as usual, $G$ is the
overall interaction strength, $R$ and $r$ are center of mass and relative
coordinates respectively, and $P(r)$ is defined as:

$P(r)=\frac{1}{(4\pi a^{2})^{3/2}e^{-\frac{r^{2}}{2a^{2}}}}.$ (15)

The overall strength in the pairing channel is multiplied by a dimensionless
factor $f$ defined as

$f=\left\\{\begin{array}[]{c}1,\quad T=1,S=0,\\\ V_{0pp},\quad T=0,S=1,\\\
0,\quad\text{ otherwise},\end{array}\right.$ (16)

where $V_{0pp}$ represents the isoscalar pairing strength. For $T=1$ channel,
the same separable pairing interaction is used as in the RH-BCS model (see
Sec. III.1). The strength parameters $G_{p(n)}$ for protons(neutrons) and the
width $a$ are defined in Ref. [106], but other parameterizations also exist,
e.g., the one constrained with the DD-PCX interaction [97]. For $T=0$ channel,
the pairing strength parameter $V_{0pp}$ has to be constrained at the level of
REDF-QRPA, e.g., by using experimental data on Gamow-Teller transitions or
$\beta$-decays. In our treatment of 2$\nu\beta\beta$ decay, this allows us to
explore the dependence of the NME on the pairing strength parameter $V_{0pp}$,
showing the variation in the possible values for the matrix elements. However,
to provide predictions on 2$\nu\beta\beta$ decay matrix elements, in this work
$V_{0pp}$ is also constrained by the experimental data on single $\beta$-decay
half-lives (Sec. III.3).

The transition matrix elements necessary for 2$\nu\beta\beta$ NME in Eq. (3)
are given by following expressions [43],

$\displaystyle\begin{split}&\bra{i}\sum_{a}\sigma_{a}\tau^{-}_{a}\ket{m}\\\
&=\sum_{pn}\langle
p||\boldsymbol{\sigma}||n\rangle\left[u_{p}v_{n}X^{m}_{pn}+v_{p}u_{n}Y^{m}_{pn}\right]\\\
\end{split}$ (17)
$\displaystyle\begin{split}&\bra{m}\sum_{a}\sigma_{a}\tau^{+}_{a}\ket{f}=\\\
&\sum_{pn}\langle
p||\boldsymbol{\sigma}||n\rangle\left[v_{p}u_{n}X^{m}_{pn}+u_{p}v_{n}Y^{m}_{pn}\right].\end{split}$
(18)

Here, the terms next to the $X$ and $Y$ amplitudes represent particle-type and
hole-type one-quasi-particle transitions, respectively [78].

###  III.3 Determination of the isoscalar pairing strength

In order to constrain the isoscalar ($T=0$) pairing strength $V_{0pp}$ we use
similar approach as suggested in Refs. [81, 80, 107], namely $V_{0pp}$ is
determined from the global fit of $\beta$-decay half-lives to experimental
data. We have used even-even nuclei in the range $8\leq Z\leq 82$ for which
experimental data on $\beta$-decay half-lives is available, where $V_{0pp}$
has the effect in the relative change of half-lives $T_{1/2}$ by more than
20$\%$, and whose half-lives are $<\text{10}^{3}$ s. In this way we have
additionally optimized the fitting procedure, because too long half-lives may
be challenging for quantitative description within the QRPA method [80]. Since
the isoscalar pairing strength $V_{0pp}$ shows isotopic dependence [80], in
the fitting procedure we use a functional form,

$V_{0pp}=V_{1}+V_{2}\left(\frac{N-Z}{A}\right)$ (19)

that appears to provide a comparable quality of the fit to $\beta$-decay half-
lives to other previously used functional forms [80]. The $\beta$-decay half-
lives are calculated using the REDF-QRPA as described in Ref. [81], including
both allowed and first-forbidden transitions. We have used the same quenched
value of axial-coupling $g_{A}=1.0$ as in this work. The fitting procedure
determines $V_{0pp}$ which reproduces the experimental $\beta$-decay half-
lives. Here we note that only those $V_{0pp}$ which yield real solutions of
the QRPA equation should be included in the fit to the proposed ansatz. Using
models with DD-ME2, DD-PC1 and DD-PCX interactions, the obtained average
values of the fitted parameters $V_{1}$ and $V_{2}$ as well as
$1\sigma-$uncertainties are given in Tab. 1. This results in rather narrow
range of values for the isoscalar pairing strength $V_{0pp}$ for the set of
nuclei considered in 2$\nu\beta\beta$ decay study in this work. The respective
values for $V_{0pp}$ with the uncertainties are given in Tab. 2. In the
following Section we refer to _optimal_ $V_{0pp}$ when using the values from
Tab. 2. We note that most of the $V_{0pp}$ values are close to 1.0, meaning
the strengths of the isoscalar and isovector pairing are very similar, which
in turns points to a ”soft” breaking of spin-isospin SU(4) symmetry in present
calculations.

Table 1: The fit parameters for the $T=0$ pairing strength functional form in Eq. (19) obtained from the optimization on $\beta$-decay half-lives for the models with DD-PCX, DD-PC1 and DD-ME2 interactions | DD-ME2 | DD-PC1 | DD-PCX
---|---|---|---
$V_{1}$ | 0.574 $\pm$ 0.338 | 0.522 $\pm$ 0.407 | 0.592 $\pm$ 0.366
$V_{2}$ | 2.301 $\pm$ 1.422 | 3.092 $\pm$ 1.768 | 2.321 $\pm$ 1.559
Table 2: Optimal values of the $T=0$ pairing strength parameter $V_{0pp}$ constrained from the $\beta$-decay half lives. |  | $V_{0pp}$ |
---|---|---|---
| DD-ME2 | DD-PC1 | DD-PCX
48Ca | 0.98 $\pm$ 0.12 | 1.04 $\pm$ 0.12 | 0.96 $\pm$ 0.11
76Ge | 0.96 $\pm$ 0.13 | 1.01 $\pm$ 0.14 | 0.94 $\pm$ 0.12
82Se | 0.99 $\pm$ 0.11 | 1.05 $\pm$ 0.12 | 0.97 $\pm$ 0.10
96Zr | 0.98 $\pm$ 0.12 | 1.04 $\pm$ 0.12 | 0.96 $\pm$ 0.11
100Mo | 0.96 $\pm$ 0.13 | 1.02 $\pm$ 0.13 | 0.94 $\pm$ 0.12
116Cd | 0.99 $\pm$ 0.11 | 1.06 $\pm$ 0.11 | 0.97 $\pm$ 0.10
124Xe | 0.89 $\pm$ 0.17 | 0.92 $\pm$ 0.18 | 0.87 $\pm$ 0.16
128Te | 1.07 $\pm$ 0.07 | 1.16 $\pm$ 0.07 | 1.05 $\pm$ 0.07
130Te | 1.03 $\pm$ 0.09 | 1.10 $\pm$ 0.09 | 1.01 $\pm$ 0.08
136Xe | 1.06 $\pm$ 0.08 | 1.14 $\pm$ 0.08 | 1.03 $\pm$ 0.07
150Nd | 1.06 $\pm$ 0.08 | 1.14 $\pm$ 0.08 | 1.03 $\pm$ 0.07

##  IV Results and discussion

By employing the framework of a relativistic theory for nuclear properties and
transitions as outlined in the previous sections, we have performed
calculations of the $2\nu\beta\beta$ decay matrix elements. In the first step,
the RH-BCS model [87] has been used to describe the ground state properties of
the initial and final nuclei involved in the decay. The RH-BCS model is
formulated in the harmonic oscillator basis and we restrict calculations to 20
oscillator shells both for protons and neutrons. The single-particle wave
functions and the corresponding occupation probabilities in the RH-BCS
quasiparticle basis are used in the REDF-QRPA to describe Gamow-Teller (GT)
transitions involved in the $2\nu\beta\beta$ decay. The REDF-QRPA calculations
are performed in two steps, for beta minus (plus) matrix elements for decays
from the initial (final) nucleus to the intermediate nucleus, which are then
used to calculate the 2$\nu\beta\beta$ decay NME. All the other quantities
that appear in the calculations, including the energies of the intermediate
states and the $Q$-values, are also taken self-consistently from the RH-BCS
and REDF-QRPA calculations. In order to assess the information on the model
dependence of the 2$\nu\beta\beta$ decay NME, three relativistic energy
density functionals are used in the study, including the density-dependent
meson-exchange interaction with the DD-ME2 parameterization [90], and density
dependent point coupling interactions DD-PC1 [96], and the more recently
established DD-PCX [97]. One of the open questions in the description of
double-beta decays is the problem of quenching of the axial-vector coupling
constant $g_{A}$, that has attracted attention in many recent studies, e.g.
see Refs. [7, 2]. As has already been mentioned above, $g_{A}=1.0$ is
systematically used in the present study.

Throughout this section we provide the absolute values of the $2\nu\beta\beta$
decay matrix elements, since only the squares of their absolute values have
physical significance. The REDF-QRPA calculations are first performed to set
the cut-off energy for the quasiparticle pairs that compose the configuration
space. To all relevant configurations, a convergence test is performed in
order to restrict the maximal two-quasiparticle energy with the condition that
the value of the NME converges with increasing energy. We illustrate the
convergence of the NME values in Fig. 1 for the case of 48Ca, where the matrix
elements are shown as a function of the maximal proton-neutron two-
quasiparticle energy $E_{2qp}^{max}$. The DD-ME2 effective interaction is used
for this demonstration. The Gamow-Teller (GT) transitions are considered. One
can observe that the NMEs converge with high accuracy at two-quasiparticle
energy $E_{2qp}^{max}\approx$ 100 MeV.

Figure 1: The NMEs for $2\nu\beta\beta$ decay based on Gamow-Teller
transitions, shown as a function of the maximal two-quasiparticle energy. The
DD-ME2 interaction is used in calculations. Matrix elements are dimensionless.

Insight into the contributions of various states to the final NMEs can be
obtained from the running sum, that includes the sum of all contributions to
the matrix element up to a specific maximal excitation energy in the
intermediate nucleus [108] which we denote as Eexc in the following
discussion.

Figures 2 and 3 show the respective running sums for the Gamow-Teller double-
beta decay transitions for 48Ca and 76Ge, displayed as functions of Eexc. The
sums are taken at two values of the $T=0$ pairing strength parameter,
$V_{0pp}=0$, representing a situation where isoscalar pairing vanishes, and
the optimal values $V_{0pp}=0.98$ and 0.96 for 48Ca and 76Ge, respectively.
One can observe that for 48Ca the hypothesis of single state dominance [109]
holds in the REDF-QRPA calculations, i.e., most of the contribution to the
NMEs mainly comes from a single low-lying state in the intermediate nucleus.
For 76Ge, many states contribute to $M_{GT}$ up to 20 MeV, and in some cases
destructive interference in their contributions to the NMEs is obtained.

Figure 2: The running sum of the GT NMEs for the $2\nu\beta\beta$ decay of
48Ca for the DD-ME2 interaction, shown as a function of maximal excitation
energy Eexc in the intermediate nucleus. The cases with and without $T=0$
pairing are shown separately.  Figure 3: The same as in Fig.2, but for 76Ge.

The same dependence is shown in Fig. 4 for other nuclei considered in this
work, from 82Se to 150Nd, using optimal values of $V_{0pp}$. We see that
generally the single-state dominance (SSD) hypothesis is approximatively
fulfilled, although in some cases more complicated structure is obtained due
to contributions from several states. Particularly noticeable are significant
cancellations when the nuclear matrix element is significantly lowered from
the value at zero isoscalar pairing due to a high value of V0pp. It is
interesting to note that presented results for some nuclei are in reasonable
agreement with recent experimental results, which show that the SSD hypothesis
holds for 100Mo [110], while for 82Se we obtain a more complex structure than
the SSD measured in Ref. [111]. Our results for the running sum in 136Xe are
at variance with the experimental results in Ref. [112]. In the same study, it
has been shown that the running sum is strongly model dependent, resulting in
considerable differences between the reported QRPA and shell model
calculations [112]. Our results also differ from other calculations of running
sums, including also contributions with a negative sign. They appear and
become more prominent with increasing V0pp, and cancel out to a large extent
an initially large NME, and the fact that we show the running sums at higher
optimal V0pp value explains in part their appearance.

Figure 4: Running sums for the GT NMEs for DD-ME2 interaction and optimal
values of $V_{0pp}$ from Tab. 2, for the $2\nu\beta\beta$ decays of
82Se–150Nd, shown as a function of maximal excitation energy Eexc in the
intermediate nucleus.

Before systematic implementation of the optimal values of $V_{0pp}$ from Tab.
2 in model calculations, it is interesting to explore the general dependence
and sensitivity of the $2\nu\beta\beta$ decay NMEs on $V_{0pp}$. In Fig. 5 the
dependence of the absolute values of nuclear matrix elements is shown for 48Ca
as a function of $V_{0pp}$. The black line in figure represents the value of
the nuclear matrix element deduced from the experimental data on
$2\nu\beta\beta$ half-lives [79, 113]. The grey band around this line denotes
an uncertainty at the level of 3$\sigma$. As one can see in Fig. 5, the NME
systematically decreases with $V_{0pp}$, and an overlap between the calculated
curve and the experimental limit for the NMEs is obtained.

Figure 5: The dependence of the GT NMEs for 2$\nu\beta\beta$ decay on the
isoscalar pairing strength $V_{0pp}$ for 48Ca, using DD-ME2 interaction, in
comparison to the result obtained from the experimental data on
$2\nu\beta\beta$ decay [79].

Here, the overlap prescription given in Eq. (7) is used. We have verified that
the result using the overlap prescription in Eq. (7) appears almost identical
to the one given in Eq. (8). Figure 6 shows the difference between the NMEs
for 48Ca obtained using these two methods for calculating the overlap, given
as a function of $V_{0pp}$. One can observe that not only the difference is
small, but it also reduces as one approaches the $V_{0pp}$ value deduced from
the experiment as the optimal one, as shown in Fig. 5. Therefore, in the
following investigation we consider only the results based on the overlap
prescription given in Eq. (7).

Figure 6: The difference between the NMEs calculated using two prescriptions
for the GT $2\nu\beta\beta$ decay matrix element, shown as a function of the
isoscalar pairing strength $V_{0pp}$ (see text for details). ’MGT,simpl.’
refers to the usual prescription for the overlap factor in Eq. (7), while
’MGT,Simkovic’ refers to the prescription in Eq.(8).

A quantity related to the NMEs that has received recent theoretical and
experimental interest [114, 115] is the function $C(r)$, representing the
contribution to the NME at a given internucleon distance $r_{12}\equiv r$. It
is related to the NME, denoted here as $\mathcal{M}$, as

$\int_{0}^{\infty}drC(r)=\mathcal{M}.$ (20)

Figure 7 shows the $C(r)$ function for 2$\nu\beta\beta$ decay of 48Ca,
obtained using the DD-ME2 interaction for the range of values of $T=0$ pairing
strength parameters $V_{0pp}$. The quantity is evaluated in the closure
approximation and is directly connected to the dimensionless NME. It needs to
be noted that, due to the computational resources necessary for the evaluation
of the $C(r)$ function, the results we report have been calculated at a lower
energy cutoff than the results for the NMEs in the present work, specifically
at 50 MeV. Comparing the $C(r)$ function to those calculated in Ref. [115] we
conclude that qualitatively the same shape is obtained, with a sharp peak
around 1 fm, converging to zero with an increase of $r$. This result is
consistent with the findings that this shape is universal for 2$\nu\beta\beta$
decay [116], although our model does not contain $S=0$ pairing in the $T=0$
channel, as well as $T=0$ channel in the ground state. The behaviour obtained
with increasing $T=0$ pairing strength $V_{0pp}$ in the residual QRPA
interaction is also consistent with previous studies [115], i.e., the $C(r)$
function reaches its highest values when pairing is not taken into account and
with the increase of the pairing strength the central peak becomes reduced. We
would note that the contribution of high inter-nucleon distances in many
nuclei is negative and serves to lower the NME; these negative contributions,
once again, become more noticeable with increasing V0pp.

Figure 7: The function $C(r_{12})$, showing the contributions to the NME at
nucleon separations $r_{12}$, for the Gamow-Teller $2\nu\beta\beta$ decay of
48Ca, for $V_{0pp}$ from 0.0 to 1.5.

The same function, for other nuclides considered in this work up to 150Nd,
evaluated using DD-ME2 interaction and optimal values of $V_{0pp}$, is shown
in Fig. 8. Similar general behavior of $C(r_{12})$ function is obtained for
all nuclei, that peaks strongly at low $r_{12}$ values. Where the NMEs are
close to zero, we can note significant cancellations as one integrates the
$C(r)$ function over the inter-nucleon distance r. 150Nd is an exception to
the general rule that the first and most prominent peak of the function is
near 1 fm, with a dominant peak at significantly higher radii. This suggests
that the results for 150Nd should be taken with less certainty than those for
other nuclides. However, we need to point out that in certain cases our C(r)
functions also show contributions at higher inter-nucleon distances than is
found in previous QRPA calculations. We note that the contribution of high
inter-nucleon distances in many nuclei is negative and serves to lower the
NME; these negative contributions, once again, become more noticeable with
increasing V0pp.

Figure 8: The function $C(r_{12})$, showing the contributions to the NME at
nucleon separations $r_{12}$, for the Gamow-Teller $2\nu\beta\beta$ decay of
76Ge-150Nd, evaluated at optimal values of $V_{0pp}$.

The GT transitions relevant for the $2\nu\beta\beta$ decay of 48Ca have been
studied experimentally in Ref. [117]. The GT- and GT+ strength distributions
in 48Sc have been measured by the ${}^{48}\text{Ca}(p,n)$ and
${}^{48}\text{Ti}(n,p)$ reactions, respectively. The integrated GT strengths
up to an excitation energy of 30 MeV in 48Sc obtained from $(p,n)$ and $(n,p)$
spectra amount $\text{B}(\text{GT}^{-})=15.3\pm 2.2$ and
$\text{B}(\text{GT}^{+})=2.8\pm 0.3$. The REDF-QRPA calculations with DD-ME2
interaction for the corresponding transitions result in
$\text{B}(\text{GT}^{-})=23.47$ and $\text{B}(\text{GT}^{+})=3.48$, thus the
experiment provides 65$\%$ of GT- strength and 80$\%$ of the GT+ strength
obtained in model calculations. Clearly, further experimental studies of GT
transitions are needed to provide more transition strength that is relevant
for double beta decays. The missing strength in measured GT spectra has been
confirmed in recent REDF-QRPA calculations for other nuclei in Ref. [67], and
studies going beyond the RPA level including couplings between single nucleon
and collective nuclear vibrations could not resolve the discrepancy between
the theoretical and experimental GT strengths [118].

In the following, the NMEs are investigated for the set of nuclides usually
considered in $2\nu\beta\beta$ decay studies because there are experimental
data available: 76Ge, 82Se, 96Zr, 100Mo, 116Cd, 128Te, 130Te, 136Xe, 150Nd,
and most recently 124Xe. Figure 9 shows the NMEs for the set of 11 nuclides
listed above, given as a function of the isoscalar pairing strength $V_{0pp}$.
The results extracted from the experimental data on $2\nu\beta\beta$ decay
half-lives are shown for comparison [113]. We conclude that the dependence of
the nuclear matrix elements on $V_{0pp}$ is qualitatively similar for all
nuclei considered, with the values of the GT based NMEs decreasing and their
slopes increasing, with increasing absolute values of $V_{0pp}$. However,
there is a significant variation in the values of $V_{0pp}$ needed to
reproduce the NMEs based on experimental data. A similar issue has been
observed in $\beta$-decay studies [119], indicating that some mass dependence
is neccessary for the optimization of the isoscalar pairing channel in the
residual interaction of the REDF-QRPA, as we have also discussed in Sec.
III.3.

Figure 9: The dependence of the NMEs for $2\nu\beta\beta$ decay on the
isoscalar pairing strength $V_{0pp}$ for $48<A<150$ nuclei that decay through
the $2\nu\beta\beta$ channel, in comparison to the result obtained from the
experimental data [79].

In order to assess the information about the model dependence of the
$2\nu\beta\beta$ decay matrix elements, we conducted a detailed comparison of
our results with a selection of previous studies. A direct comparison of the
results of the present study with those from previous calculations can be
difficult both due to different theoretical frameworks and parametrisations
used.

In Ref. [1], the pn-QRPA matrix elements are given for several $A>100$
nuclides, using a value of the $T=0$ pairing strength $g_{pp}$ = 0.7. A
comparison with the REDF-QRPA results with $V_{0pp}=$ 0 is given in Table 3.
We see that the upper limit of NMEs in the present study remains lower than
those reported by other pn-QRPA calculations.

Table 3: The NMEs for 2$\nu\beta\beta$ decay based on the REDF-QRPA (DD-ME2 interaction, $V_{0pp}=$ 0) and the NMEs based on the pn-QRPA from Ref. [1]. | REDF-QRPA | pn-QRPA [1]
---|---|---
| (DD-ME2) |
100Mo | 0.259 | 0.6560
116Cd | 0.056 | 0.2169
128Te | 0.078 | 0.1041
130Te | 0.079 | 0.1066

Further comparison can be made with the recent pn-QRPA results of Šimkovic et
al. [83], which is done in Table 4. We note that the values from Ref. [83],
except for 48Ca and 82Se, are in excellent agreement with our results, even
though the latter were obtained with vanishing isoscalar pairing. We also note
that the NMEs given in Ref. [83] are effective matrix elements obtained with a
quenched value of $g_{A}$ = 0.904. The renormalisation factor was split
according to isospin channel into $g_{pp}^{T=1}$ and $g_{pp}^{T=0}$; the
values of both parameters for each respective nucleus are given in Tab. 4.

Table 4: The NMEs for 2$\nu\beta\beta$ decay based on the REDF-QRPA (DD-ME2 interaction, $V_{0pp}=$ 0) and the NMEs based on the pn-QRPA from Ref. [83] with the pairing strength parameters $g_{pp}^{T=1}$ and $g_{pp}^{T=0}$ given in the last two columns. | REDF-QRPA | pn-QRPA [83] | g${}_{pp}^{T=1}$[83] | g${}_{pp}^{T=0}$[83]
---|---|---|---|---
| (DD-ME2) |  |  |
48Ca | 0.078 | 0.019 | 1.028 | 0.745
76Ge | 0.070 | 0.077 | 1.021 | 0.733
82Se | 0.041 | 0.071 | 1.016 | 0.737
96Zr | 0.140 | 0.162 | 0.961 | 0.739
100Mo | 0.259 | 0.306 | 0.985 | 0.799
116Cd | 0.056 | 0.059 | 0.892 | 0.877
128Te | 0.078 | 0.076 | 0.965 | 0.741
130Te | 0.079 | 0.065 | 0.963 | 0.737

Next, in Table 5 we compare the NMEs using the REDF-QRPA (DD-ME2 interaction,
$V_{0pp}=$ 0 MeV) with those of recent Interacting Boson Model (IBM)
calculations [36, 120]. The IBM results have been calculated within a closure
approximation, and to allow for a direct comparison, the values we report for
the NMEs have been obtained in the closure approximation as well. To obtain
NMEs comparable to our own from the NMEs quoted by Barea and Iachello, which
do not include the energy denominator, we divide the NMEs with an average
energy denominator as tabulated in [121]. The resulting NMEs based on the IBM
are generally larger, but comparable to the REDF-QRPA results. The inclusion
of isoscalar pairing in the REDF-QRPA residual interaction would further
increase this difference. This result suggests an advantage of our
calculations over the IBM, as the experimental values are even lower for most
nuclei considered, and, for most nuclei, are within the reach for the REDF-
QRPA with a suitable choice of $V_{0pp}$.

Further comparison is made with the Interacting Shell Model (ISM) [111], as
shown in Table 6. Since the ISM contains isoscalar pairing as a significant
effect, we compare the shell model results to those of the REDF-QRPA
calculated at the optimal pairing strength $V_{0pp}$ as given in Table 2. The
REDF-QRPA results are considerably smaller than those of the ISM.

Table 5: The NMEs for $2\nu\beta\beta$ decay based on the REDF-QRPA (DD-ME2 interaction, $V_{0pp}=$ 0), calculated in the closure approximation, and the results of the Interacting Boson Model (IBM) [120]. | REDF-QRPA | IBM [120]
---|---|---
| (DD-ME2) |
48Ca | 0.121 | 0.213
76Ge | 0.110 | 0.471
82Se | 0.071 | 0.356
96Zr | 0.125 | 0.208
100Mo | 0.210 | 0.272
116Cd | 0.056 | 0.197
128Te | 0.079 | 0.308
Table 6: The NMEs for $2\nu\beta\beta$ decay based on the REDF-QRPA (DD-ME2 interaction, optimal $V_{0pp}$), and the results of the Interacting Shell Model (ISM) (Ref. [111] unless otherwise noted).  | REDF-QRPA | ISM ([111])
---|---|---
| (DD-ME2) |
48Ca | 0.019 | 0.026
76Ge | 0.001 | 0.104
82Se | 0.001 | 0.109
124Xe[122] | 0.023 | 0.028 — 0.072
128Te [123] | 0.006 | 0.030
130Te | 0.002 | 0.061
136Xe [123] | 0.001 | 0.013

More recent PN-QRPA results can be found in Ref. [124], but only for some of
the nuclei considered in this work. The parameters for the calculation are
determined using two procedures, fitting of the ft values to an isobaric
triplet and multiplet. The comparison to the NMEs based on the REDF-QRPA with
optimal $V_{0pp}$ is given in Table 7, showing reasonable agreement for 100Mo,
while for 116Cd and 128Te the REDF-QRPA provides smaller NMEs.

Table 7: The NMEs for $2\nu\beta\beta$ decay based on the REDF-QRPA (DD-ME2 interaction, optimal $V_{0pp}$), and recent results of the proton-neutron QRPA (PN-QRPA) (Ref. [124].  | REDF-QRPA | PN-QRPA ([124]) | PN-QRPA ([124])
---|---|---|---
| (DD-ME2) | triplet | multiplet
100Mo | 0.189 | 0.153 | 0.131
116Cd | 0.038 | 0.153 | 0.160
128Te | 0.006 | 0.069 | 0.095

Finally, we compare our calculations with the results of a recent effective
theory (ET) treatment of double beta decay [125]. Effective theories concern
the behaviour of nuclei at suitably low energies, which is described in terms
of low-energy collective degrees of freedom, with any influence from high-
energy physics being encoded into the low-energy constants which are fitted to
available experimental data. The effective theory described in Ref. [125]
takes as its degrees of freedom nucleons and quadrupole phonon excitations of
an even-even spherical core the nucleons are coupled to. The comparison of
REDF-QRPA results using DD-ME2 interaction and optimal $V_{0pp}$, with those
of the ET, is given in Table 8. With the exception of 100Mo, our results are
systemically smaller, although they are comparable to the effective theory
result in the case of 116Cd.

Table 8: The NMEs for $2\nu\beta\beta$ decay based on the REDF-QRPA (DD-ME2 interaction, optimal $V_{0pp}$), and the results of the effective theory outlined in Ref [125].  | REDF-QRPA | ET ([125])
---|---|---
| (DD-ME2) |
76Ge | 0.001 | 0.054
82Se | 0.001 | 0.097
100Mo | 0.189 | 0.111
116Cd | 0.038 | 0.085
128Te | 0.006 | 0.031
130Te | 0.002 | 0.021

One of the open questions in the description of double-beta decays is
optimization of the strength parameter of the isoscalar pairing in the
residual QRPA interaction. Since this interaction channel is difficult to
constrain, and in most of the models cannot be determined based on the ground
state properties, it is necessary to explore its role in the NMEs for double-
beta decays. In the present study, the experimental data on $\beta$-decay half
lives are used to constrain the value of the isoscalar pairing strength
parameter $V_{0pp}$ (Sec. III.3), and in this way the REDF-QRPA can provide
predictions for $2\nu\beta\beta$ decay properties. Fig. 10 shows the NME for
$2\nu\beta\beta$ decay nuclides from 48Ca up to 150Nd, obtained using the
REDF-QRPA with the DD-ME2 interaction, for the range of values for $V_{0pp}$ =
0–1. The experimental data adopted from Ref. [79] are also shown for
comparison. One can observe systematic decrease of the NMEs with increasing
the value of $V_{0pp}$. For several nuclei, the experimentally determined NMEs
(with their error bars) are within the range of calculated NMEs when
considering the full range of $V_{0pp}$ values from 0 to 1. The exceptions are
76Ge, 82Se, 96Se, 116Cd, and 150Nd where the calculated NMEs are smaller (also
the experimental errors are rather small), and 96Zr where the calculated NMEs
are larger. In the case of 116Cd and 96Zr the differences between the
experimental and calculated NMEs within given range of $V_{0pp}$ are very
small.

In Table 9 we give the NMEs for the DD-ME2 interaction, using the optimal
value for the isoscalar pairing strength $V_{0pp}$, in comparison to the
values obtained from the experimental data on $2\nu\beta\beta$ decay [79]. The
inclusion of $T=0$ pairing further reduces the NMEs, thus allowing the
agreement with the experimental data where previously the calculated NMEs have
been too large. However, in many cases the isoscalar pairing reduced the NMEs
such that the final result is below the experimental range.

Figure 10: The dependence of the NME for $2\nu\beta\beta$ decay on the isoscalar pairing strength $V_{0pp}$ for $A=48$–150 nuclides. Table 9: The NMEs for $2\nu\beta\beta$ decay obtained with the DD-ME2 interaction and optimal values of $V_{0pp}$ given in Tab. 2, and NMEs obtained from the experimental data. | REDF-QRPA | exp. [126]
---|---|---
| (DD-ME2) |
48Ca | 0.019${}^{+0.014}_{-0.019}$ | 0.035 $\pm$ 0.003
76Ge | 0.001${}^{+0.024}_{-0.001}$ | 0.106 $\pm$ 0.004
82Se | 0.001${}^{+0.011}_{-0.001}$ | 0.085 $\pm$ 0.001
96Zr | 0.121${}^{+0.006}_{-0.010}$ | 0.080 $\pm$ 0.004
100Mo | 0.189${}^{+0.020}_{-0.033}$ | 0.185 $\pm$ 0.002
116Cd | 0.038${}^{+0.004}_{-0.004}$ | 0.108 $\pm$ 0.003
124Xe | 0.023${}^{+0.021}_{-0.023}$ | 0.059 $\pm$ 0.015
128Te | 0.006${}^{+0.011}_{-0.006}$ | 0.043 $\pm$ 0.003
130Te | 0.002${}^{+0.010}_{-0.002}$ | 0.0293 $\pm$ 0.0009
136Xe | 0.001${}^{+0.004}_{-0.001}$ | 0.0177 $\pm$ 0.0002
150Nd | 0.013${}^{+0.007}_{-0.005}$ | 0.055 $\pm$ 0.003

In order to assess the information about the systematic model dependence of
the NMEs for $2\nu\beta\beta$ decay in the REDF framework, we extend our
calculations also to the relativistic point coupling interactions (see Sec.
III). In Fig. 11 the NMEs are shown for 48Ca, obtained using point coupling
interactions DD-PC1 and DD-PCX, those with meson exchange interaction DD-ME2,
and the NME obtained from the experimental data on $2\nu\beta\beta$ decay
[79]. The figure shows the dependence of the NMEs on the $T=0$ pairing
strength $V_{0pp}$. The results for the two point coupling interactions
display some variations of the NMEs, though not considerable; the NMEs for the
DD-PC1 interaction are somewhat lower, e.g. at $V_{0pp}=0$ the difference in
the NMEs for DD-PC1 and DD-PCX interactions is less than $\approx 0.01$. The
NMEs for DD-ME interaction show qualitatively the same dependence on
$V_{0pp}$. Similar analysis of the sensitivity of the NMEs on the effective
interaction employed is performed for other nuclei of interest, from 76Ge
toward 128Te, as shown in Fig. 12. The NMEs for DD-PC1 interaction appear
systematically smaller than those of DD-PCX interaction, and their differences
provide the insight into the model uncertainties when using different point
coupling interactions. In the case of DD-ME2 interaction, the NMEs are smaller
than those of DD-PC1 and DD-PCX interactions, except for 76Ge.The differences
depend on the specific nucleus under consideration. For example, for 48Ca the
results for the two interactions almost coincide, while for 150Nd considerable
differences in the NMEs are obtained.

Figure 11: Comparison of the NMEs for 48Ca using the DD-PC1 and DD-PCX interactions and the value obtained from experimental data [126].  Figure 12: The dependence of the NMEs on the isoscalar pairing strength $V_{0pp}$ for the $2\nu\beta\beta$ decay for the set of nuclides in the mass range 76Ge–150Nd, for all relativistic interactions employed in this work shown in comparison to the values obtained from the experimental data [126].  Table 10: The same as Tab. 9 but for DD-PC1 and DD-PCX interactions and the corresponding optimal values of the isoscalar pairing strength $V_{0pp}$ given in Tab. 2. | REDF-QRPA | REDF-QRPA | exp. [79]
---|---|---|---
| (DD-PC1) | (DD-PCX) |
48Ca | 0.008${}^{+0.019}_{-0.008}$ | 0.024${}^{+0.016}_{-0.021}$ | 0.035 $\pm$ 0.003
76Ge | 0.003${}^{+0.021}_{-0.003}$ | 0.009${}^{+0.014}_{-0.009}$ | 0.106 $\pm$ 0.004
82Se | 0.009${}^{+0.019}_{-0.009}$ | 0.005${}^{+0.014}_{-0.005}$ | 0.085 $\pm$ 0.001
96Zr | 0.148${}^{+0.016}_{-0.027}$ | 0.160${}^{+0.007}_{-0.011}$ | 0.080 $\pm$ 0.004
100Mo | 0.188${}^{+0.028}_{-0.049}$ | 0.233${}^{+0.016}_{-0.026}$ | 0.185 $\pm$ 0.002
116Cd | 0.042${}^{+0.006}_{-0.008}$ | 0.087${}^{+0.005}_{-0.006}$ | 0.108 $\pm$ 0.003
124Xe | 0.034${}^{+0.036}_{-0.019}$ | 0.034${}^{+0.084}_{-0.034}$ | 0.059 $\pm$ 0.015
128Te | 0.009${}^{+0.015}_{-0.005}$ | 0.012${}^{+0.015}_{-0.012}$ | 0.043 $\pm$ 0.003
130Te | 0.014${}^{+0.015}_{-0.013}$ | 0.006${}^{+0.018}_{-0.006}$ | 0.0293 $\pm$ 0.0009
136Xe | 0.007${}^{+0.006}_{-0.005}$ | 0.003${}^{+0.008}_{-0.003}$ | 0.0177 $\pm$ 0.0002
150Nd | 0.023${}^{+0.010}_{-0.008}$ | 0.066${}^{+0.019}_{-0.022}$ | 0.055 $\pm$ 0.003

The optimal values for the $T=0$ pairing strength, obtained for the REDF-QRPA
with DD-PC1 and DD-PCX interactions (Table 2), can now be employed in
description of $2\nu\beta\beta$ decay NMEs. The calculated NMEs are summarized
in Table 10 for DD-PC1 and DD-PCX interactions, respectively, in comparison to
the experimental values. Although there are some variations in the NMEs when
compared to the experimental data, an overall reasonable qualitative
agreeement is obtained. This is illustrated in Fig. 13, where we summarize the
results of the present study, including the NMEs for $2\nu\beta\beta$ decay
obtained using the REDF-QRPA with DD-ME2, DD-PC1, and DD-PCX interactions with
the corresponding optimal values of isoscalar pairing strength parameters. For
comparison, the NMEs from previous studies are shown, including pn-QRPA
implementations by Suhonen [127], Pirinen [1], and Šimkovic [83], Interacting
Boson Model (IBM) [120] and the Interacting Shell Model (ISM) [111], as well
as the experimental result [79]. Clearly, the results of the REDF-QRPA are of
comparable quality than those from previous studies, though for most of nuclei
somewhat lower than those determined from the experimental data. In addition,
rather than providing just a single NME value, our study provides the insight
into systematic uncertainties due to variations in the formulation of the REDF
and parametrizations used. Additional uncertainty is accounted for due to the
isoscalar pairing strength interaction that is constrained by $\beta$-decay
half-lives.

Figure 13: Summary of the $2\nu\beta\beta$ decay NMEs from the REDF-QRPA with
DD-ME2, DD-PC1 and DD-PCX interactions with optimal $T=0$ pairing together
with respective $1\sigma$ uncertainties, compared to the calculations based on
the pn-QRPA by Deppisch and Suhonen [124], Pirinen [1], and Šimkovic [83],
Interacting Boson Model (IBM) [120], Interacting Shell Model (ISM) [111], and
the effective theory (ET) outlined in [125]. The NMEs from the experimental
data [79] are also shown.

The nuclear matrix elements calculated in this work appear rather small
compared to previous studies. As already discussed, the isoscalar and
isovector pairing strengths are very similar, indicating that the spin-isospin
SU(4) symmetry is softly broken. Weak breaking of the SU(4) symmetry could
partly explain small NMEs obtained within REDF-QRPA approach, but even for
$V_{0pp}=0$ the the NMEs appear smaller than in other approaches. Another
issue is that our treatment did not include relevant effects going beyond the
QRPA, as mentioned in the Sec.I.

In order to preform an additional sensitivity check of our calculations we
have also preformed calculations in which Q values derived from the
experimenal masses have been used. Figure 14 shows the NMEs calculated using
DD-ME2 interaction, both with experimental and calculated Q values. One can
observe that the NMEs results are nearly the same for 48Ca, 76Ge, and 82Se,
while for other nuclei the same trend with $V_{0pp}$ is obtained, but some
differences can be observed. In several cases the NMEs become smaller when
using experimental Q values instead of calculated ones. Thus, the choice of Q
value does not provide a solution to the problem of small NMEs values, except
for 96Zr. Whether the 2$\nu\beta\beta$ NMEs could be increased by an
appropriate treatment of the effects we could not include in the current work
is an interesting perspective for the further research. For example, in Ref.
[129] the structure of 76Ge and 76Se, has been studied with Gogny functionals,
indicating that the comparison with the experimental data could only be
obtained when triaxial shapes have been included, and these are precisely the
nuclei for which we obtain the NMEs that diverge the most from experimental
results. Also in the present study the systematic uncertainty due to
implementation of three different relativistic interactions is considerably
larger than the experimental uncertainty. Therefore, further investigations
and improvements of the theory frameworks are necessary in order to reduce the
theoretical uncertainty obtained from our, but also from other studies.

Figure 14: The NMEs for the DD-ME2 interaction, calculated with the Q-values
derived self-consistently in the RH-BCS framework, shown next to the NMEs for
the same interaction, but with the Q-values based on experimental data [130].

##  V Conclusion

In this work a theory framework is established for the study of
$2\nu\beta\beta$ decay nuclear matrix elements based on the relativistic
nuclear energy functional. Model calculations include two different
formulations of the effective interactions, density-dependent meson-exchange
and point coupling interactions, and three parameterizations (DD-ME2, DD-PC1,
DD-PCX) have been employed in order to assess the information on the
systematic uncertainties on the NMEs in the relativistic framework. The ground
states of nuclei involved in the decay are calculated within the relativistic
Hartree-BCS model, while nuclear transitions in the $2\nu\beta\beta$ decay are
described using the REDF-QRPA. In addition to the isovector pairing
correlations taken into account in the ground state calculations within the
RH-BCS model, the isoscalar paring channel has also been included in the
residual REDF-QRPA interaction, and its optimal values of the strength
parameter $V_{0pp}$ are constrained by the experimental data on $\beta$-decay
half lives, in order to allow predictions in double beta decay studies.

We have calculated the NMEs for a set of nuclei that undergo $2\nu\beta\beta$
decay, 48Ca, 76Ge, 82Se, 96Zr, 100Mo, 116Cd, 128Te, 130Te, 124Xe, 136Xe, and
150Nd. The dependence of the NMEs on the isoscalar paring strength $V_{0pp}$
has been investigated at the limit $V_{0pp}=0$ and using optimal $V_{0pp}$
values. The NMEs from the REDF-QRPA provide an improvement over the
interacting boson model, being closer to the values of the NMEs obtained from
the experimental data. However, when compared to the non-relativistic pn-QRPA,
interacting shell model, or effective theory, the NMEs from the present study
for most of studied nuclei, with a few exceptions, are rather small. While for
some studied nuclei the NMEs are already below experimental ones in the
$V_{0pp}=0$ limit, for several nuclei at this limit the NMEs are above
experimental values, but when introducing the optimal value of $V_{0pp}$, the
NMEs become rather low. It has been shown that different treatment of the Q
value calculation for most of studied nuclei could not increase the NMEs, and
further studies of additional effects are required to resolve this question.
Rather than providing just a single NME value like most of previous studies,
our work provides the insight into systematic uncertainties due to different
formulations of the REDF and parametrizations used. Additional uncertainty is
accounted for due to the isoscalar pairing strength interaction that is
constrained by $\beta$-decay half-lives.

This work provides an important benchmark for the future applications of the
relativistic framework in studies of neutrinoless double-beta decay. However,
we note that the present study represents our first study of double-beta
decays in the relativistic framework, and some effects have not been
considered, e.g., nuclear deformation, symmetry restoration, configuration
mixing, etc. Future improvements of the REDF based theory framework, in
particular the on-going development of the deformed REDF-QRPA for the
functionals used in this work, will allow additional improvements in modeling
nuclear double-beta decays.

##  VI Acknowledgements

We thank Deni Vale for support regarding the REDF-QRPA and useful discussions,
and Jenni Kotila for discussion on the phase space factors for double beta
decays. This work is supported by the QuantiXLie Centre of Excellence, a
project co-financed by the Croatian Government and European Union through the
European Regional Development Fund, the Competitiveness and Cohesion
Operational Programme (KK.01.1.1.01.0004).

## References

  * Pirinen and Suhonen [2015] P. Pirinen and J. Suhonen, Phys. Rev. C 91, 054309 (2015).
  * Engel and Menéndez [2017] J. Engel and J. Menéndez, Reports on Progress in Physics 80, 046301 (2017).
  * Jiao  _et al._ [2018] C. F. Jiao, M. Horoi, and A. Neacsu, Phys. Rev. C 98, 064324 (2018).
  * Haxton and Stephenson [1984] W. Haxton and G. Stephenson, Progress in Particle and Nuclear Physics 12, 409 (1984).
  * Bilenky and Giunti [2015] S. M. Bilenky and C. Giunti, International Journal of Modern Physics A 30, 1530001 (2015), https://doi.org/10.1142/S0217751X1530001X .
  * Faessler and Simkovic [1998] A. Faessler and F. Simkovic, Journal of Physics G: Nuclear and Particle Physics 24, 2139 (1998).
  * Mustonen and Engel [2013] M. T. Mustonen and J. Engel, Phys. Rev. C 87, 064302 (2013).
  * Staudt  _et al._ [1990] A. Staudt, T. Kuo, and H. Klapdor-Kleingrothaus, Physics Letters B 242, 17 (1990).
  * Primakoff and Rosen [1959] H. Primakoff and S. P. Rosen, Reports on Progress in Physics 22, 121 (1959).
  * Tomoda [1991] T. Tomoda, Reports on Progress in Physics 54, 53 (1991).
  * Suhonen and Civitarese [1998] J. Suhonen and O. Civitarese, Physics Reports 300, 123 (1998).
  * Elliott and Vogel [2002] S. R. Elliott and P. Vogel, Annual Review of Nuclear and Particle Science 52, 115 (2002), https://doi.org/10.1146/annurev.nucl.52.050102.090641 .
  * Avignone  _et al._ [2008] F. T. Avignone, S. R. Elliott, and J. Engel, Rev. Mod. Phys. 80, 481 (2008).
  * Vergados  _et al._ [2012] J. D. Vergados, H. Ejiri, and F. Šimkovic, Reports on Progress in Physics 75, 106301 (2012).
  * Saakyan [2013] R. Saakyan, Annual Review of Nuclear and Particle Science 63, 503 (2013), https://doi.org/10.1146/annurev-nucl-102711-094904 .
  * Vergados  _et al._ [2016] J. D. Vergados, H. Ejiri, and F. Šimkovic, International Journal of Modern Physics E 25, 1630007 (2016), https://doi.org/10.1142/S0218301316300071 .
  * Dolinski  _et al._ [2019] M. J. Dolinski, A. W. Poon, and W. Rodejohann, Annual Review of Nuclear and Particle Science 69, 219 (2019), https://doi.org/10.1146/annurev-nucl-101918-023407 .
  * Ejiri  _et al._ [2019] H. Ejiri, J. Suhonen, and K. Zuber, Physics Reports 797, 1 (2019), neutrino-nuclear responses for astro-neutrinos, single beta decays and double beta decays.
  * Suhonen and Civitarese [2012] J. Suhonen and O. Civitarese, Journal of Physics G: Nuclear and Particle Physics 39, 085105 (2012).
  * Terasaki and Iwata [2019] J. Terasaki and Y. Iwata, Phys. Rev. C 100, 034325 (2019).
  * Ejiri [2019] H. Ejiri, Frontiers in Physics 7, 30 (2019).
  * Simkovic  _et al._ [1997] F. Simkovic, J. Schwieger, G. Pantis, and A. Faessler, Foundations of Physics 27, 1275 (1997).
  * Šimkovic  _et al._ [2009] F. Šimkovic, A. Faessler, H. Müther, V. Rodin, and M. Stauf, Phys. Rev. C 79, 055501 (2009).
  * Rodin  _et al._ [2006] V. Rodin, A. Faessler, F. Šimkovic, and P. Vogel, Czechoslovak Journal of Physics 56, 495 (2006).
  * Stoica and Klapdor-Kleingrothaus [2001] S. Stoica and H. Klapdor-Kleingrothaus, Nuclear Physics A 694, 269 (2001).
  * Beneš  _et al._ [2006] P. Beneš, F. Šimkovic, A. Faessler, and W. Kamiński, Progress in Particle and Nuclear Physics 57, 257 (2006), international Workshop of Nuclear Physics 27th course.
  * Stoica and Klapdor-Kleingrothaus [2003] S. Stoica and H. V. Klapdor-Kleingrothaus, The European Physical Journal A - Hadrons and Nuclei 17, 529 (2003).
  * Rodríguez and Martínez-Pinedo [2010] T. R. Rodríguez and G. Martínez-Pinedo, Phys. Rev. Lett. 105, 252503 (2010).
  * Menéndez  _et al._ [2014] J. Menéndez, T. R. Rodríguez, G. Martínez-Pinedo, and A. Poves, Phys. Rev. C 90, 024311 (2014).
  * Ferreira  _et al._ [2020] V. d. S. Ferreira, A. R. Samana, F. Krmpotić, and M. Chiapparini, Phys. Rev. C 101, 044314 (2020).
  * Caurier  _et al._ [2005] E. Caurier, G. Martínez-Pinedo, F. Nowacki, A. Poves, and A. P. Zuker, Rev. Mod. Phys. 77, 427 (2005).
  * Kostensalo and Suhonen [2020] J. Kostensalo and J. Suhonen, Physics Letters B 802, 135192 (2020).
  * Suhonen and Kostensalo [2019] J. Suhonen and J. Kostensalo, Frontiers in Physics 7, 29 (2019).
  * Barea and Iachello [2009] J. Barea and F. Iachello, Phys. Rev. C 79, 044301 (2009).
  * Brown  _et al._ [2018] B. A. Brown, D. L. Fang, and M. Horoi, Journal of Physics: Conference Series 966, 012017 (2018).
  * Barea  _et al._ [2015] J. Barea, J. Kotila, and F. Iachello, Phys. Rev. C 91, 034304 (2015).
  * Jiao  _et al._ [2017] C. F. Jiao, J. Engel, and J. D. Holt, Phys. Rev. C 96, 054310 (2017).
  * Rath  _et al._ [2019] P. K. Rath, A. Kumar, R. Chandra, R. Gautam, P. K. Raina, and B. M. Dixit, International Journal of Modern Physics E 28, 1950096 (2019), https://doi.org/10.1142/S0218301319500964 .
  * Kotila  _et al._ [2010] J. Kotila, J. Suhonen, and D. S. Delion, Journal of Physics G: Nuclear and Particle Physics 37 (2010), 10.1088/0954-3899/37/1/015101.
  * Simkovic  _et al._ [1998] F. Simkovic, G. Pantis, and A. Faessler,  _Nonaccelerator new physics. Proceedings, 1st International Workshop, NANP’97, Dubna, Russia, July 7-11, 1997_ , Phys. Atom. Nucl. 61, 1218 (1998), [Yad. Fiz.61,1318(1998)], arXiv:nucl-th/9711060 [nucl-th] .
  * Ünlü  _et al._ [2018] S. Ünlü, N. Çakmak, and C. Selam, Nuclear Physics A 970, 379 (2018).
  * Suhonen [2011] J. Suhonen, Nuclear Physics. A 864 (2011), 10.1016/j.nuclphysa.2011.06.021.
  * Vogel and Zirnbauer [1986] P. Vogel and M. R. Zirnbauer, Phys. Rev. Lett. 57, 3148 (1986).
  * Engel  _et al._ [1988] J. Engel, P. Vogel, and M. R. Zirnbauer, Phys. Rev. C 37, 731 (1988).
  * Martínez-Pinedo  _et al._ [1999] G. Martínez-Pinedo, K. Langanke, and P. Vogel, Nuclear Physics A 651, 379 (1999).
  * Frauendorf and Macchiavelli [2014] S. Frauendorf and A. Macchiavelli, Progress in Particle and Nuclear Physics 78, 24 (2014).
  * Hinohara and Engel [2014] N. Hinohara and J. Engel, Phys. Rev. C 90, 031301 (2014).
  * Menéndez  _et al._ [2016] J. Menéndez, N. Hinohara, J. Engel, G. Martínez-Pinedo, and T. R. Rodríguez, Phys. Rev. C 93, 014305 (2016).
  * Song  _et al._ [2014] L. S. Song, J. M. Yao, P. Ring, and J. Meng, Phys. Rev. C 90, 054309 (2014).
  * Yao  _et al._ [2015] J. M. Yao, L. S. Song, K. Hagino, P. Ring, and J. Meng, Phys. Rev. C 91, 024316 (2015).
  * Song  _et al._ [2017] L. S. Song, J. M. Yao, P. Ring, and J. Meng, Phys. Rev. C 95, 024305 (2017).
  * Paar  _et al._ [2005] N. Paar, D. Vretenar, and P. Ring, Phys. Rev. Lett. 94, 182501 (2005).
  * Paar  _et al._ [2009] N. Paar, Y. F. Niu, D. Vretenar, and J. Meng, Phys. Rev. Lett. 103, 032502 (2009).
  * Niu  _et al._ [2013a] Y. F. Niu, Z. M. Niu, N. Paar, D. Vretenar, G. H. Wang, J. S. Bai, and J. Meng, Phys. Rev. C 88, 034308 (2013a).
  * Niu  _et al._ [2009] Y. Niu, N. Paar, D. Vretenar, and J. Meng, Physics Letters B 681, 315 (2009).
  * Khan  _et al._ [2011] E. Khan, N. Paar, and D. Vretenar, Phys. Rev. C 84, 051301 (2011).
  * Yüksel  _et al._ [2020] E. Yüksel, N. Paar, G. Colò, E. Khan, and Y. F. Niu, Phys. Rev. C 101, 044305 (2020).
  * Samana  _et al._ [2011] A. R. Samana, F. Krmpotić, N. Paar, and C. A. Bertulani, Phys. Rev. C 83, 024303 (2011).
  * Fantina  _et al._ [2012] A. F. Fantina, E. Khan, G. Colò, N. Paar, and D. Vretenar, Phys. Rev. C 86, 035805 (2012).
  * Vale  _et al._ [2016] D. Vale, T. Rauscher, and N. Paar, Journal of Cosmology and Astroparticle Physics 2016, 007 (2016).
  * Petković  _et al._ [2019] J. Petković, T. Marketin, G. Martínez-Pinedo, and N. Paar, Journal of Physics G: Nuclear and Particle Physics 46, 085103 (2019).
  * Roca-Maza and Paar [2018] X. Roca-Maza and N. Paar, Progress in Particle and Nuclear Physics 101, 96 (2018).
  * Mondal  _et al._ [2016] C. Mondal, B. K. Agrawal, M. Centelles, G. Colò, X. Roca-Maza, N. Paar, X. Viñas, S. K. Singh, and S. K. Patra, Phys. Rev. C 93, 064303 (2016).
  * Paar  _et al._ [2014] N. Paar, C. C. Moustakidis, T. Marketin, D. Vretenar, and G. A. Lalazissis, Phys. Rev. C 90, 011304 (2014).
  * Dobaczewski  _et al._ [2014] J. Dobaczewski, W. Nazarewicz, and P.-G. Reinhard, Journal of Physics G: Nuclear and Particle Physics 41, 074001 (2014).
  * Paar  _et al._ [2004] N. Paar, T. Nikšić, D. Vretenar, and P. Ring, Phys. Rev. C 69, 054303 (2004).
  * Vale  _et al._ [2021] D. Vale, Y. F. Niu, and N. Paar, Phys. Rev. C 103, 064307 (2021).
  * Ravlić  _et al._ [2021a] A. Ravlić, Y. F. Niu, T. Nikšić, N. Paar, and P. Ring, Phys. Rev. C 104, 064302 (2021a).
  * Fang  _et al._ [2011] D.-L. Fang, A. Faessler, V. Rodin, and F. Simkovic, Physical Review C 83 (2011), 10.1103/PhysRevC.83.034320.
  * Hyvarinen and Suhonen [2016] J. Hyvarinen and J. Suhonen, Advances in High Energy Physics 2016 (2016), 10.1155/2016/4714829.
  * Delion  _et al._ [2019] D. S. Delion, A. Dumitrescu, and J. Suhonen, Physical Review C 100 (2019), 10.1103/PhysRevC.100.024331.
  * Fang  _et al._ [2018] D.-L. Fang, A. Faessler, and F. Simkovic, Physical Review C 97 (2018), 10.1103/PhysRevC.97.045503.
  * Raduta and Raduta [2013] A. A. Raduta and C. M. Raduta, Journal of Physics: Conference Series 413, 012014 (2013).
  * Šimkovic  _et al._ [2018a] F. Šimkovic, R. Dvornický, D. c. v. Štefánik, and A. Faessler, Phys. Rev. C 97, 034315 (2018a).
  * Robledo  _et al._ [2018] L. M. Robledo, T. R. Rodríguez, and R. R. Rodríguez-Guzmán, Journal of Physics G: Nuclear and Particle Physics 46, 013001 (2018).
  * Šimkovic  _et al._ [2013] F. Šimkovic, V. Rodin, A. Faessler, and P. Vogel, Phys. Rev. C 87, 045501 (2013).
  * Stoica and Mirea [2019] S. Stoica and M. Mirea, Frontiers in Physics 7, 12 (2019).
  * Suhonen [2007] J. Suhonen,  _From Nucleons to Nucleus_, Theoretical and Mathematical Physics (Springer, Berlin, Germany, 2007).
  * Barabash [2015a] A. Barabash, Nuclear Physics A 935, 52 (2015a).
  * Marketin  _et al._ [2016a] T. Marketin, L. Huther, and G. Martínez-Pinedo, Phys. Rev. C 93, 025805 (2016a).
  * Ravlić  _et al._ [2021b] A. Ravlić, E. Yüksel, Y. F. Niu, and N. Paar, Phys. Rev. C 104, 054318 (2021b).
  * Ravlić  _et al._ [2020] A. Ravlić, E. Yüksel, Y. F. Niu, G. Colò, E. Khan, and N. Paar, Phys. Rev. C 102, 065804 (2020).
  * Šimkovic  _et al._ [2018b] F. Šimkovic, A. Smetana, and P. Vogel, Phys. Rev. C 98, 064325 (2018b).
  * de Oliveira  _et al._ [2015] L. de Oliveira, A. R. Samana, F. Krmpotic, A. E. Mariano, and C. A. Barbero, Journal of Physics: Conference Series 630, 012048 (2015).
  * Vogel [2012] P. Vogel, Journal of Physics G: Nuclear and Particle Physics 39, 124002 (2012).
  * Aunola and Suhonen [1996] M. Aunola and J. Suhonen, Nuclear Physics A 602, 133 (1996).
  * Nikšić  _et al._ [2014] T. Nikšić, N. Paar, D. Vretenar, and P. Ring, Computer Physics Communications 185, 1808 (2014).
  * Kohn and Sham [1965] W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
  * Kohn [1999] W. Kohn, Rev. Mod. Phys. 71, 1253 (1999).
  * Lalazissis  _et al._ [2005] G. A. Lalazissis, T. Nikšić, D. Vretenar, and P. Ring, Phys. Rev. C 71, 024312 (2005).
  * Hofmann  _et al._ [2001] F. Hofmann, C. M. Keil, and H. Lenske, Phys. Rev. C 64, 034314 (2001).
  * Typel and Wolter [1999] S. Typel and H. Wolter, Nuclear Physics A 656, 331 (1999).
  * Fuchs  _et al._ [1995] C. Fuchs, H. Lenske, and H. H. Wolter, Phys. Rev. C 52, 3043 (1995).
  * Paar  _et al._ [2007] N. Paar, D. Vretenar, E. Khan, and G. Colò, Reports on Progress in Physics 70, 691 (2007).
  * Paar  _et al._ [2015] N. Paar, T. Marketin, D. Vale, and D. Vretenar, International Journal of Modern Physics E 24, 1541004 (2015), https://doi.org/10.1142/S0218301315410049 .
  * Nikšić  _et al._ [2008] T. Nikšić, D. Vretenar, and P. Ring, Phys. Rev. C 78, 034318 (2008).
  * Yüksel  _et al._ [2019] E. Yüksel, T. Marketin, and N. Paar, Phys. Rev. C 99, 034318 (2019).
  * Bonche  _et al._ [2005] P. Bonche, H. Flocard, and P. Heenen, Computer Physics Communications 171, 49 (2005).
  * Reinhard  _et al._ [2021] P.-G. Reinhard, B. Schuetrumpf, and J. Maruhn, Computer Physics Communications 258, 107603 (2021).
  * Tian  _et al._ [2009a] Y. Tian, Z.-y. Ma, and P. Ring, Phys. Rev. C 80, 024313 (2009a).
  * Bai  _et al._ [2014] C. L. Bai, H. Sagawa, G. Colò, Y. Fujita, H. Q. Zhang, X. Z. Zhang, and F. R. Xu, Phys. Rev. C 90, 054335 (2014).
  * Toivanen and Suhonen [1995] J. Toivanen and J. Suhonen, Phys. Rev. Lett. 75, 410 (1995).
  * Paar  _et al._ [2003] N. Paar, P. Ring, T. Nikšić, and D. Vretenar, Phys. Rev. C 67, 034312 (2003).
  * Daoutidis and Ring [2011] I. Daoutidis and P. Ring, Phys. Rev. C 83, 044303 (2011).
  * Nikšić  _et al._ [2013] T. Nikšić, N. Kralj, T. Tutiš, D. Vretenar, and P. Ring, Phys. Rev. C 88, 044327 (2013).
  * Tian  _et al._ [2009b] Y. Tian, Z.-y. Ma, and P. Ring, Phys. Rev. C 79, 064301 (2009b).
  * Niu  _et al._ [2013b] Z. Niu, Y. Niu, H. Liang, W. Long, T. Nikšić, D. Vretenar, and J. Meng, Physics Letters B 723, 172 (2013b).
  * de Guerra  _et al._ [2012] E. M. de Guerra, O. Moreno, P. Sarriguren, and M. Ramon, Journal of Physics: Conference Series 366, 012011 (2012).
  * Simkovic  _et al._ [2001] F. Simkovic, P. Domin, and S. V. Semenov, Journal of Physics G: Nuclear and Particle Physics 27, 2233 (2001).
  * Azzolini  _et al._ [2019] O. Azzolini, J. W. Beeman, F. Bellini, M. Beretta, M. Biassoni, C. Brofferio, C. Bucci, S. Capelli, L. Cardani, P. Carniti, N. Casali, D. Chiesa, M. Clemenza, O. Cremonesi, A. Cruciani, I. Dafinei, S. Di Domizio, F. Ferroni, L. Gironi, A. Giuliani, P. Gorla, C. Gotti, G. Keppel, J. Kotila, M. Martinez, S. Nagorny, M. Nastasi, S. Nisi, C. Nones, D. Orlandi, L. Pagnanini, M. Pallavicini, L. Pattavina, M. Pavan, G. Pessina, V. Pettinacci, S. Pirro, S. Pozzi, E. Previtali, A. Puiu, C. Rusconi, K. Schäffner, C. Tomei, M. Vignati, and A. Zolotarova, Phys. Rev. Lett. 123, 262501 (2019).
  * Coraggio  _et al._ [2019] L. Coraggio, L. De Angelis, T. Fukui, A. Gargano, N. Itaco, and F. Nowacki, Phys. Rev. C 100, 014316 (2019).
  * Gando  _et al._ [2019] A. Gando, Y. Gando, T. Hachiya, M. Ha Minh, S. Hayashida, Y. Honda, K. Hosokawa, H. Ikeda, K. Inoue, K. Ishidoshiro, Y. Kamei, K. Kamizawa, T. Kinoshita, M. Koga, S. Matsuda, T. Mitsui, K. Nakamura, A. Ono, N. Ota, S. Otsuka, H. Ozaki, Y. Shibukawa, I. Shimizu, Y. Shirahata, J. Shirai, T. Sato, K. Soma, A. Suzuki, A. Takeuchi, K. Tamae, K. Ueshima, H. Watanabe, D. Chernyak, A. Kozlov, S. Obara, S. Yoshida, Y. Takemoto, S. Umehara, K. Fushimi, S. Hirata, B. E. Berger, B. K. Fujikawa, J. G. Learned, J. Maricic, L. A. Winslow, Y. Efremenko, H. J. Karwowski, D. M. Markoff, W. Tornow, T. O’Donnell, J. A. Detwiler, S. Enomoto, M. P. Decowski, J. Menéndez, R. Dvornický, and F. Šimkovic (KamLAND-Zen Collaboration), Phys. Rev. Lett. 122, 192501 (2019).
  * Barabash [2015b] A. Barabash, Physics Procedia 74, 416 (2015b), fundamental Research in Particle Physics and Cosmophysics.
  * Shimizu  _et al._ [2018] N. Shimizu, J. Menéndez, and K. Yako, Phys. Rev. Lett. 120, 142502 (2018).
  * Šimkovic  _et al._ [2008] F. Šimkovic, A. Faessler, V. Rodin, P. Vogel, and J. Engel, Phys. Rev. C 77, 045503 (2008).
  * Menéndez  _et al._ [2009] J. Menéndez, A. Poves, E. Caurier, and F. Nowacki, Nuclear Physics A 818, 139 (2009).
  * Yako  _et al._ [2009] K. Yako, M. Sasano, K. Miki, H. Sakai, M. Dozono, D. Frekers, M. B. Greenfield, K. Hatanaka, E. Ihara, M. Kato, T. Kawabata, H. Kuboki, Y. Maeda, H. Matsubara, K. Muto, S. Noji, H. Okamura, T. H. Okabe, S. Sakaguchi, Y. Sakemi, Y. Sasamoto, K. Sekiguchi, Y. Shimizu, K. Suda, Y. Tameshige, A. Tamii, T. Uesaka, T. Wakasa, and H. Zheng, Phys. Rev. Lett. 103, 012503 (2009).
  * Niu  _et al._ [2014] Y. F. Niu, G. Colò, and E. Vigezzi, Phys. Rev. C 90, 054328 (2014).
  * Marketin  _et al._ [2016b] T. Marketin, L. Huther, and G. Martínez-Pinedo, Phys. Rev. C 93, 025805 (2016b).
  * Iachello  _et al._ [2015] F. Iachello, J. Barea, and J. Kotila, Nuclear and Particle Physics Proceedings 265-266, 25 (2015), proceedings of the Neutrino Oscillation Workshop.
  * Kotila and Iachello [2012] J. Kotila and F. Iachello, Phys. Rev. C 85, 034316 (2012).
  * Coello Pérez  _et al._ [2019] E. Coello Pérez, J. Menéndez, and A. Schwenk, Physics Letters B 797, 134885 (2019).
  * Caurier  _et al._ [2012] E. Caurier, F. Nowacki, and A. Poves, Physics Letters B 711, 62 (2012).
  * Deppisch and Suhonen [2016] F. F. Deppisch and J. Suhonen, Phys. Rev. C 94, 055501 (2016).
  * Coello Pérez  _et al._ [2018] E. A. Coello Pérez, J. Menéndez, and A. Schwenk, Phys. Rev. C 98, 045501 (2018).
  * Barabash [2020] A. Barabash, Universe 6 (2020), 10.3390/universe6100159.
  * Suhonen [2005] J. Suhonen, Physics Letters B 607, 87 (2005).
  * Barea  _et al._ [2013] J. Barea, J. Kotila, and F. Iachello, Phys. Rev. C 87, 014315 (2013).
  * Rodríguez [2017] T. R. Rodríguez, Journal of Physics G: Nuclear and Particle Physics 44, 034002 (2017).
  * Wang  _et al._ [2021] M. Wang, W. Huang, F. Kondev, G. Audi, and S. Naimi, Chinese Physics C 45, 030003 (2021).