Isoscalar pairing interaction for the quasiparticle random-phase approximation approach to double-$\beta$ and $\beta$ decays

TL;DR

The paper refines a method to accurately determine nuclear matrix elements for neutrinoless double-beta decay, reducing uncertainties in key parameters, and applies it to specific isotopes with mixed success.

Contribution

It introduces a new approach to precisely determine the isoscalar pairing interaction and axial-vector coupling, improving the reliability of double-beta decay calculations.

Findings

Successful application to $^{136}$Xe with self-check validation.

Inability to accurately predict for $^{130}$Te.

Satisfactory results for $eta$ decay of $^{138}$Xe.

Abstract

We have proposed in a series of previous papers a method to determine the effective axial-vector current coupling and the strength of the isoscalar proton-neutron pairing interaction for calculating the nuclear matrix elements of the neutrinoless double-$\beta$ decay by the quasiparticle random-phase approximation. The combination of these two parameters have had an uncertainty in this approach, but now this uncertainty is removed. In this paper, we apply our method to the neutrinoless double-$\beta$ decays of $^{136}$Xe and $^{130}$Te and predict the nuclear matrix elements and reduced half-lives. Our calculation is tested first by a self-check method using the two-neutrino double-$\beta$ decay, and this test ensures the application of our method to $^{136}$Xe. It turns out, however, that our method is not successful in $^{130}$Te. Further test is made for our calculation of the $\beta$ decay of $^{138}$Xe, and a satisfactory result is obtained.

Tables

Table 1. Table 1: Properties of the HFB solutions of 136 Xe, 136 Ba, 130 Te, and 130 Xe. Δ p subscript Δ 𝑝 \Delta_{p} and Δ n subscript Δ 𝑛 \Delta_{n} are the average pairing gaps of the protons and neutrons, respectively, and β 𝛽 \beta denotes the quadrupole deformation parameter associated with the quadrupole moment.

Nucleus	${\Delta }_p$ (MeV)	${\Delta }_n$ (MeV)	$\beta$
136Xe	$1.341$	$0$	$0$
136Ba	$1.641$	$1.158$	$0$
130Te	$1.442$	$1.359$	$0$
130Xe	$1.692$	$1.438$	$0.112$

Table 2. Table 2: Strengths of proton-proton ( G p ​ p subscript 𝐺 𝑝 𝑝 G_{pp} ), neutron-neutron ( G n ​ n subscript 𝐺 𝑛 𝑛 G_{nn} ), isovector proton-neutron ( G p ​ n IV superscript subscript 𝐺 𝑝 𝑛 IV G_{pn}^{\mathrm{IV}} ), and isoscalar proton-neutron ( G p ​ n IS superscript subscript 𝐺 𝑝 𝑛 IS G_{pn}^{\mathrm{IS}} ) pairing interactions. The unit is MeV fm 3 for all strengths.

Nucleus	$G_{pp}$	$G_{nn}$	$G_{pn}^{\mathrm{IV}}$	$G_{pn}^{\mathrm{IS}}$
136Xe	$-194.3$	$-179.9$	$-187.1$	$-55.0$
136Ba	$-200.5$	$-189.3$	$-194.9$	$-55.0$
130Te	$-219.8$	$-179.9$	$-199.85$	$-50.0$
130Xe	$-219.8$	$-181.8$	$-200.8$	$-50.0$

Table 3. Table 3: Choice of set of intermediate-state energies, NMEs of the 2 ​ ν ​ β ​ β 2 𝜈 𝛽 𝛽 2\nu\beta\beta decay, its components, g A subscript 𝑔 𝐴 g_{A} , the calculated half-life in 2 ​ ν ​ β ​ β 2 𝜈 𝛽 𝛽 2\nu\beta\beta decay T 1 / 2 ( 2 ​ ν ) superscript subscript 𝑇 1 2 2 𝜈 T_{1/2}^{(2\nu)} , and the corresponding experimental half-life T 1 / 2 ( 2 ​ ν ) superscript subscript 𝑇 1 2 2 𝜈 T_{1/2}^{(2\nu)} (exp) Barabash ( 2013 ) . M ( 2 ​ ν ) superscript 𝑀 2 𝜈 M^{(2\nu)} implies M ( 2 ​ ν ) ​ ( I ) superscript 𝑀 2 𝜈 𝐼 M^{(2\nu)}(I) or M ( 2 ​ ν ) ​ ( F ) superscript 𝑀 2 𝜈 𝐹 M^{(2\nu)}(F) , and the same abbreviation is applied to other NMEs and T 1 / 2 ( 2 ​ ν ) superscript subscript 𝑇 1 2 2 𝜈 T_{1/2}^{(2\nu)} .

Initial nucleus	Intermediate-	$M^{\left(2\nu \right)}$	$\frac{M_{GT}^{\left(2\nu \right)}}{{\mu }_0}$	$\frac{M_F^{\left(2\nu \right)}}{{\mu }_{0F}}$	$g_A$	$T_{1/2}^{\left(2\nu \right)}$	$T_{1/2}^{\left(2\nu \right)}$(exp)
	state energy					(${10}^{21}$yr)	(${10}^{21}$yr)
136Xe	$E_{aK}^I$	$0.0748$	$0.0673$	$-0.00181$	$0.49$	$2.164$	$2.20\pm 0.06$
	$E_{aK}^F$	$0.0744$	$0.0666$	$-0.00188$	$0.49$	$2.188$
