Nuclear ab initio calculations of 6He $\beta$-decay for beyond the Standard Model studies

Precision measurements of $\beta$-decay observables offer the possibility to search for deviations from the Standard Model. A possible discovery of such deviations requires accompanying first-principles calculations. Here we compute the nuclear structure corrections for the $\beta$-decay of $^6$He which is of central interest in several experimental efforts. We employ the impulse approximation together with wave functions calculated using the ab initio no-core shell model with potentials based on chiral effective field theory. We use these state-of-the-art calculations to give a novel and comprehensive analysis of theoretical uncertainties. We find that nuclear corrections, which we compute within the sensitivity of future experiments, create significant deviation from the naive Gamow-Teller predictions, making their accurate assessment essential in searches for physics beyond the Standard Model.

The Standard Model (SM) is very successful in explaining the vast majority of observed phenomena in particle physics. Nevertheless, several important questions remain unanswered and the search for evidence of Beyond the Standard Model (BSM) physics is one of the heralds of contemporary particle physics [1]. In particular, recent years have brought advances in the sensitivity of β-decay studies, and several high precision experimental efforts [2,3] have been deployed, as a "precision frontier" to search for BSM physics-alternative to the high-energy frontier represented by the Large Hadron Collider (LHC) [4]. β-decay observables are sensitive to interference of currents of SM particles and hypothetical BSM physics. Such couplings are proportional to v/Λ 2 , with v ≈ 174 GeV, the SM vacuum expectation value, and Λ the new physics energy Email addresses: christian.forssen@chalmers.se (Christian Forssén), doron.gazit@mail.huji.ac.il (Doron Gazit) scale. This entails that a difference between SM theoretical predictions and experiment that can be inferred as a result of a ∼ 10 −4 coupling between SM and BSM physics would suggest new physics at a scale that is out of the reach of current particle accelerators.
However, discovering such minute deviations from the SM predictions demands also high-precision theoretical calculations. In the case of β-decays this challenge entails nuclear structure calculations with quantified uncertainties. The field of ab initio nuclear structure calculations has significantly evolved in the last two decades. Concomitantly, model-independent nuclear interactions based on chiral effective field theory (χEFT) [5,6] have been developed. A showcase example of this progress, relevant for studies of weak interactions with nuclei, is the recent accurate calculations of β-decay rates for nuclei with masses up to A = 100 [7].
Our motivation to study 6 He is twofold: First, 6 He is a light nucleus for which ab initio many-body calculations are numerically tractable, and thus can test feasibility to reach the precision needed by experiments to constrain BSM signatures. Second, 6 He is being studied in several ongoing, or soon to be initiated, experimental campaigns at Laboratoire de Physique Corpusculaire de CAEN (France) [3], National Superconducting Cyclotron Laboratory (USA) [9], the University of Washington CENPA (USA) by the He6-CRES Collaboration [10], and at SARAF accelerator (Israel) [11]. The SARAF experiment is focusing on the angular correlation between the emitted β-particles, while the other aforementioned experiments will measure the β-electron energy spectrum. As the experiments aim for a per-mil level of accuracy, precise calculations based on the SM are needed.
The β-decay spectrum of a pure GT transition takes the simple form: dω ∝ 1 + a βν β ·ν + b F me E [12]. Here, β = k E , while m e , E and k are the mass, energy and momentum of the emitted β-electron, respectively, and ν = νν is the momentum of the emitted anti-neutrino. a βν is the angular correlation coefficient between the emitted electron and anti-neutrino, and b F is the so-called Fierz interference term. The V − A structure of the weak interaction within the SM entails that a βν = − 1 3 and b F = 0 for pure GT transitions [2], neglecting other effects. These values can be modified in the presence of BSM physics, but also by nuclear-structure corrections. Thus, disregarding nuclear structure effects might lead to wrong interpretation of the experiments. Notice that b F is linear in BSM couplings such that constraining b F < 10 −3 in a GT transition yields a sensitivity to BSM currents characterized by tensor coupling T < 1.5 · 10 −4 , or new physics at a scale Λ > 14 TeV [3,4,13] . On the other hand, a βν is proportional to 2 T [14,12], thus demanding higher precision to reach the same BSM constraints. The extraction of both these parameters from 6 He β-decay demands two experimental settings as b F is obtained from spectral measurements in which the angular correlation term vanishes. This is in contrast to the case of first-forbidden unique transitions, where a simultaneous extraction of a βν and b F is possible from the β-energy spectrum [15].
Currently, the state-of-the-art measurement of 6 He angular correlations is a recoil ion energy measurement from 1963 [16], resulting in a βν value consistent with − 1 3 up to the experimental error of 0.9%. Corrections to this result were offered over the years by adding radiative corrections [17], and updating the 6 He shake-off probability [18] and Q-value [19,20]. The present work is, to our knowledge, the first consistent calculation of nuclear-structure related corrections to these observables, taking into account the full nuclear dynamics, shown as important already in Ref. [21], and using χEFT to quantify systematic uncertainties in the nuclear modeling. Nuclear-structure effects are particularly important since they entail a finite b 1 + β − F value which in turn has been shown to distort the extraction of a βν [22]. In the following, we use ab initio calculations of nuclear wave functions and the weak transition matrix elements to predict experimentally relevant observables. We focus on corrections related to nuclear structure, and use a recent β-decay formalism [23], to augment the point prediction with a theoretical uncertainty estimate.
The full expression for 6 He β − -decay differential distribution within the SM-including the leading shape and recoil corrections, i.e., next to leading order (NLO) in GTtakes the following form: where E 0 is the maximal electron energy, k = | k|, L A 1 is the reduced matrix element between the initial-and finalstate wave functions of the rank-1 spherical tensor longitudinal operator of the axial current (proportional to the GT operator). In general, we use the superscript A(V ) to denote axial-(polar-)vector contribution to the weak nuclear current. Furthermore, F − (Z f , E) is the Fermi function (calculated here according to [24]), which takes into account the deformation of the β-particle wave function due to the long-range electromagnetic interaction with the nucleus, and C corr represents other corrections which do not originate purely in the weak matrix element, such as radiative corrections, finite-mass and electrostatic finite-size effects, and atomic effects. These corrections are assumed to be well known [25] and are not taken into account in the following calculations as they do not affect the observables.
Taking recoil and shape corrections into account, the β − ν correlation becomes Similarly, an m e /E spectral behavior appears, similar to the BSM-induced Fierz interference term, via a nuclearstructure dependent factor M J Jl1 is the vector spherical harmonic, and the spher- of the multipole operators, which are expansions of the hadronic currents and charges within the nucleus. At the low energies characterizing nuclear β-decays, this dynamics, microscopically governed by quantum chromodynamics (QCD), can be effectively reduced into a field theory of nucleons, pions and short-range interactions by the use of χEFT. This results in a consistent expansion governed by a small parameter EFT , which dictates the accuracy of the theory. Below we estimate that EFT 0.15 for the present study. A detailed derivation of the powercounting of electro-weak operators in χEFT can be found in Ref. [29] and references therein. Weak probes generally interact with currents of ever growing clusters of nucleons. However, within χEFT, interactions with currents of bigger clusters are suppressed. For weak magnetism,M V 1 , the two-body current part is suppressed by EFT compared to the leading-order single-nucleon current, while theL A 1 and C A 1 two-body current terms are associated with the next order, 2 EFT . We calculate the needed multipole operators within the so-called impulse approximation, i.e., single-nucleon currents weakly interacting with the β particles, while neglecting two-and higher-body currents. In this approximation, the three nuclear operatorsL A ,Ĉ A andM V appearing in Eqs. (1) and (4) can be expressed in terms of four basic multipole operatorsΣ ,Ω ,∆, andΣ [30] (see Appendix A for definitions) aŝ Here, J (M J ) is the multipole rank (projection), m N is the nucleon mass, r j (τ + j ) is the jth nucleon position vector (isospin-raising operator). Note that the sum runs over the A nucleons (not to be confused with the A labeling axial quantities). In (5), µ ≈ 4.706 is the nucleon isovector magnetic moment; [32] are the hadronic vector, axialvector and pseudo-scalar charges, which correspond to the nucleon form factors at q = 0. In general, the nucleon form factors include momentum-transfer corrections proportional to q 2 µ [2]. These, however, are ∼ 1 6 2 qr , and thus smaller than the needed precision.
As aforementioned, the correction ∆E c to the energy transfer E 0 is the difference between the Coulomb energies of the final and initial states of the nucleus, ∆E c ≡ 6 Li 1 + gs |V c | 6 Li 1 + gs − 6 He 0 + gs |V c | 6 He 0 + gs , where V c de-3 notes the full Coulomb potential operator. Quantum Monte Carlo calculations in Ref. [33] present the Coulomb energy difference ∆E c = 0.85 (3) MeV, which is the value that we use in our calculations. This result is consistent with the experimental value for the Coulomb displacement energy between a pair of isobaric analog levels, which is given by [34], where M Z> (M Z< ) is the atomic mass of the higher (lower) Z member of the analog pair, and ∆ nH is the neutron-hydrogen mass difference.
Wave functions of 6 He and 6 Li and the many-body matrix elements of the multipole operators in (4) are obtained within the ab initio no-core shell model (NCSM) [35,36,37] using χEFT interactions as the only input. We utilized two chiral interactions in this work, namely NNLO opt [38] and NNLO sat [39]. The former was constructed from χEFT at the NNLO order with inclusion of only the NN terms. This interaction reproduces reasonably well the experimental binding energies (∼ 5%) and radii for A = 3, 4 nuclei, as well as for the A = 6 systems that are relevant for this work [40]. The NNLO sat interaction is also constructed at the NNLO order of χEFT-but includes 3N forces, and is more accurate for heavier systems [41,42,43,44].
The present calculations are performed using a Slater determinant A-nucleon harmonic-oscillator (HO) basis in the M -scheme. The basis is characterized by the HO frequency Ω and contains up to N max HO excitations above the lowest Pauli-principle-allowed configuration. We apply the standard procedure of introducing one-body transition densities to compute matrix elements of one-body operators between initial-and final-state NCSM wave functions as The operator matrix elements |α| Ô J ( r) |β| , reduced in the angular momentum, are evaluated between HO states which depend on the coordinate r and are labeled by their nonmagnetic quantum numbers |α| (|β|). In Eq. (6),ã |β|,mj = (−1) j β −m β a |β|,−mj , with a † α and a β the creation and annihilation operators for the single-particle HO states |α and |β , respectively, coupled to the angular momentum J. In the present case we have |Ψ i = | 6 He 0 + gs 1 , |Ψ f = | 6 Li 1 + gs 0 , and J = 1. However, the single-particle coordinates r j and r in Eq. (6) are measured with respect to the center of the HO potential. Thus, these matrix elements clearly contain contributions from spurious center-of-mass (CM) motion. Fortunately, the exact factorization of NCSM wave functions into a product of the physical intrinsic eigenstate and a CM state in the 0 Ω excitation makes it possible to remove the effect of the CM completely. Therefore, we  introduce a translationally-invariant one-body density depending on coordinates and momenta measured from the CM of the nucleus, e.g., ξ = − A/(A − 1)( r − R CM ). This density is obtained as a direct generalization of the radial translationally-invariant density [45] by considering dependence on nucleon spins. Specifically, we replace the standard density (6) by with the M J matrix and further details given in Ref. [46]. The "one-body" HO states |a(b) depend on the Jacobi coordinate ξ as opposed to the single-particle HO states |α(β) that depend on single-particle coordinates r.
Results for the nuclear matrix elements of the onebody basic multipole operatorsΣ ,Ω ,∆, andΣ are shown in Fig. 1. These matrix elements are then used to construct the nuclear structure input,L A 1 ,Ĉ A 1 ,M V 1 , as in Eq. (5). The convergence in terms of the basis frequency Ω and the model space truncation N max is well controlled, but it is clear that results depend somewhat on the nuclear Hamiltonian. We also note that the translationallyinvariant one-body density, Eq. (7), and the standard onebody density, Eq. (6), give the same many-body matrix elements at q = 0 for theΣ ,Σ , and∆ operators while the many-body matrix elements ofΩ differ. In particular, the spurious center-of-mass component of the wave functions contaminates the matrix elements when the gradient in the first term ofΩ is applied on the wave function. With an increasing q, all the operators become contaminated by spurious center-of-mass contributions although the effect is quite small forΣ ,Σ , and∆, i.e., it is not visible on the resolution scale of Fig. 1.
The 6 He→ 6 Li nuclear matrix elements that appear in Eq. (4) are shown in Fig. 2. We note that results are indeed sensitive to the removal of spurious CM components as performed in this work. The technology for evaluating these q-dependent multipole matrix elements with NCSM wave functions was developed in Ref. [47] based on the work in [48]. We study the convergence as a function of model space parameters and the dependence on the nuclear Hamiltonian as shown in Fig. 2. These nuclearstructure uncertainties are propagated to the final BSMrelated observables considered in this work, see the dark filled bands in Fig. 3, and shown to be small. Overall, the most sophisticated description is achieved by the N max =12 NNLO sat calculation with the correction of the CM effect and the HO frequency of 20 MeV. At that frequency, the 6 He and 6 Li g.s. energies are at their minimum.
As aforementioned, the lack of two-body currents in the multipole operators leads to an absence in the theory, dominated by a small parameter EFT , where theM V 1 two-body current part is characterized by EFT , whileL A 1 andĈ A 1 two-body current terms are proportional to 2 EFT . To verify this, and estimate these EFT uncertainties better, we make use of other observables, where higher-order EFT calculations were compared to experiment, namely the 6 Li magnetic moment and M1 transition, and the 6 He half-life. According to Ref. [49] ( [50]), the two-body current, which is an NLO correction for this operator, has a vanishing contribution to the 6 Li magnetic moment, and a 20% (10%) contribution to the 6 Li(0 + →1 + ) B(M1) transition. As B(M1) contains the squared matrix element, this entails at most a 10% two-body-current contribution for M V 1 . Additionally, we compared ourL A 1 calculations with the empirical GT operator, |GT 6 He | expt = 2.161(5), calculated from the 6 He half-life [51], and found a deviation of 1.5% inL A 1 , consistent with the fact that these corrections are of higher order in EFT counting thanM V 1 . This consistency allows us to conservatively estimate that   (4). The width of the dark filled bands shows the variation with the employed nuclear Hamiltonian and NCSM model space parameters for HO frequency Ω = 16, 20, 24 MeV, Nmax = 8, 10, 12 (10,12,14) using the NNLOsat (NNLOopt) interaction, using translationally-invariant one-body densities. The width of the light filled band shows the total estimated theory error.
As shown in Fig. 3b(d), we find up to 1% (2%) corrections to the β spectrum (angular correlation), consistent with the a priori estimates based on the small parameters of the problem. However, these corrections depend on the electron kinetic energy, thus extracting a βν requires an energy-weighted average, adhering to the particular exper-imental setup. Here, we exemplify the important effect of this procedure, by using an average of a βν weighted by the spectrum dω 1 + β − dE . In this example, the total correction to a βν due to nuclear structure is i.e., a 7 per-mil correction to the SM a GT βν = − 1 3 . This, however, is a naive value, as one should keep in mind the (often neglected) dependence of the measured a βν value on the b F -analogous term detailed below.
Such a term with a similar spectral behavior as the Fierz interference can be extracted from the corrected spectrum, and our calculations indicate that it is non-zero This result, with an uncertainty of ∼ 10 −4 , is vital for ongoing experiments, aiming for a per-mil level of precision.
In order to extract the β − ν correlation coefficient a βν , one notices that the spectral shape suggests that a measured βν = a βν 1+bF me E [22], resulting in the following relation: where me E = 0.28536 (10). However, a realistic measurement cannot probe directly the correlation between the neutrino and the β particle. For example, in the 1963 experiment [16], the recoil ion energy spectrum was studied, resulting in a different effect. The effect for 6 He is given by a measured βν = a βν + 0.127 b F [22], so the measured value (including radiative corrections [17], influence of the updated shake-off probability [18] and Q-value [19,20]) a measured βν + δ rad,so,Q = −0.3324 (30) should be modified to Thus, the extracted a βν depends on corrections that imitate the spectral dependence of the Fierz term (suppressed by a numerical factor of about 0.1). Importantly, this indirectly induces a linear dependence of this observable upon BSM corrections, beyond the naive quadratic dependence of a βν . Consequently, ∼ 10 −4 experimental precision on this observable would entail tighter BSM constraints [52].
Summarizing, we have used a χEFT framework combined with the ab initio NCSM to analyze the nuclearstructure related corrections to 6 He β-decay observables. In particular, we have studied the angular correlation coefficient and a nuclear structure term with an inverse energy spectral dependence, imitating a Fierz interference term. Our analysis uses the existence of small parameters, originating mainly in the low-energy regime characterizing β-decays, to quantify the relevant theoretical uncertainties. We find that the induced m e /E behavior, that can be wrongly interpreted as a result of Fierz interference between SM and BSM currents, is significantly different than the naive SM value of zero. Our theoretical prediction comes with less than 15% uncertainty. Furthermore, 0.2 per-mil bounds were found for SM nuclear structure effects correcting the angular correlation coefficient. Albeit these are smaller than the current experimental uncertainty, future angular correlations measurements of 6 He decay, aimed at reducing the current error by one order of magnitude, should use these bounds to check for BSM signatures, due to the indirect dependence of the angular correlations on the Fierz term. These results increase significantly the potential to correctly check the SM, as well as pin-pointing possible deviations from it.

