Coupled-cluster calculations of neutrinoless double-beta decay in 48 Ca
S. J. Novario, P. Gysbers, J. Engel, G. Hagen, G. R. Jansen, T. D. Morris, P. Navrátil, T. Papenbrock, S. Quaglioni
CCoupled-cluster calculations of neutrinoless double-beta decay in Ca S. J. Novario,
1, 2
P. Gysbers,
3, 4
J. Engel, G. Hagen,
2, 1, 3
G. R. Jansen,
6, 2
T. D. Morris, P. Navr´atil, T. Papenbrock,
1, 2 and S. Quaglioni Department of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA TRIUMF, 4004 Wesbrook Mall, Vancouver BC, V6T 2A3, Canada Department of Physics and Astronomy, University of British Columbia, Vancouver BC, V6T 1Z1, Canada Department of Physics, University of North Carolina, Chapel Hill, NC 27514, USA National Center for Computational Sciences, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Lawrence Livermore National Laboratory, P.O. Box 808, L-414, Livermore, California 94551, USA
We use coupled-cluster theory and nuclear interactions from chiral effective field theory to computethe nuclear matrix element for the neutrinoless double-beta decay of Ca. Benchmarks with theno-core shell model in several light nuclei and in the traditional shell model in the pf shell informus about the accuracy of our approach. For Ca we find a relatively small matrix element. Wealso compute the nuclear matrix element for the two-neutrino double-beta decay of Ca and findagreement with data when using a quenching factor deduced from two-body currents in the recentab-initio calculation of the Ikeda sum-rule in Ca [Gysbers et al. , Nature Physics , 428-431(2019)]. Introduction and main result.—
Neutrinoless double-beta (0 νββ ) decay is a hypothesized electroweak processin which a nucleus undergoes two simultaneous beta de-cays but emits no neutrinos [1]. The observation of thislepton-number violating process would identify the neu-trino as a Majorana particle (i.e. as its own antiparti-cle) [2] and provide insights into both the origin of neu-trino mass [3, 4] and the matter-antimatter asymmetryin the universe [5]. Experimentalists are working intentlyto observe the decay all over the world; current lower lim-its on the lifetime are about 10 y [6–10], and sensitivitywill be improved by two orders of magnitude in the com-ing years.Essential for planning and interpreting these experi-ments are nuclear matrix elements (NMEs) that relatethe decay lifetime to the Majorana neutrino mass scaleand other measures of lepton-number violation. Unfor-tunately, these matrix elements are not well known andcannot be measured. Computations based on differentmodels and techniques lead to numbers that differ byfactors of three to five (see Ref. [11] for a recent re-view). Compounding these theoretical challenges is therecent discovery that, within chiral effective field theory(EFT) [12–15], the standard long-range 0 νββ decay oper-ator must be supplemented by an equally important zero-range (contact) operator of unknown strength [16]. Ef-forts to compute the strengths of this contact term fromquantum chromodynamics (QCD) [17, 18] and attemptsto better understand its impact are underway [19].The task theorists face at present is to provide moreaccurate computations of 0 νββ NMEs, including thoseassociated with contact operators, and quantify their un-certainties. In this Letter, we employ the coupled-clustermethod to perform first-principle computations of thematrix element that links the 0 νββ lifetime of Ca withthe Majorana neutrino mass scale. Among the dozen CC S D CC S D T - I M S R G + G C M R S M Q R P A I B M E D F S M ( p f ) S M ( M B P T ) S M ( s dp f ) . . . . . . . M ν SphericalDeformed
FIG. 1. (Color online) Comparison of the NME for the 0 νββ decay of Ca, calculated within various approaches (see textfor details). The coupled-cluster results use both the CCSDand CCSDT-1 approximations with both the spherical anddeformed reference states. For IMSRG+GCM, the doublebars show the effects of uncertainty in model-space size; other-wise they show those of uncertainty in short-range correlationfunctions. or so candidate nuclei for 0 νββ decay experiments [20], Ca stands out for its fairly simple structure, making itamenable for an accurate description based on chiral EFTand state-of-the-art many-body methods [21]. By vary-ing the details of our calculations, we will estimate theuncertainty of our prediction. To gauge the quality of ourapproach we also compute the two-neutrino double-betadecay of Ca and compare with data. Our results will di-rectly inform 0 νββ decay experiments that use Ca [22]and serve as an important stepping stone towards theaccurate prediction of NMEs in Ge,
Te, and
Xe,which are candidate isotopes of the next-generation 0 νββ decay experiments. a r X i v : . [ nu c l - t h ] A ug Figure 1 shows several recent results for the NME gov-erning the 0 νββ decay Ca → Ti and compares themwith those of this work. The coupled cluster results ob-tained here, with both the CCSD and CCSDT-1 approxi-mations (explained below), display uncertainties from de-tails of the computational approach. They are comparedto the very recent ab initio results from the in-mediumsimilarity group renormalization method with the gen-erator coordinator method (IMSRG+GCM) [23], a real-istic shell-model (RSM) [24], the quasi-particle randomphase approximation (QRPA) [25], the interacting bo-son model (IBM) [26], various energy-density functionals(EDF) [27, 28], and several more phenomenological shellmodel (SM) calculations. The latter either limit them-selves to the pf -shell [29, 30], include perturbative cor-rections from outside of the pf -shell [31], or are set in the sdpf shell-model space [32]. We see that the ab initio re-sults of this work and of Ref. [23] are consistent with eachother and with the most recent work [33]. Our result, inthe CCSDT-1 approximation, is 0 . ≤ M ν ≤ . Method.—
We employ the intrinsic Hamiltonian H = (cid:88) i 16, where N i = 2 n i + l i are single-particle energies. The oscillatorbasis has a frequency (cid:126) Ω = 16 MeV and we find thatworking within a model space with N i = 10 is sufficientto produce converged results.Following Refs. [35, 36], we transform the Hamiltonianfrom the spherical oscillator basis to a natural-orbitalbasis by diagonalizing the one-body density matrix. Wedenote the resulting reference state, i.e. the product stateconstructed from the A single-particle states with largestoccupation numbers, by | Φ (cid:105) and the Hamiltonian that isnormal-ordered with respect to this non-trivial vacuumby H N . We retain N N N forces at the normal-orderedtwo-body level [37, 38].Coupled-cluster theory [39–45] is based on thesimilarity-transformed Hamiltonian, H N = e − ˆ T H N e ˆ T .The cluster operator ˆ T is a sum of particle-hole (ph)excitations from the reference | Φ (cid:105) and commonly trun-cated at the two-particle two-hole (2 p –2 h ) or 3 p –3 h level.The amplitudes in ˆ T are chosen so that the referencestate | Φ (cid:105) becomes the right ground state of H N . Be-cause H N is non-Hermitian, the left ground state is (cid:104) Φ | (1 + ˆΛ), where ˆΛ is a de-excitation operator with re-spect to the reference [44, 45]. In this paper, we work atthe leading-order approximation to coupled-cluster with singles-doubles-and-triples excitations (CCSDT), knownas CCSDT-1 [46, 47]. To make the computation feasible,we truncate the 3 p –3 h amplitudes by imposing a cut onthe product of occupation probabilities n a for three par-ticles above the Fermi surface, n a n b n c ≥ E , and for threeholes below the Fermi surface, (1 − n i )(1 − n j )(1 − n k ) ≥E . This truncation favors orbitals near the Fermi sur-face. The limits are large enough so that all CCSDT-1 results presented below are stable against changes inthem.We are interested in computing | M ν | = (cid:104) Ψ I | ˆ O † ν | Ψ F (cid:105)(cid:104) Ψ F | ˆ O ν | Ψ I (cid:105) , where ˆ O ν is the 0 νββ operator and Ψ I and Ψ F denote the ground statesof the initial and final nuclei, respectively. Withincoupled-cluster theory, we can structure the calculationin two ways. In a first approach, we can use the rightand left ground states of Ca ( | Φ (cid:105) and (cid:104) Φ | (1 + ˆΛ),respectively) to compute | M ν | = (cid:104) Φ | (1 + ˆΛ) O † ν ˆ R | Φ (cid:105)(cid:104) Φ | ˆ LO ν | Φ (cid:105) . (2)In this case, we use equation-of-motion coupled-cluster(EOM-CC) techniques [44, 48–53] to represent the rightand left Ti ground states (denoted by ˆ R | Φ (cid:105) and (cid:104) Φ | ˆ L ,respectively) by generalized excited states of Ca withtwo more protons and two less neutrons [54, 55]. Here,we also work in the CCSDT-1 approximation. In Eq. (2) O ν ≡ e − ˆ T ˆ O ν e ˆ T is the similarity-transformed 0 νββ op-erator.In an alternative approach, we can decouple the groundstate of the final nucleus, i.e. take | Φ (cid:105) as a reference rightground state for Ti (with (cid:104) Φ | (1 + ˆΛ) its left groundstate), and target the initial nucleus Ca with EOM-CC. This procedure leads to the expression | M ν | = (cid:104) Φ | ˆ LO † ν | Φ (cid:105)(cid:104) Φ | (1 + ˆΛ) O ν ˆ R | Φ (cid:105) , (3)where the Ca right and left ground states ( ˆ R | Φ (cid:105) and (cid:104) Φ | ˆ L , respectively) are represented by generalized ex-cited states of Ti. Because the two approaches areidentical only when the cluster operators are not trun-cated, the difference between them is a measure of thetruncation effects. As the ground state of Ca is spher-ical, the first procedure allows us to exploit rotationalsymmetry. By contrast, starting from Ti introduces adeformed (though axially symmetric) reference state.In chiral EFT, the 0 νββ operator is organized intoa systematically improvable expansion similarly to thenuclear forces [56]. The lowest-order contributions tothe 0 νββ operator are a long-range Majorana neutrinopotential that can be divided into three components,Gamow-Teller (GT), Fermi (F), and tensor (T), thatcontain different combinations of spin operators, withˆ O ν = ˆ O GT0 ν + ˆ O F0 ν + ˆ O T0 ν . The corresponding two-body matrix elements, as is conventional, are taken fromRef. [57], which adds form factors to the leading andnext-to-leading operators. We use the closure approxi-mation (which is sufficiently accurate [29]), with closureenergies E cl = 5 MeV for all benchmarks in light nucleiand 7 . 72 MeV for the decay Ca → Ti.The NME for the 2 νββ is similar to the 0 νββ caseexcept the two-body operator is replaced by a double ap-plication of the one-body Gamow-Teller operator, στ − ,with an explicit summation over the intermediate 1 + states between them, | M ν | = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:88) µ (cid:104) Ti | στ − | + µ (cid:105)(cid:104) + µ | στ − | Ca (cid:105) ∆ E µ + Q ββ / (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) . (4)The direct computation of the matrix element (4) wouldrequire several tens of states in the intermediate nucleus Sc and several hundred Lanczos iterations, making itunfeasible in our large model space. We note that theGreen’s function at the center of this matrix elementcan be computed efficiently using the Lanczos (contin-ued fraction) method starting from a 1 + pivot state [58–62]. We generate Lanczos coefficients ( a i , b i and a ∗ i , b ∗ i )from a non-symmetric Lanczos algorithm using the 1 + subspace of H N and rewrite Eq. (4) as a continued frac-tion [58]. This computation typically requires about 10-20 Lanczos iterations. With the similarity-transformedoperator, O = στ − , and the pivot states (cid:104) ν F | = (cid:104) Φ | LO , | ν I (cid:105) = O | Φ (cid:105) , (cid:104) ν I | = (cid:104) Φ | (1 + ˆΛ) O † , and | ν F (cid:105) = O † R | Φ (cid:105) ,the NME becomes | M ν | = (cid:104) ν F | ν I (cid:105) a + Q ββ − b a + ··· (cid:104) ν I | ν F (cid:105) a ∗ + Q ββ − ( b ∗ ) a ∗ + ··· . (5) Benchmarks.— To gauge the quality of our coupled-cluster computations we benchmark with the more exactno-core shell model (NCSM) [63–65] by computing 0 νββ matrix elements in light nuclei. Although the 0 νββ de-cay of these isotopes are energetically forbidden or wouldbe swamped by successive single- β decays in an experi-ment, the benchmarks still have theoretical value. Fig-ure 2 shows the 0 νββ matrix elements of the GT, F,and T operators for the transitions He → Be, He → Be, He → Be, C → O, and O → Ne. The coupled-cluster results are shown in pairs, with both the initialand final state as the reference. For each pair, the first(second) point shows the CCSD (CCSDT-1) approxima-tion; these two points are connected by dotted lines. Thevertical error bars indicate the change of the matrix ele-ment as the model space is increased from N max = 8 to N max = 10. The NCSM results are shown in the thirdcolumn, and their error bars indicate uncertainties fromextrapolation to infinite model spaces. The shaded bandsare simply to facilitate comparison.The NMEs in the mirror-symmetric cases He → Beand C → O depend very little (within about 1%) onthe choice of the initial or final nucleus as the referencestate, a result that is consistent with the weak charge-symmetry breaking of the chiral interaction. For the A = | Φ i = | O i | Φ i = | Ne i NCSM M ν O → Ne N LO (EM) λ =1 . | Φ i = | C i | Φ i = | O i NCSM C → O . / . | Φ i = | He i | Φ i = | Be i NCSM -0.20.00.20.40.6 He → Be . / . | Φ i = | He i | Φ i = | Be i NCSM He → Be . / . | Φ i = | He i | Φ i = | Be i NCSM He → Be . / . CCSD CCSDT-1GTFT FIG. 2. (Color online) Comparison of the 0 νββ NME in sev-eral light nuclei computed with the coupled cluster methodand the no-core shell model. The first two columns corre-spond to different choices for the coupled-cluster referencestate, and results from the CCSD and CCSDT-1 approxima-tions are shown in each. The error bars indicate the uncer-tainties coming from variations with model-space size. 14 transition between doubly closed-shell nuclei, coupled-cluster theory and NCSM results agree within about 3%.The small contributions of triples correlations ( < He → Be, even though these nuclei areonly semi-magic. The case of He → Be is slightly morechallenging, with a doubly closed-shell initial nucleus anda partially closed-shell final nucleus.The cases of He → Be and O → Ne are more chal- Ca → Ti Ti → Cr Cr → Fe M ν Deformed CC Ca → Ti Spherical CC Deformed CCCCSDT-1(GXPF1a)CCSDT-1(KB3G)FCI(GXPF1a)FCI(KB3G) FIG. 3. (Color online) Comparison in several pf -shell of the0 νββ NMEs between CCSDT-1 and exact shell-model calcu-lations, with the GXPF1A and KB3G interactions. In Ca, Ti, and Cr we use a deformed reference state in the ini-tial nucleus, while for the decay Ca → Ti we use referencestates in both nuclei. lenging still, because the final nuclei are truly open-shellsystems. Adding triples correlations to the spherical re-sults induces a ∼ 50% change in the first case and wors-ens the agreement with NCSM in the second, suggestingthe need for more particle-hole excitations. Once again,however, using the deformed final state as the referenceleads to results that are both consistent with the NCSMand converged at the CCSDT-1 level. Thus, the coupled-cluster results are more accurate when the open-shell (ordeformed) nucleus is taken as the reference, and theyagree within smaller model-space uncertainties with theNCSM benchmarks.To benchmark calculations in heavier nuclei, we com-pare coupled-cluster results with exact results from shell-model calculations in the pf shell [66] with the phe-nomenological interactions GXPF1A [67] and KB3G [68].Figure 3 shows CCSDT-1 results and compares them toexact results for the transitions Ca → Ti, Ti → Cr, Cr → Fe, and Ca → Ti. In the first of these, thevalence shell contains only two nucleons and the problemis thus exactly solvable with CCSD. The next two casesare particularly challenging because initial and final nu-clei are open-shell systems. Here, coupled-cluster NMEsare significantly smaller than their exact counterparts.For the most relevant case, Ca → Ti, coupled-clusterresults are ∼ 15% lower (higher) than the benchmarkswhen Ti ( Ca) serves as the reference, i.e. when weuse a deformed (spherical) reference state.The benchmark calculations thus suggest that the twoapproaches (with a spherical Ca or a deformed Ti asthe reference state) allow us to bracket the NME. Theresult from the first approach exceeds the exact NMEbecause the imposition of spherical symmetry increasesthe overlap of the initial and final wave functions. Thesecond result underestimates the exact NME, probablybecause the deformations of the intial and final states are Q ββ / M ν Quenchingfactor 0 . Ca → Ti . / . E ( Ca) = 10 − E ( Ti) = 10 − E ( Sc) − − − − N max (6)N max (8)CCExp. FIG. 4. (Color online) The NME for the 2 νββ decay Ca → Ti computed with the 1.8/2.0 (EM) interaction asa function of the Q -value, Q ββ , and the 3 p –3 h truncationused to calculate Sc, E , at N max = 10. The results for N max = 6 , . deduced from two-bodycurrents [77]. quite different. Presumably, symmetry projection wouldincrease this result to some extent.Although the coupling strength of the leading-ordercontact potential in the 0 νββ operator is unknown [16,17, 19], we can estimate its effect by applying the coupledcluster methods discussed above with the addition of acontact term, V c ( r ) = gδ ( r ) τ (1) − τ (2) − , to the operator,ˆ O ν . Using a coupling strength of g = ± . ≤ M ν ≤ . Two-neutrino double-beta decay of Ca.— The 2 νββ decay of Ca was accurately predicted by Caurier et al. [69] before its observation [70–72]. Subsequent authorsstudied this decay further [73–75], and evaluations can befound in Refs. [20, 76]. We compute the matrix elementfor the 2 νββ decay of Ca with the 1.8/2.0 (EM) inter-action and the Lanczos continued fraction method. Weemploy a spherical Ca natural-orbital basis and con-verge our results with respect to N max and the numberof 3 p –3 h configurations included in the wave functions of Ca, Ti, and Sc. The results are also converged withrespect to the number of Lanczos iterations used in thecontinued fraction (5).Figure 4 shows the NME for the 2 νββ decay of Ca,computed in the CCSDT-1 approximation, as a func-tion the Q -value Q ββ , with different curves represent-ing both the N max convergence and E convergence ofthe intermediate nucleus Sc. The converged result, M ν = 0 . ± . Q -value, Q ββ / . 13 MeV, which is closeto the theoretical result Q ββ / . 10 MeV, i.e. the dif-ference between the ground-state energies of Ca and Ti computed from the corresponding reference states.The uncertainty in our result represents the error fromthe different convergence criteria.Multiplying our matrix element with the a quenchingfactor q = 0 . deduced from two-body currents ina recent coupled-cluster computation of the Ikeda sum-rule in Ca [77], we obtain q M ν = 0 . ± . M ν =0 . ± . 003 [76, 78]. The Ikeda sum-rule includes allfinal 1 + states in Sc and is similar to Eq. (4). In afuture work we will investigate the role of momentumdependent two-body currents on this decay. We verifiedour methods by performing two 2 νββ benchmarks, of Ca in the pf -shell and of C in a full no-core modelspace, shown in the supplementary material in Figures 10and 11, respectively. The former is compared with exactdiagonalization, and the latter with the NCSM. Conclusions.— Using interactions from chiral EFT andthe coupled-cluster method, we computed the nuclearmatrix elements for 0 νββ -decay of Ca → Ti and founda relatively small value. The uncertainties stem fromthe treatment of nuclear deformation and are supportedby extensive benchmarks. We also calculated the 2 νββ -decay of Ca → Ti and reproduced the experimentalvalue after including the ab-initio quenching factor fromtwo-body currents of the Ikeda sum-rule in Ca.We thank A. Belley, V. Cirigliano, J. de Vries, H. Herg-ert, J. D. Holt, M. Horoi, J. Men´endez, C. G. Payne, S.R. Stroberg, A. Walker-Loud, and J. M. Yao, for use-ful discussions. This work was supported by the Of-fice of Nuclear Physics, U.S. Department of Energy, un-der Grants DE-FG02-96ER40963, DE-FG02-97ER41019DE-SC0008499 (NUCLEI SciDAC collaboration), theField Work Proposal ERKBP57 at Oak Ridge NationalLaboratory (ORNL) and SCW1579 at Lawrence Liv-ermore National Laboratory (LLNL), the National Re-search Council of Canada, and NSERC, under GrantsSAPIN-2016-00033 and PGSD3-535536-2019. TRIUMFreceives federal funding via a contribution agreementwith the National Research Council of Canada. Thiswork was prepared in part by LLNL under ContractNo. DE-AC52-07NA27344. Computer time was providedby the Innovative and Novel Computational Impact onTheory and Experiment (INCITE) program. This re-search used resources of the Oak Ridge Leadership Com-puting Facility located at ORNL, which is supported bythe Office of Science of the Department of Energy underContract No. DE-AC05-00OR22725. [1] W. H. Furry, “On transition probabilities in double beta-disintegration,” Phys. Rev. , 1184–1193 (1939).[2] J. Schechter and J. W. F. Valle, “Neutrinoless double- β decay in su(2) × u(1) theories,” Phys. Rev. D , 2951–2954 (1982).[3] Peter Minkowski, “ µ → e γ at a rate of one out of 109muon decays?” Physics Letters B , 421 – 428 (1977). [4] Rabindra N. Mohapatra and Goran Senjanovi´c, “Neu-trino mass and spontaneous parity nonconservation,”Phys. Rev. Lett. , 912–915 (1980).[5] Sacha Davidson, Enrico Nardi, and Yosef Nir, “Lepto-genesis,” Physics Reports , 105 – 177 (2008).[6] A. Gando, Y. Gando, T. Hachiya, A. Hayashi,S. Hayashida, H. Ikeda, K. Inoue, K. Ishidoshiro,Y. Karino, M. Koga, S. Matsuda, T. Mitsui, K. Naka-mura, S. Obara, T. Oura, H. Ozaki, I. Shimizu, Y. Shi-rahata, J. Shirai, A. Suzuki, T. Takai, K. Tamae,Y. Teraoka, K. Ueshima, H. Watanabe, A. Kozlov,Y. Takemoto, S. Yoshida, K. Fushimi, T. I. Banks,B. E. Berger, B. K. Fujikawa, T. O’Donnell, L. A.Winslow, Y. Efremenko, H. J. Karwowski, D. M. Markoff,W. Tornow, J. A. Detwiler, S. Enomoto, and M. P. De-cowski (KamLAND-Zen Collaboration), “Search for ma-jorana neutrinos near the inverted mass hierarchy regionwith kamland-zen,” Phys. Rev. Lett. , 082503 (2016).[7] C. E. Aalseth, N. Abgrall, E. Aguayo, S. I. Alvis, M. Am-man, I. J. Arnquist, F. T. Avignone, H. O. Back, A. S.Barabash, P. S. Barbeau, C. J. Barton, P. J. Barton, F. E.Bertrand, T. Bode, B. Bos, M. Boswell, A. W. Bradley,R. L. Brodzinski, V. Brudanin, M. Busch, M. Buuck,A. S. Caldwell, T. S. Caldwell, Y-D. Chan, C. D.Christofferson, P.-H. Chu, J. I. Collar, D. C. Combs,R. J. Cooper, C. Cuesta, J. A. Detwiler, P. J. Doe, J. A.Dunmore, Yu. Efremenko, H. Ejiri, S. R. Elliott, J. E.Fast, P. Finnerty, F. M. Fraenkle, Z. Fu, B. K. Fujikawa,E. Fuller, A. Galindo-Uribarri, V. M. Gehman, T. Gilliss,G. K. Giovanetti, J. Goett, M. P. Green, J. Gruszko,I. S. Guinn, V. E. Guiseppe, A. L. Hallin, C. R. Haufe,L. Hehn, R. Henning, E. W. Hoppe, T. W. Hossbach,M. A. Howe, B. R. Jasinski, R. A. Johnson, K. J. Keeter,J. D. Kephart, M. F. Kidd, A. Knecht, S. I. Konovalov,R. T. Kouzes, B. D. LaFerriere, J. Leon, K. T. Lesko,L. E. Leviner, J. C. Loach, A. M. Lopez, P. N. Luke,J. MacMullin, S. MacMullin, M. G. Marino, R. D. Mar-tin, R. Massarczyk, A. B. McDonald, D.-M. Mei, S. J.Meijer, J. H. Merriman, S. Mertens, H. S. Miley, M. L.Miller, J. Myslik, J. L. Orrell, C. O’Shaughnessy, G. Oth-man, N. R. Overman, G. Perumpilly, W. Pettus, D. G.Phillips, A. W. P. Poon, K. Pushkin, D. C. Radford,J. Rager, J. H. Reeves, A. L. Reine, K. Rielage, R. G. H.Robertson, M. C. Ronquest, N. W. Ruof, A. G. Schu-bert, B. Shanks, M. Shirchenko, K. J. Snavely, N. Snyder,D. Steele, A. M. Suriano, D. Tedeschi, W. Tornow, J. E.Trimble, R. L. Varner, S. Vasilyev, K. Vetter, K. Vor-ren, B. R. White, J. F. Wilkerson, C. Wiseman, W. Xu,E. Yakushev, H. Yaver, A. R. Young, C. H. Yu, V. Yuma-tov, I. Zhitnikov, B. X. Zhu, and S. Zimmermann (Ma-jorana Collaboration), “Search for neutrinoless double- β decay in Ge with the majorana demonstrator,” Phys.Rev. Lett. , 132502 (2018).[8] M. Agostini, A. M. Bakalyarov, M. Balata, I. Bara-banov, L. Baudis, C. Bauer, E. Bellotti, S. Bel-ogurov, A. Bettini, L. Bezrukov, J. Biernat, T. Bode,D. Borowicz, V. Brudanin, R. Brugnera, A. Cald-well, C. Cattadori, A. Chernogorov, T. Comellato,V. D’Andrea, E. V. Demidova, N. Di Marco, A. Domula,E. Doroshkevich, V. Egorov, R. Falkenstein, A. Gan-gapshev, A. Garfagnini, P. Grabmayr, V. Gurentsov,K. Gusev, J. Hakenm¨uller, A. Hegai, M. Heisel, S. Hem-mer, R. Hiller, W. Hofmann, M. Hult, L. V. In-zhechik, J. Janicsk´o Cs´athy, J. Jochum, M. Junker, V. Kazalov, Y. Kermaidic, T. Kihm, I. V. Kirpichnikov,A. Kirsch, A. Kish, A. Klimenko, R. Kneißl, K. T.Kn¨opfle, O. Kochetov, V. N. Kornoukhov, V. V. Kuzmi-nov, M. Laubenstein, A. Lazzaro, M. Lindner, I. Lippi,A. Lubashevskiy, B. Lubsandorzhiev, G. Lutter, C. Ma-colino, B. Majorovits, W. Maneschg, M. Miloradovic,R. Mingazheva, M. Misiaszek, P. Moseev, I. Nemchenok,K. Panas, L. Pandola, K. Pelczar, L. Pertoldi, A. Pul-lia, C. Ransom, S. Riboldi, N. Rumyantseva, C. Sada,F. Salamida, C. Schmitt, B. Schneider, S. Sch¨onert, A-K.Sch¨utz, O. Schulz, B. Schwingenheuer, O. Selivanenko,E. Shevchik, M. Shirchenko, H. Simgen, A. Smolnikov,L. Stanco, L. Vanhoefer, A. A. Vasenko, A. Veresnikova,K. von Sturm, V. Wagner, A. Wegmann, T. Wester,C. Wiesinger, M. Wojcik, E. Yanovich, I. Zhitnikov, S. V.Zhukov, D. Zinatulina, A. Zschocke, A. J. Zsigmond,K. Zuber, and G. Zuzel (GERDA Collaboration), “Im-proved limit on neutrinoless double- β decay of Ge fromgerda phase ii,” Phys. Rev. Lett. , 132503 (2018).[9] J. B. Albert, G. Anton, I. Badhrees, P. S. Bar-beau, R. Bayerlein, D. Beck, V. Belov, M. Breiden-bach, T. Brunner, G. F. Cao, W. R. Cen, C. Cham-bers, B. Cleveland, M. Coon, A. Craycraft, W. Cree,T. Daniels, M. Danilov, S. J. Daugherty, J. Daugh-hetee, J. Davis, S. Delaquis, A. Der Mesrobian-Kabakian,R. DeVoe, T. Didberidze, J. Dilling, A. Dolgolenko,M. J. Dolinski, W. Fairbank, J. Farine, S. Feyzbakhsh,P. Fierlinger, D. Fudenberg, R. Gornea, K. Graham,G. Gratta, C. Hall, E. V. Hansen, J. Hoessl, P. Huf-schmidt, M. Hughes, A. Jamil, M. J. Jewell, A. John-son, S. Johnston, A. Karelin, L. J. Kaufman, T. Koffas,S. Kravitz, R. Kr¨ucken, A. Kuchenkov, K. S. Kumar,Y. Lan, D. S. Leonard, G. S. Li, S. Li, C. Licciardi,Y. H. Lin, R. MacLellan, T. Michel, B. Mong, D. Moore,K. Murray, R. Nelson, O. Njoya, A. Odian, I. Ostrovskiy,A. Piepke, A. Pocar, F. Reti`ere, P. C. Rowson, J. J. Rus-sell, S. Schmidt, A. Schubert, D. Sinclair, V. Stekhanov,M. Tarka, T. Tolba, R. Tsang, P. Vogel, J.-L. Vuilleu-mier, M. Wagenpfeil, A. Waite, T. Walton, M. Weber,L. J. Wen, U. Wichoski, G. Wrede, L. Yang, Y.-R. Yen,O. Ya. Zeldovich, J. Zettlemoyer, and T. Ziegler (EXO-200 Collaboration), “Search for neutrinoless double-betadecay with the upgraded exo-200 detector,” Phys. Rev.Lett. , 072701 (2018).[10] D. Q. Adams, C. Alduino, K. Alfonso, F. T. Avignone,O. Azzolini, G. Bari, F. Bellini, G. Benato, M. Bi-assoni, A. Branca, C. Brofferio, C. Bucci, A. Cami-nata, A. Campani, L. Canonica, X. G. Cao, S. Capelli,L. Cappelli, L. Cardani, P. Carniti, N. Casali, D. Chiesa,N. Chott, M. Clemenza, S. Copello, C. Cosmelli, O. Cre-monesi, R. J. Creswick, A. D’Addabbo, D. D’Aguanno,I. Dafinei, C. J. Davis, S. Dell’Oro, S. Di Domizio,V. Domp`e, D. Q. Fang, G. Fantini, M. Faverzani, E. Ferri,F. Ferroni, E. Fiorini, M. A. Franceschi, S. J. Freed-man, B. K. Fujikawa, A. Giachero, L. Gironi, A. Giu-liani, P. Gorla, C. Gotti, T. D. Gutierrez, K. Han,K. M. Heeger, R. G. Huang, H. Z. Huang, J. John-ston, G. Keppel, Yu. G. Kolomensky, C. Ligi, Y. G. Ma,L. Ma, L. Marini, R. H. Maruyama, Y. Mei, N. Moggi,S. Morganti, T. Napolitano, M. Nastasi, J. Nikkel,C. Nones, E. B. Norman, V. Novati, A. Nucciotti, I. Nu-tini, T. O’Donnell, J. L. Ouellet, C. E. Pagliarone,L. Pagnanini, M. Pallavicini, L. Pattavina, M. Pavan,G. Pessina, V. Pettinacci, C. Pira, S. Pirro, S. Pozzi, E. Previtali, A. Puiu, C. Rosenfeld, C. Rusconi, M. Sakai,S. Sangiorgio, B. Schmidt, N. D. Scielzo, V. Sharma,V. Singh, M. Sisti, D. Speller, P. T. Surukuchi, L. Taf-farello, F. Terranova, C. Tomei, M. Vignati, S. L. Wa-gaarachchi, B. S. Wang, B. Welliver, J. Wilson, K. Wil-son, L. A. Winslow, L. Zanotti, S. Zimmermann, andS. Zucchelli (CUORE Collaboration), “Improved limiton neutrinoless double-beta decay in Te with cuore,”Phys. Rev. Lett. , 122501 (2020).[11] Jonathan Engel and Javier Men´endez, “Status and futureof nuclear matrix elements for neutrinoless double-betadecay: a review,” Reports on Progress in Physics ,046301 (2017).[12] U. van Kolck, “Few-nucleon forces from chiral La-grangians,” Phys. Rev. C , 2932–2941 (1994).[13] P. F. Bedaque and U. van Kolck, “Effective field theoryfor few-nucleon systems,” Annual Review of Nuclear andParticle Science , 339–396 (2002), nucl-th/0203055.[14] E. Epelbaum, H.-W. Hammer, and Ulf-G. Meißner,“Modern theory of nuclear forces,” Rev. Mod. Phys. ,1773–1825 (2009).[15] R. Machleidt and D.R. Entem, “Chiral effective field the-ory and nuclear forces,” Physics Reports , 1 – 75(2011).[16] Vincenzo Cirigliano, Wouter Dekens, Jordy de Vries,Michael L. Graesser, Emanuele Mereghetti, Saori Pas-tore, and Ubirajara van Kolck, “New leading contribu-tion to neutrinoless double- β decay,” Phys. Rev. Lett. , 202001 (2018).[17] Vincenzo Cirigliano, William Detmold, Amy Nicholson,and Phiala Shanahan, “Lattice qcd inputs for nucleardouble beta decay,” Progress in Particle and NuclearPhysics , 103771 (2020).[18] Vincenzo Cirigliano, William Detmold, Amy Nicholson,and Phiala Shanahan, “Lattice qcd inputs for nucleardouble beta decay,” Progress in Particle and NuclearPhysics , 103771 (2020).[19] V. Cirigliano, W. Dekens, J. de Vries, M. L. Graesser,E. Mereghetti, S. Pastore, M. Piarulli, U. van Kolck, andR. B. Wiringa, “Renormalized approach to neutrinolessdouble- β decay,” Phys. Rev. C , 055504 (2019).[20] A.S. Barabash, “Average and recommended half-life val-ues for two-neutrino double beta decay,” Nuclear PhysicsA , 52 – 64 (2015).[21] G. Hagen, A. Ekstr¨om, C. Forss´en, G. R. Jansen,W. Nazarewicz, T. Papenbrock, K. A. Wendt, S. Bacca,N. Barnea, B. Carlsson, C. Drischler, K. Hebeler,M. Hjorth-Jensen, M. Miorelli, G. Orlandini, A. Schwenk,and J. Simonis, “Neutron and weak-charge distributionsof the Ca nucleus,” Nature Physics , 186 (2016).[22] K. Tetsuno, S. Ajimura, K. Akutagawa, T. Bat-purev, W. M. Chan, K. Fushimi, R. Hazama, T. Iida,Y. Ikeyama, B. T. Khai, T. Kishimoto, K. K. Lee, X. Li,K. Matsuoka, K. Matsuoka, K. Mizukoshi, Y. Mori,K. Nakajima, P. Noithong, M. Nomachi, I. Ogawa,H. Ohsumi, K. Ozawa, K. Shimizu, M. Shokati, F. Soberi,K. Suzuki, Y. Takemoto, Y. Takihira, Y. Tamagawa,M. Tozawa, V. T. T. Trang, S. Umehara, K. Yamamoto,S. Yoshida, I Kim, D H Kwon, H L Kim, H J Lee, M KLee, and Y H Kim, “Status of ca double beta decaysearch and its future prospect in CANDLES,” Journal ofPhysics: Conference Series , 012132 (2020).[23] J. M. Yao, B. Bally, J. Engel, R. Wirth, T. R. Rodr´ıguez,and H. Hergert, “Ab initio treatment of collective corre- lations and the neutrinoless double beta decay of Ca,”Phys. Rev. Lett. , 232501 (2020).[24] L. Coraggio, A. Gargano, N. Itaco, R. Mancino, andF. Nowacki, “Calculation of the neutrinoless double- β decay matrix element within the realistic shell model,”Phys. Rev. C , 044315 (2020).[25] Fedor ˇSimkovic, Vadim Rodin, Amand Faessler, andPetr Vogel, “0 νββ and 2 νββ nuclear matrix elements,quasiparticle random-phase approximation, and isospinsymmetry restoration,” Phys. Rev. C , 045501 (2013).[26] J. Barea, J. Kotila, and F. Iachello, “0 νββ and 2 νββ nu-clear matrix elements in the interacting boson model withisospin restoration,” Phys. Rev. C , 034304 (2015).[27] Nuria L´opez Vaquero, Tom´as R. Rodr´ıguez, and J. LuisEgido, “Shape and pairing fluctuation effects on neutri-noless double beta decay nuclear matrix elements,” Phys.Rev. Lett. , 142501 (2013).[28] J. M. Yao, L. S. Song, K. Hagino, P. Ring, and J. Meng,“Systematic study of nuclear matrix elements in neutri-noless double- β decay with a beyond-mean-field covari-ant density functional theory,” Phys. Rev. C , 024316(2015).[29] R. A. Sen’kov and M. Horoi, “Neutrinoless double- β de-cay of ca in the shell model: Closure versus nonclosureapproximation,” Phys. Rev. C , 064312 (2013).[30] J. Men´endez, A. Poves, E. Caurier, and F. Nowacki,“Disassembling the nuclear matrix elements of the neu-trinoless ββ decay,” Nuclear Physics A , 139 – 151(2009).[31] A. A. Kwiatkowski, T. Brunner, J. D. Holt, A. Chaud-huri, U. Chowdhury, M. Eibach, J. Engel, A. T. Gallant,A. Grossheim, M. Horoi, A. Lennarz, T. D. Macdonald,M. R. Pearson, B. E. Schultz, M. C. Simon, R. A. Senkov,V. V. Simon, K. Zuber, and J. Dilling, “New determina-tion of double- β -decay properties in Ca: High-precision Q ββ -value measurement and improved nuclear matrix el-ement calculations,” Phys. Rev. C , 045502 (2014).[32] Y. Iwata, N. Shimizu, T. Otsuka, Y. Utsuno,J. Men´endez, M. Honma, and T. Abe, “Large-scale shell-model analysis of the neutrinoless ββ decay of Ca,”Phys. Rev. Lett. , 112502 (2016).[33] A. Belley, C. G. Payne, S. R. Stroberg, T. Miyagi, andJ. D. Holt, “Ab initio neutrinoless double-beta decay ma-trix elements for 48Ca, 76Ge, and 82Se,” arXiv e-prints, arXiv:2008.06588 (2020), arXiv:2008.06588 [nucl-th].[34] K. Hebeler, S. K. Bogner, R. J. Furnstahl, A. Nogga,and A. Schwenk, “Improved nuclear matter calculationsfrom chiral low-momentum interactions,” Phys. Rev. C , 031301 (2011).[35] Alexander Tichai, Julius M¨uller, Klaus Vobig, andRobert Roth, “Natural orbitals for ab initio no-core shellmodel calculations,” Phys. Rev. C , 034321 (2019).[36] S. J. Novario, G. Hagen, G. R. Jansen, and T. Papen-brock, “Charge radii of exotic neon and magnesium iso-topes,” (2020), arXiv:2007.06684 [nucl-th].[37] G. Hagen, T. Papenbrock, D. J. Dean, A. Schwenk,A. Nogga, M. W(cid:32)loch, and P. Piecuch, “Coupled-clustertheory for three-body Hamiltonians,” Phys. Rev. C ,034302 (2007).[38] Robert Roth, Sven Binder, Klaus Vobig, Angelo Calci,Joachim Langhammer, and Petr Navr´atil, “Medium-Mass Nuclei with Normal-Ordered Chiral NN +3 N In-teractions,” Phys. Rev. Lett. , 052501 (2012).[39] F. Coester, “Bound states of a many-particle system,” Nuclear Physics , 421 – 424 (1958).[40] F. Coester and H. K¨ummel, “Short-range correlations innuclear wave functions,” Nuclear Physics , 477 – 485(1960).[41] Jiˇr´ı ˇC´ıˇzek, “On the Correlation Problem in Atomic andMolecular Systems. Calculation of Wavefunction Com-ponents in Ursell-Type Expansion Using Quantum-FieldTheoretical Methods,” J. Chem. Phys. , 4256–4266(1966).[42] Jiˇr´ı ˇC´ıˇzek, “On the Use of the Cluster Expansion and theTechnique of Diagrams in Calculations of Correlation Ef-fects in Atoms and Molecules,” in Advances in ChemicalPhysics (John Wiley & Sons, Inc., 2007) pp. 35–89.[43] H. K¨ummel, K. H. L¨uhrmann, and J. G. Zabolitzky,“Many-fermion theory in expS- (or coupled cluster)form,” Physics Reports , 1 – 63 (1978).[44] Rodney J. Bartlett and Monika Musia(cid:32)l, “Coupled-clustertheory in quantum chemistry,” Rev. Mod. Phys. , 291–352 (2007).[45] G Hagen, T Papenbrock, M Hjorth-Jensen, and D JDean, “Coupled-cluster computations of atomic nuclei,”Rep. Prog. Phys. , 096302 (2014).[46] John D. Watts, J¨urgen Gauss, and Rodney J. Bartlett,“Coupled-cluster methods with noniterative triple exci-tations for restricted open-shell hartree-fock and othergeneral single determinant reference functions. energiesand analytical gradients,” J. Chem. Phys. , 8718–8733(1993).[47] J. D. Watts and R. J. Bartlett, “Economical triple excita-tion equation-of-motion coupled-cluster methods for ex-citation energies,” Chem. Phys. Lett. , 81 – 87 (1995).[48] I. Shavitt and R. J. Bartlett, Many-body Methods inChemistry and Physics (Cambridge University Press,Cambridge UK, 2009).[49] G. R. Jansen, M. Hjorth-Jensen, G. Hagen, and T. Pa-penbrock, “Toward open-shell nuclei with coupled-clustertheory,” Phys. Rev. C , 054306 (2011).[50] G. R. Jansen, “Spherical coupled-cluster theory for open-shell nuclei,” Phys. Rev. C , 024305 (2013).[51] G. Hagen, M. Hjorth-Jensen, G. R. Jansen, R. Mach-leidt, and T. Papenbrock, “Continuum effects and three-nucleon forces in neutron-rich oxygen isotopes,” Phys.Rev. Lett. , 242501 (2012).[52] G. Hagen, M. Hjorth-Jensen, G. R. Jansen, R. Mach-leidt, and T. Papenbrock, “Evolution of shell structurein neutron-rich calcium isotopes,” Phys. Rev. Lett. ,032502 (2012).[53] Sven Binder, Joachim Langhammer, Angelo Calci, PetrNavr´atil, and Robert Roth, “Ab initio calculations ofmedium-mass nuclei with explicit chiral 3 n interactions,”Phys. Rev. C , 021303 (2013).[54] C. G. Payne, S. Bacca, G. Hagen, W. G. Jiang, andT. Papenbrock, “Coherent elastic neutrino-nucleus scat-tering on Ar from first principles,” Phys. Rev. C ,061304 (2019).[55] H. N. Liu, A. Obertelli, P. Doornenbal, C. A.Bertulani, G. Hagen, J. D. Holt, G. R. Jansen,T. D. Morris, A. Schwenk, R. Stroberg, N. Achouri,H. Baba, F. Browne, D. Calvet, F. Chˆateau, S. Chen,N. Chiga, A. Corsi, M. L. Cort´es, A. Delbart, J.-M.Gheller, A. Giganon, A. Gillibert, C. Hilaire, T. Isobe,T. Kobayashi, Y. Kubota, V. Lapoux, T. Motobayashi,I. Murray, H. Otsu, V. Panin, N. Paul, W. Rodriguez,H. Sakurai, M. Sasano, D. Steppenbeck, L. Stuhl, Y. L. Sun, Y. Togano, T. Uesaka, K. Wimmer, K. Yoneda,O. Aktas, T. Aumann, L. X. Chung, F. Flavigny, S. Fran-choo, I. Gaˇspari´c, R.-B. Gerst, J. Gibelin, K. I. Hahn,D. Kim, T. Koiwai, Y. Kondo, P. Koseoglou, J. Lee,C. Lehr, B. D. Linh, T. Lokotko, M. MacCormick,K. Moschner, T. Nakamura, S. Y. Park, D. Rossi,E. Sahin, D. Sohler, P.-A. S¨oderstr¨om, S. Takeuchi,H. T¨ornqvist, V. Vaquero, V. Wagner, S. Wang,V. Werner, X. Xu, H. Yamada, D. Yan, Z. Yang, M. Ya-suda, and L. Zanetti, “How robust is the n = 34 subshellclosure? first spectroscopy of Ar,” Phys. Rev. Lett. , 072502 (2019).[56] Vincenzo Cirigliano, Wouter Dekens, EmanueleMereghetti, and Andr´e Walker-Loud, “Neutrino-less double- β decay in effective field theory: Thelight-majorana neutrino-exchange mechanism,” Phys.Rev. C , 065501 (2018).[57] Fedor ˇSimkovic, Amand Faessler, Vadim Rodin, Petr Vo-gel, and Jonathan Engel, “Anatomy of the 0 νββ nuclearmatrix elements,” Phys. Rev. C , 045503 (2008).[58] J. Engel, W. C. Haxton, and P. Vogel, “Effective sum-mation over intermediate states in double-beta decay,”Phys. Rev. C , R2153–R2157 (1992).[59] Wick C. Haxton, Kenneth M. Nollett, and Kathryn M.Zurek, “Piecewise moments method: Generalized lanczostechnique for nuclear response surfaces,” Phys. Rev. C , 065501 (2005).[60] M. A. Marchisio, N. Barnea, W. Leidemann, and G. Or-landini, “Efficient method for lorentz integral transformsof reaction cross sections,” Few-Body Systems , 259–276 (2003).[61] M. Miorelli, S. Bacca, N. Barnea, G. Hagen, G. R.Jansen, G. Orlandini, and T. Papenbrock, “Electricdipole polarizability from first principles calculations,”Phys. Rev. C , 034317 (2016).[62] J. Rotureau, P. Danielewicz, G. Hagen, F. M. Nunes, andT. Papenbrock, “Optical potential from first principles,”Phys. Rev. C , 024315 (2017).[63] P. Navr´atil, J. P. Vary, and B. R. Barrett, “Large-basisab initio no-core shell model and its application to C ,”Phys. Rev. C , 054311 (2000).[64] Petr Navr´atil, Sofia Quaglioni, Ionel Stetcu, andBruce R. Barrett, “Recent developments in no-core shell-model calculations,” Journal of Physics G: Nuclear andParticle Physics , 083101 (2009).[65] Bruce R. Barrett, Petr Navr´atil, and James P. Vary, “Abinitio no core shell model,” Prog. Part. Nucl. Phys. ,131 – 181 (2013).[66] Javier Men´endez, Tom´as R. Rodr´ıguez, GabrielMart´ınez-Pinedo, and Alfredo Poves, “Correlations andneutrinoless ββ decay nuclear matrix elements of pf -shellnuclei,” Phys. Rev. C , 024311 (2014).[67] M. Honma, T. Otsuka, B. A. Brown, and T. Mizusaki,“New effective interaction for pf -shell nuclei and its im-plications for the stability of the n = z = 28 closed core,”Phys. Rev. C , 034335 (2004).[68] A. Poves, J. S´anchez-Solano, E. Caurier, and F. Nowacki,“Shell model study of the isobaric chains a=50, a=51 anda=52,” Nuclear Physics A , 157 – 198 (2001).[69] E. Caurier, A. Poves, and A. P. Zuker, “A full 0 (cid:126) ω de-scription of the 2 νββ decay of 48ca,” Physics Letters B , 13 – 17 (1990).[70] A. Balysh, A. De Silva, V. I. Lebedev, K. Lou, M. K.Moe, M. A. Nelson, A. Piepke, A. Pronskiy, M. A. Vient, and P. Vogel, “Double beta decay of Ca,” Phys. Rev.Lett. , 5186–5189 (1996).[71] V. B. Brudanin, N. I. Rukhadze, Ch. Briancon, V. G.Egorov, V. E. Kovalenko, A. Kovalik, A. V. Sala-matin, I. ˇStekl, V. V. Tsoupko-Sitnikov, Ts. Vylov, andP. ˇCerm´ak, “Search for double beta decay of 48ca in thetgv experiment,” Physics Letters B , 63 – 68 (2000).[72] R. Arnold, C. Augier, A. M. Bakalyarov, J. D. Baker,A. S. Barabash, A. Basharina-Freshville, S. Blondel,S. Blot, M. Bongrand, V. Brudanin, J. Busto, A. J.Caffrey, S. Calvez, M. Cascella, C. Cerna, J. P. Ce-sar, A. Chapon, E. Chauveau, A. Chopra, D. Duches-neau, D. Durand, V. Egorov, G. Eurin, J. J. Evans,L. Fajt, D. Filosofov, R. Flack, X. Garrido, H. G´omez,B. Guillon, P. Guzowski, R. Hod´ak, A. Huber, P. Hu-bert, C. Hugon, S. Jullian, A. Klimenko, O. Ko-chetov, S. I. Konovalov, V. Kovalenko, D. Lalanne,K. Lang, V. I. Lebedev, Y. Lemi`ere, T. Le Noblet,Z. Liptak, X. R. Liu, P. Loaiza, G. Lutter, F. Mame-dov, C. Marquet, F. Mauger, B. Morgan, J. Mott,I. Nemchenok, M. Nomachi, F. Nova, F. Nowacki,H. Ohsumi, R. B. Pahlka, F. Perrot, F. Piquemal,P. Povinec, P. Pˇridal, Y. A. Ramachers, A. Remoto, J. L.Reyss, B. Richards, C. L. Riddle, E. Rukhadze, N. I.Rukhadze, R. Saakyan, R. Salazar, X. Sarazin, Yu. Shi-tov, L. Simard, F. ˇSimkovic, A. Smetana, K. Smolek,A. Smolnikov, S. S¨oldner-Rembold, B. Soul´e, I. ˇStekl,J. Suhonen, C. S. Sutton, G. Szklarz, J. Thomas,V. Timkin, S. Torre, Vl. I. Tretyak, V. I. Tretyak, V. I.Umatov, I. Vanushin, C. Vilela, V. Vorobel, D. Waters,S. V. Zhukov, and A. ˇZukauskas (NEMO-3 Collabora-tion), “Measurement of the double-beta decay half-lifeand search for the neutrinoless double-beta decay of Cawith the nemo-3 detector,” Phys. Rev. D , 112008(2016).[73] M. Horoi, S. Stoica, and B. A. Brown, “Shell-model cal-culations of two-neutrino double- β decay rates of Cawith the gxpf1a interaction,” Phys. Rev. C , 034303(2007).[74] C. M. Raduta, A. A. Raduta, and I. I. Ursu, “New theo-retical results for 2 νββ decay within a fully renormalizedproton-neutron random-phase approximation approachwith the gauge symmetry restored,” Phys. Rev. C ,064322 (2011).[75] Mihai Horoi, “Shell model analysis of competing contri-butions to the double- β decay of ca,” Phys. Rev. C ,014320 (2013).[76] A. S. Barabash, “Average and recommended half-life val-ues for two-neutrino double beta decay: Upgrade-2019,”AIP Conference Proceedings , 020002 (2019).[77] P. Gysbers, G. Hagen, J. D. Holt, G. R. Jansen,T. D. Morris, P. Navr´atil, T. Papenbrock, S. Quaglioni,A. Schwenk, S. R. Stroberg, and K. A. Wendt, “Discrep-ancy between experimental and theoretical beta -decayrates resolved from first principles,” Nature Physics ,428–431 (2019).[78] Sabin Stoica and Mihail Mirea, “New calculations forphase space factors involved in double- β decay,” Phys.Rev. C , 037303 (2013).[79] Miroslav Urban, Jozef Noga, Samuel J. Cole, and Rod-ney J. Bartlett, “Towards a full ccsdt model for electroncorrelation,” The Journal of Chemical Physics , 4041–4046 (1985), https://doi.org/10.1063/1.449067.[80] Jozef Noga, Rodney J. Bartlett, and Miroslav Urban, “Towards a full ccsdt model for electron correlation.ccsdt-n models,” Chemical Physics Letters , 126 –132 (1987).[81] C. F. Jiao, M. Horoi, and A. Neacsu, “Neutrino-less double- β decay of Sn, Te, and Xe in thehamiltonian-based generator-coordinate method,” Phys.Rev. C , 064324 (2018). Supplemental MaterialBenchmarks for energies of light nuclei We also computed the ground-state energies for thebenchmark nuclei He → Be, He → Be, He → Be, C → O, and O → Ne. Figure 5 shows the resultsfrom coupled-cluster CCSD and CCSDT-1 computationsand compares them to data for the 1.8/2.0 (EM) inter-action. We remind the reader that this interaction yieldsaccurate binding energies across the lower half of thenuclear chart. As indicated, the coupled-cluster resultsused both the initial and final nuclei as reference states.While deformed reference states were sufficient to matchthe NCSM results for the 0 νββ nuclear matrix elementsshown in Fig. 2 of the main text, the ground-state ener-gies are underbound by a few MeV which are expectedto be obtained when restoring the broken spherical sym-metry. Benchmarks for spectra in the pf shell In addition to the pf -shell benchmark results for thenuclear matrix element shown in Fig. 3 of the main text,we also calculated the low-lying spectrum in Ca and Ti with a spherical Ca Hartree-Fock basis. The re-sults are shown in Figs. 6 and 7, respectively, using theEOM-CCSD, EOM-CCSDT-1, and the EOM-CCSDT-3approximations [48, 79, 80].With triples contributions included, the spectra ofboth, Ca and Ti, agree with the exact diagonal-ization. The closed-shell nucleus Ca is well-describedalready in the EOM-CCSD approximation. As the re-stricted model space does not allow for any 3 p –3 h con-figurations, the spectrum does not change in the EOM-CCSDT-1 or EOM-CCSDT-3 approximations. The nu-cleus Ti is computed with the double-charge exchangeEOM-CC. Despite the quality of both these spectra, thenuclear matrix element (shown in Fig. 3 of the main text)deviates by about ∼ 15% from the exact result. Thisreflects the sensitivity of the this matrix element withrespect to the spectra of the initial and final nuclei [81]. Spectrum of Ti Because of the strong correlation between the accu-racy of the 0 νββ NME and the quality of the excitationspectra of the initial and final nuclei, we calculate theexcitation spectrum of Ti with the double-charge ex-change EOM-CCSDT-3 approximation using a spherical Ca Hartree-Fock basis. The spectrum for the 1.8/2.0(EM) interaction is shown in Fig. 8 and compared withexperiment. | Φ i = | O i | Φ i = | Ne i Exp. -180-170-160-150-140 E ( M e V ) O → Ne . / . | Φ i = | C i | Φ i = | O i Exp. -100-90-80 C → O . / . | Φ i = | He i | Φ i = | Be i Exp. -60-40-20 He → Be . / . | Φ i = | He i | Φ i = | Be i Exp. -60-40-20 He → Be . / . | Φ i = | He i | Φ i = | Be i Exp. -30-25-20-15 He → Be . / . CCSD CCSDT-1InitialFinal FIG. 5. (Color online) Comparison of the ground-state ener-gies for the several light nuclei involved in our 0 νββ bench-mark calculations with their experimental values. The firsttwo columns indicate which nucleus was taken as the refer-ence state, and results from the CCSD and CCSDT-1 approx-imations are shown. The error bars indicate the uncertaintieswith respect to the model-space size. The compressed 2 + and 4 + states of the 1.8/2.0 (EM)spectrum show that the triples correlations in a sphericalbasis are insufficient to represent the deformed nucleusand motivates the need for the deformed coupled-clusterapproach.1 EOM-CCSD EOM-CCSDT-1 EOM-CCSDT-3 FCI − − − − − − E ( M e V ) + -73.660 + -68.300 + -66.382 + -69.922 + -66.952 + -66.804 + -69.414 + -67.714 + -66.113 + -69.131 + -65.77 + -73.66 + -68.31 + -66.38 + -69.92 + -66.95 + -66.80 + -69.41 + -67.71 + -66.11 + -69.13 + -65.77 + -73.66 + -68.31 + -66.38 + -69.92 + -66.95 + -66.80 + -69.41 + -67.71 + -66.11 + -69.13 + -65.77 + -73.66 + -68.39 + -66.85 + -69.93 + -67.19 + -67.07 + -69.40 + -67.72 + -66.29 + -69.13 + -66.21 Ca (GXPF1A) FIG. 6. (Color online) Energies of low-lying states in Cawith respect to the ground states using the GXPF1A in-teractions in the pf-shell compared with full diagonalization(FCI). These results use the EOM-CCSD, EOM-CCSDT-1,and EOM-CCSDT-3 approximations with a spherical CaHartree-Fock reference state (see text for details). The triplesapproximations don’t add any binding energy because thereare no 3 p –3 h configurations for Ca in the pf -shell. EOM-CCSD EOM-CCSDT-1 EOM-CCSDT-3 FCI − − − − − − − − − E ( M e V ) + -15.280 + -10.842 + -14.212 + -13.302 + -12.054 + -13.314 + -12.554 + -11.313 + -12.411 + -11.79 + -16.59 + -12.70 + -15.57 + -14.57 + -13.51 + -14.62 + -13.98 + -13.78 + -13.40 + -16.65 + -12.69 + -15.62 + -14.61 + -13.53 + -14.66 + -14.01 + -13.80 + -13.42 + -17.35 + -16.34 + -15.17 + -14.03 + -15.18 + -14.42 + -14.39 + -13.96 Ti (GXPF1A) FIG. 7. (Color online) Energies of low-lying states in Tiwith respect to the ground states using the GXPF1A in-teractions in the pf-shell compared with full diagonalization(FCI). These results use the EOM-CCSD, EOM-CCSDT-1,and EOM-CCSDT-3 approximations with a spherical CaHartree-Fock reference state (see text for details). Additional νββ decay material The convergence of the NME for the 2 νββ decay of Ca with respect to the the 3 p –3 h truncation, E , is com-puted for the initial nucleus, Ca, the final nucleus, Ti,and the intermediate nucleus, Sc, successively. The lat-ter is shown in Figure 4, and the former two are shownin Figure 9. These calculations utilize the CCSDT-1 ap-proximation in a spherical Ca natural orbital basis withthe 1.8/2.0 (EM) interaction. Not shown is the conver-gence with respect to the number of iterations used in EM (1.8/2.0) Exp.0 . . . . . . . . E ( M e V ) + + + + + + + + + + Ti (EOM-CCSDT-1) FIG. 8. (Color online) Energies of low-lying states in Tiwith respect to the ground states using the 1.8/2.0 (EM)interaction compared with experiment. These results usethe EOM-CCSDT-3 approximation with a spherical CaHartree-Fock reference state (see text for details). the Lanczos (continued fraction) method. Our final re-sults need only 20 Lanczos iterations which convergesvery rapidly and does not contribute to the uncertainty. Q ββ / M ν E ( Ca) = 10 − E ( Sc) = 10 − Exp. E ( Ti) − − − − Ca → Ti . / . E ( Ti) = 10 − E ( Sc) = 10 − Exp. E ( Ca) − − − − FIG. 9. (Color online) The NME for the 2 νββ decay Ca → Ti computed with the Lanczos method and the1.8/2.0 (EM) interaction as a function of the double-betadecay Q-value, Q ββ , and the 3 p –3 h truncation, E , used tocalculate Ca (top) and Ti (bottom). The results use theCCSDT-1 approximation and N max = 10. The experimentalNME and Q-value are are shown in vertical and horizontalbands, respectively. +2 +1 0 -2 -4 -6 -20 -200 ∆(MeV) . . . . . . . . M ν Ca → TiGXPF1a FCICCSDT-1 (Sum)CCSDT-1 (Lanczos) FIG. 10. (Color online) Comparison of the nuclear matrixelement for the 2 νββ decay of Ca → Ti between CCSDT-1and exact shell-model calculation, computed in the pf -shellwith the GXPF1A interaction as a function of the gap be-tween the f / and p / shells, where ∆ = 0 corresponds tothe phenomenological value of the GXPF1A interaction. TheCC results are shown using both the explicit summation fromEq. (4) and the Lanczos method of Eq. (5). To benchmark our 2 νββ decay results of Ca, wecompare NMEs computed with coupled-cluster in theCCSDT-1 approximation with exact results from shell-model calculations in the pf shell phenomenological in-teraction GXPF1A [67]. Figure 10 shows the NME as afunction of the the gap between the f / and p / shells,∆. The original GXPF1A interaction is given by ∆ = 0,and ∆ → −∞ minimizes any correlations, which essen-tially makes the exact shell-model method equivalent tothe approximate coupled-cluster method. Additionally,the Lanczos method is compared to the explicit sumover intermediate 1 + states in Sc as in Eq. (4). Forthese results, the Lanczos method used only 20 itera-tions while the summation used 60 intermediate stateswhich required ∼ 300 iterations. These results confirmthe validity of the Lanczos method and the validity ofthe coupled-cluster method for the 2 νββ NME when im-portant correlations are included.We perform an additional benchmark for the ficticious2 νββ decay of C → O by comparing our results to theno-core shell model in a full model space using the 1.8/2.0(EM) interaction. Both methods use the Lanczos contin-ued fraction method and are converged with respect to N max . Given the relatively small size of the calculations,the CCSDT-1 results include all 3 p –3 h configurations.Additionally, the coupled cluster results are computedin a spherical C natural orbital basis. These results,shown in Figure 9, once again bolster the validity of theLanczos method applied within coupled cluster theory,and shows the importance of including 3 p –3 h configura- tions in these calculations.The shapes of the curves in Figure 11 capture the spec-tra of 1 + states relative to the 1 + ground state in N.The absolute position of these curves is with respect to Q ββ / M ν C → O . / . CCSDCCSDT-1NCSM FIG. 11. (Color online) Comparison of the NME for the 2 νββ decay of C → O computed with the no-core shell model andcoupled cluster at both the CCSD and CCSDT-1 approxima-tions. All results use the Lanczos continued fraction method,and the CCSDT-1 results include all 3 p –3 h configurations andare converged at N max = 10. Q ββ / M ν C → O . / . CCSDCCSDT-1NCSM