130Te	$E_{aK}^I$	$0.133$	$0.132$	$-0.00029$	$0.48$	$0.697$	$0.69\pm 0.13$
	$E_{aK}^F$	$0.0864$	$0.0856$	$-0.00028$	$0.60$	$0.677$

Table 4. Table 4: NMEs and R 1 / 2 ( 0 ​ ν ) subscript superscript 𝑅 0 𝜈 1 2 R^{(0\nu)}_{1/2} of the 0 ​ ν ​ β ​ β 0 𝜈 𝛽 𝛽 0\nu\beta\beta decay of 136 Xe and 130 Te with the g A subscript 𝑔 𝐴 g_{A} obtained by fitting the experimental T 1 / 2 ( 2 ​ ν ) subscript superscript 𝑇 2 𝜈 1 2 T^{(2\nu)}_{1/2} .

Initial	$M^{\left(0\nu \right)}$	$M_{GT}^{\left(0\nu \right)}$	$M_F^{\left(0\nu \right)}$	$g_A$	$R_{1/2}^{\left(0\nu \right)}$
nucleus					(${10}^{13}$MeV${}^2$yr)
136Xe	$5.040$	$3.094$	$-0.467$	$0.49$	$1.223$
130Te	$6.691$	$3.613$	$-0.709$	$0.48$	$0.772$
	$5.400$			$0.63$	$0.400$

Table 5. Table 5: M ( 0 ​ ν ) superscript 𝑀 0 𝜈 M^{(0\nu)} and R 1 / 2 ( 0 ​ ν ) superscript subscript 𝑅 1 2 0 𝜈 R_{1/2}^{(0\nu)} for g A = 1.27 subscript 𝑔 𝐴 1.27 g_{A}=1.27 . The closure approximation is used for M ( 0 ​ ν ) superscript 𝑀 0 𝜈 M^{(0\nu)} , thus we have only one NME and R 1 / 2 ( 0 ​ ν ) superscript subscript 𝑅 1 2 0 𝜈 R_{1/2}^{(0\nu)} for 130 Te.

Initial	$M^{\left(0\nu \right)}$	$R_{1/2}^{\left(0\nu \right)}$
nucleus		(${10}^{13}$MeV${}^2$yr)
136Xe	$3.109$	$0.075$
130Te	$3.725$	$0.054$

Equations

$$M^{(2\nu )}(I)$$

$$M^{(2\nu )}(F)$$

$$\frac{M_{GT}^{(2\nu )}(I)}{{\mu }_0}$$

$$\frac{M_{GT}^{(2\nu )}(F)}{{\mu }_0}$$

$$\frac{M_F^{(2\nu )}(I)}{{\mu }_{0F}}$$

$$\frac{M_F^{(2\nu )}(F)}{{\mu }_{0F}}$$

$${\mu }_a(I)$$

$${\mu }_a(F)$$

$$M^{(2\nu )}(I)\approx M^{(2\nu )}(F).$$

$$T_{1/2}^{(2\nu )}$$

$$M^{(0\nu )}$$

$$M_{GT}^{(0\nu )}$$

$$M_F^{(0\nu )}$$

$$V_{p{p}^{\prime },n{n}^{\prime }}^{GT(0\nu )}({\bar{E}}_a)$$

$$V_{p{p}^{\prime },n{n}^{\prime }}^{F(0\nu )}({\bar{E}}_a)=⟨p{p}^{\prime }|h_{+}(r_{12},{\bar{E}}_a){\tau }^{-}(1){\tau }^{-}(2)|n{n}^{\prime }⟩.$$

$$h_{+}(r_{12},{\bar{E}}_a)$$

$${\bar{\mu }}_a$$

$$\mathrm{si}(x)$$

$$\mathrm{ci}(x)$$

$$R_{1/2}^{(0\nu )}=\frac{m_e^2{c}^4}{G_{0\nu }{g}_A^4|M^{(0\nu )}{|}^2},$$

$$⟨m_{\nu }{⟩}^2=\frac{R_{1/2}^{(0\nu )}}{T_{1/2}^{(0\nu )}},$$

$$85.145-1.138=84.007,$$

$$B_{\mathrm{GT}}=g_A^2\sum _{K=0,\pm 1}|M_{\mathrm{GT}}^K{|}^2,$$

$$\frac{M_{GT\text{-}3}^{(2\nu )}}{{\mu }_0}$$

$${\xi }_{31}^{(2\nu )}=\frac{M_{GT\text{-}3}^{(2\nu )}}{M_{GT}^{(2\nu )}},$$

$$O_{\mathrm{id}}^{-}(0^{+})$$

$$O_{\mathrm{id}}^{-}(0^{-})$$

$$O_{\mathrm{Jid}}^{+}(J^{\pi })$$

$$O_{\mathrm{Jid}}^{+}(0^{-})$$

$$O_{\mathrm{Jid}}^{+}(J^{\pi })$$

Full text

I Introduction

Neutrino physics is in an important era. The neutrino-oscillation experiments have significant progress in terms of precision and accuracy, and nowadays there is a realistic possibility that the hierarchy of the neutrino mass is clarified in the future Neu . The experiments of the neutrinoless double-$\beta$ ($0\nu\beta\beta$) decay also progress, and the upper limit of the effective neutrino mass (Majorana neutrino mass) is lowered significantly compared to many years ago, e.g., Gando et al. (2016).