Acknowledgments
This work was initiated as a result of the stimulating environment at the ECT* workshop "Precise beta decay calculations for searches for new physics" in Trento. We wish to acknowledge the support of the ISF grant no

Appendices Appendix A. Nuclear multipole operators
The four basic operators from SM electroweak theory that appear in Eq. (5) in the main text are defined as [30] with σ j being the Pauli spin matrices associated with nucleon j. Furthermore, M JM J (q r j ) = j J (qr j )Y JM J (r j ) and ) are the spherical harmonics (vector spherical harmonics), and J (M J ) is the multipole rank (projection).
When evaluating the one-body-like Jacobi-coordinate matrix elements appearing in Eq. (7) in the main text we first carry out the gradients in the parenthesis ofΣ ,Σ (see, e.g., Refs. [30,53]) and then replace r by − A−1 A ξ in all the operators, and, in addition, we replace the gradients (momenta) inΩ and∆ by − A−1 A ∇ ξ . We note that one-body matrix elements of the seven basic multipole operators for electroweak processes can be carried out analytically in the HO basis as demonstrated in Refs. [30,53]. In Ref. [53], a Mathematica script is provided for the calculation of the matrix elements. These results can be readily applied to calculate the matrix elements of the translationally-invariant versions of the operators we use here. In the analytic results, e.g., Eqs. (17)- (19) in Ref. [53], (i) the q is replaced by − A−1 A q, (ii) the matrix elements (18) and (19) in Ref. [53] are multiplied by one more factor of − A−1 A due to the gradient (momentum) acting on the wave function, and, finally, (iii) yet another factor of − A−1 A is applied to terms with 1/q, i.e., Σ ,Σ , and∆ (A.1), to compensate for the extra scaling in step (i).

Appendix B. Comparison to literature
We would like to compare our results to the ones presented in the original 1975 Calaprice calculation [54], which is following the notation of Holstein and Treiman. According to that notation, there are three nuclear form factors needed to describe the beta-dacay transition to first order in recoil: the Gamow-Teller c, the weak magnetism b, and the induced tensor d. These can be connected, at leading order, to the matrix elements we calculated in Eq. (6) (in the main text), through the following leading-order relations [28]: We found that c 1 ∼ = 2.85 (14), in agreement with c ∼ = 2.75 of Calaprice. Also the magnetic form factor b ∼ = 66.7 (8.3) we obtained is in agreement with b = 69.0 (1.0) experimental value that Calaprice presents. As we mentioned, the uncertainty of the magnetic multipole is large because it is dominated by EFT . However, in δ1 + β − a the contribution of the magnetic multipole is averaged out. Last, unlike the Calaprice values, our calculations result in a negative value for d I Ac1 . Using the leading-order relation from Eq. (B.1), we obtain d I Ac1 ∼ = −0.45 (4). If we use the exact same operator as Calaprice, then we get d I Ac1 ∼ = −0.29. This value itself is somehow between the theoretical value 0.12, and the experimental value 2.0 (1.5) presented by Calaprice.