The nuclear matrix element (NME) of the $0\nu\beta\beta$ decay plays a crucial role for the determination of the effective neutrino mass Engel and Menéndez (2017). This NME and the phase-space factor are necessary for the determination, and these theoretical quantities cannot be confirmed directly experimentally. It seems fair to write that the calculation of the NME is more difficult than the phase-space factor because all of the candidate nuclei of the $0\nu\beta\beta$ decay are heavy, thus, approximation is essential for obtaining the nuclear wave functions. For this reason, an important problem has been how the reliability of the NME calculation can be shown, and if possible, how the calculation can be improved.

Several theoretical methods or models have been applied for obtaining the nuclear wave functions to calculate the NMEs, and those NMEs are in the range of a factor of 2$-$3 depending on the method Engel and Menéndez (2017). The quasiparticle random-phase approximation (QRPA) approach has a long history of the application to the $0\nu\beta\beta$ decay, e.g., Engel et al. (1988); Suhonen et al. (1990, 1991); Suhonen (1993); Suhonen et al. (1997); Suhonen and Civitarese (1998); Faessler and Šimkovic (1998); Šimkovic et al. (1998, 1999); Bobyk et al. (2001); Pacearescu et al. (2003); Šimkovic et al. (2004); Álvarez-Rodríguez et al. (2004); Suhonen (2005); Rodin et al. (2006); Kortelainen and Suhonen (2007); Šimkovic et al. (2008); Moreno et al. (2009); Šimkovic et al. (2009); Fang et al. (2011); Šimkovic et al. (2013); Mustonen and Engel (2013); Fang et al. (2015); Hyvärinen and Suhonen (2015); Suhonen (2017a). Those applications clarified that a crucial point in the calculation procedure is to determine the effective axial-vector current coupling $g_{A}$ and the strength of the isoscalar, necessarily proton-neutron, pairing interaction, e.g., Suhonen (2017a). Currently phenomenological value is necessary for the former parameter in any approach to the $\beta\beta$ decay. The value of the latter parameter is not well-established in the approaches using phenomenological energy-density functional because the proton-neutron pairing gap is not clear in the experimental data. It was pointed out Suhonen (2017a) that the different combinations of the two parameters result in similar NMEs as long as the experimental half-lives in the two-neutrino double-$\beta$ ($2\nu\beta\beta$) decays are reproduced. It would be even better, if this uncertainty of the combination were removed; the calculation of the NME is a matter of reliability.

We have proposed Terasaki (2016) a method to solve this problem by introducing a mathematical identity for determining the strength of the isoscalar pairing interaction; no additional symmetry is imposed. The effective $g_{A}$ is determined so as to reproduce the experimental $2\nu\beta\beta$ half-life; this is as usual in the QRPA approach, e.g., Šimkovic et al. (2013). This method is clear-cut theoretically, however, the obtained value of $g_{A}$ is $\simeq$ 0.5 much smaller than the usual values $\simeq$ 1.0, e.g., Šimkovic et al. (2013). It would not be surprising, if the appropriate effective $g_{A}$ depends on the approximation method. Considering the current situation of many $\beta\beta$ studies Engel and Menéndez (2017), it is necessary in our method to examine if our small $g_{A}$ does not cause any problem. A useful check is the consistency with the $\beta$ decay. As a matter of course, any of the candidate nuclei of the $0\nu\beta\beta$ decay used by the experiments do not have the $\beta$ decay, thus to our knowledge, there has been no theoretical paper on the $\beta\beta$ decay discussing the $\beta$ decay of the parent nucleus simultaneously. In this case, one can study the $\beta$ decay of nuclei close to the candidate nuclei in the nuclear chart assuming that the change in the effective $g_{A}$ is negligible.

The QRPA is a method to obtain the transition from the ground state to the excited states, e.g., Ring and Schuck (1980); Blaizot and Ripka (1985). When the pairs of the proton and neutron quasiparticles are used as the building blocks to express the transition, i.e., the proton-neutron QRPA (pnQRPA), the intermediate states of the $\beta\beta$ decay can be represented in two ways; one way is the pnQRPA based on the ground state of the mother nucleus, and another way is to use that of the grand-daughter nucleus. This structure of the dual intermediate-state spaces is a unique feature of the QRPA approach to the $\beta\beta$ decay. Because the QRPA is an approximation, the intermediate state of one space is not identical to a state of another space. We have proposed Terasaki (2018) a simple method of self-check of the validity of the dual intermediate-state spaces by comparing the $2\nu\beta\beta$ NME with the energies of the intermediate states obtained from the ground state of the mother nucleus and that with those energies obtained from the grand-daughter nucleus. The two results are close, if the QRPA is a good approximation. We do this test for 136Xe $\rightarrow$ 136Ba and 130Te $\rightarrow$ 130Xe in this paper. In addition, we compare the spectra of the intermediate states obtained by the two pnQRPA calculations. This comparison shows how the two intermediate spaces are close. The purpose of this paper is to show the results of the above tests of the methodology and discuss the validity of our approach.

Section II is a concise description of the calculation methods of the Hartree-Fock-Bogoliubov (HFB) approximation, QRPA, and NMEs of the $\beta\beta$ decays. The calculation results of the NMEs are shown in Sec. III, and several tests of the calculation are discussed. Section IV is the summary.

II Calculation method

II.1 HFB calculation

We proceed according to the original QRPA theory Ring and Schuck (1980); Blaizot and Ripka (1985). That is, the HFB calculation is performed at the beginning, and the subsequent QRPA calculation is performed using the HFB ground state and the same Hamiltonian as used in the HFB calculation. The HFB code developed according to Refs. Terán et al. (2003); Blazkiewicz et al. (2005); Oberacker et al. (2007) is used. The quasiparticle wave functions are represented in the two-dimensional cylindrical B-spline mesh with the vanishing boundary condition at the edge of the cylindrical box. The parity and the $z$ component of the angular momentum $j_{z}$ are good quantum numbers. Only $z>0$ is treated numerically, and the calculation with respect to the angle around the $z$ axis is processed analytically. These two symmetries are always conserved throughout our calculations below. The maximum $z$ is 20 fm, and the maximum of the radius variable $\rho$ of the direction perpendicular to the $z$ direction is also 20 fm. The number of the B-spline mesh points is 42 for each coordinate variable with a non-uniform distribution.

For the Hamiltonian, we use the Skyrme interaction (energy density functional) with the parameter set SkM∗ Bartel et al. (1982) and the volume (density independent) contact pairing interactions. The effective single-particle energy range of applying the pairing interaction is from the lowest level up to 60 MeV above the continuum threshold. The nuclear densities and pairing tensors are calculated in this range. The strengths of the pairing interactions for the like-particles are determined so as to reproduce by the HFB calculation the pairing gaps obtained from the experimental masses through the three-point formula Bohr and Mottelson (1969). The properties of the HFB ground-state solutions of 136Xe, 136Ba, 130Te, and 130Xe are summarized in Table 1.

II.2 QRPA calculation

After the HFB solution is obtained, the canonical single-particle Ring and Schuck (1980) wave functions are obtained by diagonalization of the density. The QRPA equation for the computation is constructed in the so-called matrix formulation Ring and Schuck (1980) with the canonical-quasiparticle basis. The dimension of the two-quasiparticle space to define the QRPA equation is truncated by using parameters based on the occupation probabilities of the canonical single-particle states Terasaki and Engel (2010); the physically relevant states are used by this scheme. The dimension is $\simeq$70000 for $K\leq 1$ and less than 40000 for $K\geq 2$ in the like-particle QRPA (lpQRPA) calculations. $K$ is the $z$ component of the nuclear angular momentum, and the actual computation is performed only for $K\geq 0$ on the basis of the time-reversal invariance of the ground state. The larger space is used for $K\leq 1$ for the separation of the spurious states from the real states. For the proton-neutron QRPA (pnQRPA), the dimension is always less than 40000. The solutions of $K\leq 8$ were obtained.

II.3 Calculation of $\bm{\beta\beta}$ nuclear matrix element

Our unique procedure Terasaki (2016) to obtain the NME of the $\beta\beta$ decay is as follows: first, the Gamow-Teller (GT) component of the $0\nu\beta\beta$ NME is calculated using the virtual two-particle transfer path [$(Z,N)\rightarrow(Z,N-2)\rightarrow(Z+2,N-2)$, where $Z$ is the proton number, and $N$ is the neutron number], which is possible under the closure approximation. The lpQRPA is used for obtaining the intermediate states. Our numerical solutions of the pnQRPA (see below) do not have a complex-energy solution, thus, the HFB ground state is not in the proton-neutron pair condensation. Therefore, the proton-neutron pairing interaction does not contribute to the lpQRPA calculation. The spurious states inherent to the lpQRPA are not included in the calculation of the NME.

The $0\nu\beta\beta$ NME is also obtained by the original $\beta\beta$ path [$(Z,N)\rightarrow(Z+1,N-1)\rightarrow(Z+2,N-2)$] using the pnQRPA. The average of the proton-proton and neutron-neutron pairing interactions is used for the strength of the isovector proton-neutron pairing interaction assuming the isospin invariance of the system. The strength of the isoscalar pairing interaction is determined so as to reproduce the GT component of the $0\nu\beta\beta$ NME obtained by the virtual-path calculation. The equivalence of the two decay paths is a constraint on the effective interactions for the QRPA. The assumption for this step is that the effective interactions other than the isoscalar pairing interaction are established. The strengths of the pairing interactions are summarized in Table 2. The Fermi component of the NME is controlled mainly by the isovector proton-neutron pairing interaction. The tensor component Doi et al. (1985) of the NME is neglected throughout our paper because it is known that this component is relatively small.

Next, the $2\nu\beta\beta$ NME is calculated by the pnQRPA with those strengths of the pairing interactions, and the effective $g_{A}$ is determined so as to reproduce the experimental half-life in the $2\nu\beta\beta$ decay. Finally, the $0\nu\beta\beta$ NME is calculated with this $g_{A}$ and the pnQRPA solutions already obtained. This $g_{A}$ is a choice, and the result with the bare value $g_{A}=1.27$ is also shown later. We checked the convergence of the result with respect to the dimension of the single-particle spaces (see below). Therefore, our calculation is less uncertain than calculations without this check.

III Calculation results and tests

III.1 Self-check using $\bm{2\nu\beta\beta}$ decay

Two NMEs of the $2\nu\beta\beta$ decay, $M^{(2\nu)}(I)$ and $M^{(2\nu)}(F)$, are calculated by111The notation of Ref. Doi et al. (1985) is used. It is also used for the $0\nu\beta\beta$ NME.

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

$|I\rangle$ and $|F\rangle$ are the ground states of the initial and final nuclei, respectively. The intermediate states are denoted by $|a^{K}_{I}\rangle$ (obtained from $|I\rangle$) and $|a^{K}_{F}\rangle$ (obtained from $|F\rangle$). Those used for the Fermi NME [Eqs. (LABEL:eq:M2vF(I)) and (LABEL:eq:M2vF(F))] are the $K=0$ states, and the $K$ is omitted. The operator changing a neutron to proton is denoted by $\tau^{-}$, and $\sigma_{K}$ is the spherical component of the spin-Pauli matrix operator; $\tau^{-}\sigma_{K}$ in the equations of the GT component [Eqs. (3) and (4)] is a one-body operator, and $\tau^{-}$ in Eqs. (LABEL:eq:M2vF(I)) and (LABEL:eq:M2vF(F)) is another one-body operator. $E^{I}_{aK}$ and $E^{F}_{aK}$ are the energies of the intermediate states obtained from $|I\rangle$ and $|F\rangle$, respectively. $\bar{M}$ is the mean value of the masses of the initial and final nuclei, and $m_{e}c^{2}$ is the electron mass. The parameter $g_{V}$ is the vector-current coupling, which is always equal to 1 in our calculations. The overlap $\langle a_{F}^{K}|a_{I}^{K}\rangle$ is calculated according to the method of Refs. Terasaki (2012, 2013, 2015); the many-body correlations are included by an expansion-truncation. The only difference between $M^{(2\nu)}(I)$ and $M^{(2\nu)}(F)$ is which energy ($E^{I}_{aK}$ or $E^{F}_{aK}$) is used in the energy denominator. Usually, the mean value of $E^{I}_{aK}$ and $E^{F}_{aK}$ is used. If the higher-order many-body effects beyond the QRPA is small, we will have

[TABLE]

This is the check point of the validity of the QRPA approach to the $2\nu\beta\beta$ decay.

The half-life in this decay is calculated by

[TABLE]

where $G_{2\nu}^{(0)}$ is the phase-space factor, and we refer to the values in Ref. Kotila and Iachello (2012). Table 3 shows the results of calculation of these decays of 136Xe and 130Te. For 136Xe $\rightarrow$ 136Ba, the same value of $g_{A}$ can be used for the two calculations to reproduce the experimental half-life, thus, the QRPA is good for this decay instance. On the other hand, the $g_{A}$ for fitting the $T^{(2\nu)}_{1/2}$ depends on the choice of the set of the intermediate-state energy for 130Te $\rightarrow$ 130Xe. This decay instance turns out to be an example for which the QRPA is not good. The absolute values of $M^{(2\nu)}_{F}/\mu_{0F}$ are much smaller than those of $M^{(2\nu)}_{GT}/\mu_{0}$ reflecting on the approximate isospin invariance of the isovector pairing interaction.

III.2 Spectrum of intermediate nucleus

The calculated spectra of the intermediate nucleus 136Cs are presented in Fig. 1 to further show the validity of the QRPA approach to the decay of 136Xe. The left spectrum is obtained by the pnQRPA calculation based on 136Xe, and the right one is obtained based on 136Ba. We used the transition strength for identifying the intermediate nucleus; see Appendix A. Overall, the two results are close to each other. States of angular momentum $J\leq 3$ and $\pi=-$ were not found in the shown energy region. There is a major gap of levels around 1 MeV for any $J^{\pi}$. The correspondence of the levels can also be found in many parts between the two calculations. There are two major differences. There is no $J^{\pi}=4^{+}$ level in the right figure. This is because there is no transition probability from the ground state of 136Ba. The levels around 3 MeV in the right lower figure are more tightly gathered than those in the left lower figure. The satisfactory result of the check using the $2\nu\beta\beta$ decay is endorsed by these spectra.

The latest status of the observed states nnd is illustrated in Fig. 2. A gap is seen between 1.0 and 1.9 MeV, and no low-$J$ negative-parity level is seen except for one to which 2- is assigned tentatively. The $1^{+}$ levels are observed through the $\beta$ decays, and the high-$J$ levels higher than 2.0 MeV are observed through fissions. Apparently, what states are observed depends on the experimental method. Thus, it is speculated that not all states have been observed in the discrete-energy region. The features of the calculated spectrum seem to be seen in the experimental data, although this agreement is not conclusive. The ground state is a $7^{-}$ state in our calculation, but it is actually a $5^{+}$ state in the experimental data. It is noted that the Skyrme parameter set SkM∗ is constructed so as to reproduce the properties of the even-even nuclear ground states on average in a broad region of the nuclear chart with emphasis on the doubly-magic nuclei, and fission barrier is also taken into account. Properties of odd-odd nuclei are not taken into account at all.

We also show the calculated spectrum of 130I in Fig. 3. Panels a and b show the levels obtained from 130Te and those obtained from 130Xe, respectively. The level density is clearly different between the two results. The origin of this difference is that the HFB ground state of 130Te is spherical, but that of 130Xe is deformed (see Table 1). The QRPA has a mathematical property that the broken symmetry in the HFB ground state is restored, e.g., Ring and Schuck (1980), however, this property is limited in the QRPA order. In our experience with the like-particle excitations of deformed nuclei Terasaki and Engel (2010, 2011), $J$ is clear in the states with large transition strengths, but mixture is seen in the states with small transition strengths. The deformation implies disappearance of the degeneracy of spherical states. Thus, Fig. 3 indicates that those states obtained from 130Xe are deformed. Therefore, the two most major values of $J$ are assigned to each level in panel b; for the method to assign $J$ to the deformed states, see Appendix B. The $\beta\beta$ decay of 130Te is disadvantageous to the QRPA approach due to the difference between the two spectra.

Figure 4 illustrates most of the experimental levels compiled in Ref. nnd . The experimental spectrum is closer to that obtained from 130Xe. Thus, we speculate that the experimental data indicate the deformation effect, assuming that the many-particle-many-hole correlations do not have significant effects to increase the level density in the low-energy region.

III.3 $\bm{0\nu\beta\beta}$ decay

The $0\nu\beta\beta$-decay NME $M^{(0\nu)}$ is calculated by

[TABLE]

[TABLE]

[TABLE]

where $p$ and $p^{\prime}$ denote the protons, and $n$ and $n^{\prime}$ denote the neutrons. Operators $c^{\dagger}_{i}$ and $c_{i}$ denote the creation and annihilation operators of particle $i$, respectively. The two-body transition matrix elements are defined by

[TABLE]

[TABLE]

The argument of $\bm{\sigma}$ and $\tau^{-}$ distinguishes the two particles that the operators act on. We use the following equation for the neutrino potential:

[TABLE]

[TABLE]

This neutrino potential is derived by neglecting the effective neutrino mass compared to the major momentum transfer by the propagating neutrino Doi et al. (1985). $R$ is the root-mean-square radius of the nucleus, $r_{12}$ is the distance between two particles, and $\bar{E}_{a}$ is the average energy introduced in the closure approximation. $R=1.2A^{1/3}$ fm, with the mass number $A$, and $\bar{\mu}_{a}=18.51$ are used in our calculations. In Eq. (16), functions

[TABLE]

[TABLE]

are used. Our simple neutrino potential does not include the effects of the dipole form factors Šimkovic et al. (1999), which reduces the $0\nu\beta\beta$ NME by $20-30$ %, the short-range correlations Šimkovic et al. (2009) also reducing the NME, and other subleading terms. Our $0\nu\beta\beta$ NME could be overestimated by up to, or at least, 50 %.

We calculate reduced half-life in the $0\nu\beta\beta$ decay

[TABLE]

where $G_{0\nu}$ is the phase-space factor of the $0\nu\beta\beta$ decay. $R^{(0\nu)}_{1/2}$ is the quantity necessary for determining the effective neutrino mass $\langle m_{\nu}\rangle$ as seen from

[TABLE]

with the half-life in the $0\nu\beta\beta$ decay $T^{(0\nu)}_{1/2}$. The correct approximations should give the same $R^{(0\nu)}_{1/2}$, if different effective $g_{A}$ values are used, because $\langle m_{\nu}\rangle$ and $T^{(0\nu)}_{1/2}$ are unique.222 $M^{(0\nu)\prime}=[g_{A}^{2}/(g_{A}^{\mathrm{bare}})^{2}]M^{(0\nu)}$ is also used sometimes in the literatures, e.g., Šimkovic et al. (2008), with the bare value of $g_{A}^{\mathrm{bare}}\approx 1.27$.

Our calculated NMEs and $R^{(0\nu)}_{1/2}$ are summarized in Tables 4 and 5, and $R^{(0\nu)}_{1/2}$’s of different groups are compared in Fig. 5 ($g_{A}$ = 1.25$-$1.27) together with those of 48Ca that we calculated previously. It is usual to use the bare value of $g_{A}$ around 1.26 for comparison of different calculations, e.g., Engel and Menéndez (2017). Our values, solid red circles in Fig. 5, are always low but are not the lowest in the distribution and are close to those of IBM-2. The possible overestimation of the $0\nu\beta\beta$ NME by 50 % mentioned before corresponds to the possible underestimation of $R^{(0\nu)}_{1/2}$ by a factor of 1/2.25. The method of Ref. Mustonen and Engel (2013) to calculate the QRPA wave functions is close to our method, however, the corresponding $R^{(0\nu)}_{1/2}$’s are rather different in 130Te and 136Xe. It is speculated that the difference comes from the different Skyrme interactions used in those calculations.

The IBM-2 Barea et al. (2015) was used with two values of $g_{A}$: the unquenched value (the authors use 1.269, see Fig. 5), and their effective values of “maximal quenching” (the result is not shown in the figure). These effective values are determined by their simple $A$-dependent formula obtained from the $2\nu\beta\beta$ decays; they are 0.528 (130Te) and 0.524 (136Xe). These values are close to our effective $g_{A}$ in Table 4.

The authors of Ref. Pirinen and Suhonen (2015) picked up triplets of nuclei connected by the $\beta^{-}$ or $\beta^{+}$ decay or the electron capture from the mass region of A=100$$-$$134 and found that the geometric means of the two calculated NMEs for the “left-middle” and “right-middle” combinations of nuclei in those triplets are rather independent of the strength of the pairing interaction (scaled by $g_{pp}$ in their notation), if $g_{pp}$ is not large. On the basis of this interesting discovery, they determined the effective $g_{A}$ so as to reproduce the experimental mean values, and a simple fitting formula of $g_{A}$ was obtained as a function of $A$. If this formula can be applied to $A=136$, their effective value is 0.833. The $0\nu\beta\beta$ NME is not shown in that paper.

An early example of an effective $g_{A}$ value smaller than one is shown in Ref. Faessler et al. (2008). The conclusion of this reference is that $g_{A}<1$ is necessary for reproducing the experimental lifetimes to both the single-charge-change and the $\beta\beta$ decays by the pnQRPA. An effective $g_{A}$ value of 0.39 for 128Te was also included in the first version of that reference Faessler et al. (2007).

Recently, the authors of Ref. Šimkovic et al. (2018a) suggested a condition to determine $g_{pp}$ using the $2\nu\beta\beta$ NME with the closure approximation and the spin and isospin symmetries. They obtained a quenching factor to $g_{A}$ of 0.712 using the averaged ratio of the $2\nu\beta\beta$ NME calculated without the closure approximation and the experimental one with respect to $A$. Their effective $g_{A}$ value is 1.00, which is one of the values that they have used; for their previous calculations, see Ref. Šimkovic et al. (2018b) and references therein. For comprehensive recent review on the effective $g_{A}$, see Refs. Suhonen (2017b); Engel and Menéndez (2017).

Figure 6 shows the convergence of the NME of 136Xe $\rightarrow$ 136Ba with respect to the single-particle space. We use very large two-quasiparticle spaces for the QRPA calculations, and the test of the convergence was made by varying the two-single-particle spaces used for calculating $V^{GT(0\nu)}_{pp^{\prime},nn^{\prime}}(\bar{E}_{a})$ and $V^{F(0\nu)}_{pp^{\prime},nn^{\prime}}(\bar{E}_{a})$. The well-converged results were used for the calculation of the NMEs. The behavior of the convergence of 136Xe is less smooth than that of 130Te (Fig. 7) because the ground states of 136Xe and 136Ba are spherical, and the neutrons of 136Xe are unpaired; see Table 1.

III.4 Charge-change reaction

There are experimental data of 136Xe(3He,t)136Cs Puppe et al. (2011); Frekers et al. (2017), that can be used for test of the charge-change transition density $\langle a^{K}_{I}|c_{p^{\prime}}^{\dagger}c_{n^{\prime}}|I\rangle$ [see Eqs. (12) and (13)] independently of $g_{A}$. The experimental and calculated GT- strengths are shown in Fig. 8 for a low-energy region in which the GT--strength data are obtained. The calculated one is the summation with respect to $K$ of the squared transition matrix element $|M_{\mathrm{GT}}^{K}|^{2}$ of the GT- operator $\sigma\tau^{-}$, and the data are obtained from the cross section. The two results are similar in terms of the energy dependence, however, the calculated values are larger than the data by one order of magnitude. The GT sum-rule value of our calculation is

[TABLE]

where the first term is the GT- component, and the second term is the GT+ component (136Xe $\rightarrow$ 136I). The sum rule is satisfied well, i.e., our value is very close to $3(N-Z)$.

There was a similar problem in 48Ca $\rightarrow$ 48Sc and 48Ti $\rightarrow$ 48Sc Yako et al. (2009). This problem was solved Terasaki (2018) by using the transition operator consisting of the GT and isovector spin-monopole operators phenomenologically. This idea can be tested, if the GT-strength data including the giant-resonance region are obtained. Actually the yield data are obtained up to 30 MeV. For the shell-model approach to these data, see Ref. Neacsu and Horoi (2015).

We show in Fig. 9 the GT- strength obtained from the reaction of 130Te(3He,t)130I Puppe et al. (2012) and that obtained by our calculation. Contrary to the GT- strength of 136Xe$\rightarrow$136Cs, we have much fewer $1^{+}$ states than the data have. The deformation effect discussed in Sec. III.2 is also reflected here. (The data show more $1^{+}$ states than the compilation of Ref. nnd .) The total strength up to 5 MeV is 1.5 Puppe et al. (2012) and 11.4 in our calculation, thus, the same argument as that for 136Xe$\rightarrow$136Cs is applied here; the experimental strength distribution in a very large energy region is necessary for discussing the possible necessity of the isovector spin-monopole operator. Our GT sum-rule value is 77.951, and the exact value is 78.

III.5 $\bm{\beta}$ decay

We calculated the probability of the $\beta$ decay of 138Xe for testing our $g_{A}$. The calculated spectrum and logft are compared to those of the experimental data nnd for $1^{+}$ states in Fig. 10. The level distribution of the data is slightly denser and lower than the corresponding calculated levels. We obtained a level around $-$2.9 MeV, however, there is no corresponding experimental level. Actually there is a level of $0^{-}$ or $1^{-}$ around that energy in the data. We plot

[TABLE]

in Figs. 12 (our calculation) and 12 (experiment). The summation of the calculated $B_{\mathrm{GT}}$ between $-1.2$ MeV and $-0.3$ MeV is 0.319, and the corresponding experimental value is 0.281; the former is 14 % larger than the latter. Thus, our $g_{A}$ fitting the $T^{(2\nu)}_{1/2}$ of 136Xe is consistent with the $\beta$ decay of a nearby nucleus.

The GT- spectrum of 132Te $\rightarrow$ 132I is shown in Fig. 13 with the logft values. Only one $1^{+}$ final state is observed experimentally nnd . The logft by our calculation is smaller than the data by $\simeq$ 0.5$-$0.8, and the corresponding $B_{\mathrm{GT}}$ ratio is 3.2$-$6.3. The quality of this calculation is not as good as that for 138Xe.

III.6 Higher-order effect in 2$\bm{\nu\beta\beta}$ decay

Very recently, a higher-order component of the $2\nu\beta\beta$ NME Šimkovic et al. (2018b)

[TABLE]

was suggested, and the experimental value of the ratio

[TABLE]

was extracted for 136Xe by fitting the $2\nu\beta\beta$ spectrum Gando et al. (2019). In Eq. (22), $\mu_{a}$ refers to Eq. (7) or (8); there is no difference for 136Xe. Our value is $\xi^{(2\nu)}_{31}=0.040$, which is within the experimentally reported range of $\xi^{(2\nu)}_{31}=-0.26^{+0.31}_{-0.25}$ determined by $1\sigma$ for exclusion. This comparison also enhances the validity of our calculation. Because this discussion is very new, $M^{(2\nu)}_{GT\textrm{-}3}/\mu_{0}$ is not included in the $2\nu\beta\beta$ NME in Sec. III.1.

IV Summary

In this paper, we have calculated the NMEs of the $\beta\beta$ decays of 136Xe and 130Te, from which $R^{(0\nu)}_{1/2}$’s were derived. Several tests have been conducted to investigate the reliability of our calculation of $R^{(0\nu)}_{1/2}$:

the self-check of the dual intermediate states using the $2\nu\beta\beta$ decay. 2. 2.

The convergence of the NME with respect to the dimension of the single-particle space. 3. 3.

The comparison of the two spectra of the intermediate nucleus obtained by the pnQRPA and comparison of them with the experimental data. 4. 4.

The GT sum rule. 5. 5.

The comparison of the GT- strength with the data of the charge-change reaction. 6. 6.

The comparison of the $\beta$ decay spectrum and $B_{\mathrm{GT}}$ with the experimental data. 7. 7.

Test using a new quantity expressing a higher-order effect in the $2\nu\beta\beta$ NME.

No problem was found for 136Xe concerning our methodology. On the other hand, it turned out that our approach was not successful for 130Te.

A key point of our procedure is the strength of the isoscalar pairing interaction. Our value is not very large, so that the pnQRPA solutions are not close to the unstable region. Actually, this value makes the effective $g_{A}$ smaller. Thus, we investigated if this $g_{A}$ was consistent with the $\beta$ decay, and the consistency was obtained. The effective $g_{A}$ depends on the approximation. To our knowledge, the spectrum of the intermediate nucleus by the pnQRPA has been shown for the first time. The comparison of the calculated GT- strength with the experimental data (136Xe) is encouraging in terms of the energy dependence, however, a problem remains about how the experimental absolute value can be reproduced; a possibility is to include the isovector spin-monopole operator in the transition operator. The GT--strength data in a broader energy region are necessary for the test of the calculation. That transition operator is not obvious a priori in the charge-change reaction.

Acknowledgements.

The numerical calculations of this paper were performed by the K computer at RIKEN Center for Computational Science through the program of High Performance Computing Infrastructure in 2017B (hp170288) and 2018B (hp180232). Computer Oakforest-PACS operated by Joint Center for Advanced High Performance Computing was also used through Multidisciplinary Cooperative Research Program of Center for Computational Sciences, University of Tsukuba in 2018 (xg18i006). Oakforest-PACS has also been used through High Performance Computing Infrastructure (hp190001, 2019). This study is supported by European Regional Development Fund, Project “Engineering applications of microworld physics” (No. CZ.02.1.01/0.0/0.0/16_019/0000766).

Appendix A Identification of nucleus in QRPA solutions

The QRPA solutions include multiple nuclei, if the HFB ground state is paired. In the low-energy region, the identification of the solutions corresponding to a specific nucleus is possible only approximately. We introduce the auxiliary transition operator $O^{-}_{\mathrm{id}}(J^{\pi})$ for this identification:

[TABLE]

[TABLE]

and those for other $J^{\pi}$ in the analogous way. We calculate the transition strengths of these operators for the QRPA solutions based on 136Xe and regard those having relatively large strengths as the states of 136Cs. The $J^{\pi}$ of the state is identified with that of the transition operator of the large strength. The expectation value of $r^{8}$ is calculated with respect to the HFB ground state. This $r^{8}$ is used for covering a large energy region, as seen by expressing it using the creation and annihilation operators of the harmonic oscillator. The exponent of 8 may be reasonable because $K$ values up to 8 are treated in the NME calculation of the $\beta\beta$ decay.

For identifying 136Cs obtained from 136Ba, the same operators but for $\tau^{+}=\tau^{-\dagger}$ instead of $\tau^{-}$ are used. In neutron-rich nuclei, the $p\rightarrow n$ transition is more rare compared to the $n\rightarrow p$ transition. Therefore, all QRPA solutions that have the finite transition strength of $p\rightarrow n$ were picked up.

Appendix B Assignment of $\bm{J}$ to deformed QRPA solutions

We modify slightly the auxiliary transition operator in Appendix A for assigning major $J$ to the deformed states obtained from 130Xe (Sec. III.2). If $J^{\pi}$ is of the unnatural parity,

[TABLE]

is used. ${\cal S}(J^{\pi})$ is the factor normalizing the sum of the transition strengths for different $J^{\pi}$’s. The transition strengths of this operator are calculated for all QRPA solutions ($K\leq 8$) with the specified parity in a low-energy region, and the most and next-most major $J$’s are noted in panel b of Fig. 3. This method was confirmed to give the correct $J$ for the spherical states of 130Te; the correct $J$ is seen by the degeneracy. This check is not trivial because our calculation uses the single-particle wave functions in the cylindrical box with the vanishing boundary condition. For the natural-parity states,

[TABLE]

is used.