Eta Decay and Muonic Puzzles

New physics motivated by muonic puzzles (proton radius and muon $g-2$ discrepancies) is studied. Using a light scalar boson $\phi$, assuming Yukawa interactions, accounts for these muonic puzzles simultaneously. Our previous work limits the existence of such a scalar boson's mass $m_\phi$ from about 160 keV to 60 MeV. We improve this result by including the influence of all of the possible particles that couple to the $\phi$ in computing the decay rate. Doing this involves including the strong interaction physics, involving quarks, necessary to compute the $\eta\pi\phi$ vertex function. The Nambu-Jona-Lasinio model, which accounts for the spontaneous symmetry breaking that yields the constituent mass is employed to represent the relevant strong-interaction physics. We use the $\eta\pi\phi$ vertex function to reanalyze the electron beam dump experiments. The result is that the allowed range of $m_\phi$ lies between about 160 keV and 3.5 MeV. This narrow range represents an inviting target for ruling out or discovering this scalar boson. A possible UV completion of our phenomenological model is discussed.


I. INTRODUCTION
The proton charge radius measured using the Lamb shift in muonic hydrogen, r p = 0.84087 (39) fm [1,2], differs from the CODATA average obtained from hydrogen spectroscopy and e − p scattering, r p = 0.8751 (61) fm [3], by more than 5σ. Although the discrepancy may arise from subtle lepton-nucleon non-perturbative effects within the Standard Model (SM), or experimental uncertainties [4,5], it could also be a signal of new physics involving a violation of lepton universality [6,7]. The muon anomalous magnetic moment provides another potential signal of new physics [8]. The BNL measurement [9] differs from the SM prediction by more than 3σ, ∆a µ = a exp µ − a th µ = 287(80) × 10 −11 [10,11].
A new scalar boson φ (we have concluded in our previous work [12] that other spin-0 and spin-1 bosons are ruled out), which couples to the muon and proton could explain both the proton radius and (g − 2) µ puzzles [13,14]. Phenomenological motivation for such a scalar boson as a Higgs portal in the dark sector has been considered theoretically [15][16][17][18][19][20][21][22] and experimentally [23][24][25][26][27]. We investigate the couplings of this boson to SM fermions, ψ, which appear as Yukawa terms in the effective Lagrangian, L ⊃ e f φψ f ψ f , where f = g f /e, e is the electric charge of the proton, and f is the flavor index. Other authors have pursued this idea, but made further assumptions relating the couplings to different particle species, mass range, etc. We make no a priori assumptions regarding signs or magnitudes of the coupling constants. The Lamb shift in muonic hydrogen fixes µ and p to have the same sign. Without loss of generality, we take both µ and p to be positive, and e and n are allowed to have either sign.
Coupling a single scalar to up and down quarks in an effective Lagrangian at the MeV scale is not consistent with the SU(2) L × U(1) Y gauge invariance of the Standard Model above the electroweak symmetry breaking scale. An ultraviolet (UV) * mestelqure@gmail.com † icloet@anl.gov ‡ miller@phys.washington.edu completion is needed. As is also well known [28], while it is difficult to create a viable model of dark scalars with masses in the MeV range, interesting attempts have been made [20,22]. Electron beam dump experiments have been aimed at searching for new particles [18,[29][30][31][32]. The typical setup of an electron beam dump experiment involves a beam stopped by a large amount of material. The ensuing interactions could produce new particles via a bremsstrahlung-like process. Such particles would pass through a shield region and decay. These new particles can be detected by their decay products, electron and/or photon pairs, measured by the detector downstream of the decay region. In our previous work, it was assumed that the new particle only couples to electrons. In our simple model, considering the φ couplings to other SM particles could dramatically change the exclusion range.
It is worthwhile to study the production of a new scalar boson by eta decay. This is because there are no selection rules preventing φ emission and possible complications involving strangeness are absent in many channels. We will show that η → π 0ψ f ψ f and η → π 0 γγ decay channels are particularly useful. Eta decay to the π 0ψ f ψ f final state is forbidden at tree level in the SM by charge conjugation symmetry, but it is allowed by a virtual φ emission. Eta decay to π 0 γγ is observed, and the existence of the φ may open up new channels, whose decay rate should not exceed the observed value. We will use the Nambu-Jona-Lasinio (NJL) model [33][34][35][36], a chiral effective theory of QCD exhibiting dynamical chiral symmetry breaking, to provide the strong-interaction input necessary to predict these decay rates. The NJL model satisfies the soft-pion theorems making it an ideal tool with which to determine the coupling of a scalar boson to the Goldstone bosons. Therefore, using the current η decay data, we will significantly improve the constraints on the new scalar boson.
A recent experiment extracts the proton radius to be r p = 0.8335(95) fm [37] by measuring the 2S − 4P transition frequency in electronic hydrogen. This result agrees with the previous muonic hydrogen experiments [1,2] but is more than 3 standard deviations away from the CODATA value [3] that is dominated by many previous hydrogen spectroscopy experiments. Three possible scenarios can immediately be envisioned: 1. the proton radius puzzle is solved, 2. it is too early to use the new experiment as a replacement for many others, 3. new physics may coexist with the CODATA value and the new experiment result.
It is tempting to accept the first scenario, however it defies the results of decades of the electron-proton scattering experiments. On the other hand, the preliminary nuclei radii from laser spectroscopy of µ 4 He + and µ 3 He + [38] agree with the electron-nucleus experiments [39,40]. The PRad experiment [41][42][43] may shed some light on this direction. For the second approach, one may argue that the measurement of 2S − 4P transition frequency is very difficult because of quantum interference effects that involve the details of the experimental setup. It is desirable to have a second experiment on regular hydrogen for this transition. Moreover, a more recent electron hydrogen experiment [44] on the 1 S − 3 S transition finds a radius in agreement with the CODATA value and earlier hydrogen spectroscopy measurements. For the third approach, the true value of proton radius may lie within 3 standard deviations of the new experiments and the old CODATA value, and the muonic hydrogen experiments still signal new physics. In other words, the existence of a new scalar meson may not conflict with any of the experiments. We examine the latter two possibilities here. This paper is organized as follows: Sec. II discusses the Lagrangian, introducing φ couplings to u and d quarks. A possible UV completion is discussed. The ηπ 0 φ vertex is discussed in Sect. III . Sec. IV presents the φ and η decay rates. Sec. V revisits the beam dump experiments. Sec. VI and Sec. VII show the new exclusion region obtained by different η decay channels. Sec. VIII discusses third scenario which the new physics coexists with the new regular hydrogen experiments and the old CODATA value. A conclusion is given in Sec. IX.

II. LAGRANGIAN
In our previous work, the Lagrangian involved interactions between the φ and nucleons. This is not sufficient to study effects involving mesons. Coupling between the φ and quarks is examined here.

A. φ couplings to u and d quarks
Here we use a simplified Lagrangian including the new boson φ in the mostly plus metric: where f is the flavor index, f = g f /e, e is the electric charge, and ψ f is the fermion field (quarks and leptons) in the SM. The couplings to the neutron, n , and proton p are given by  Figure 1. Exclusion plot for p (shaded region is excluded). The solid black, dotted red, and dashed blue lines are from our previous work [12,45] corresponding to combining muonic hydrogen [1][2][3] with muon g − 2 experiments [9][10][11], the binding energy difference of 3 He and 3 H [46][47][48][49][50][51][52][53], and the binding energy of nuclear matter per nucleon [54]. The thick yellow solid curve is from the preliminary muonic 3 He ion laser spectroscopy experiment [38][39][40]55] combining with the n constraint in Fig. 2. The vertical line indicates the allowed mass range obtained in Fig. 11.  Since n can take either sign, we present n as a ratio to p . The solid black, dotted red, and dashed blue lines are from our previous work [12,45] corresponding to the low energy scattering of neutron on 208 Pb [56], the preliminary muonic 4 He ion laser spectroscopy experiment [38,40,57], and the laser spectroscopy experiment of muonic deuterium [58][59][60]. The vertical line indicates the allowed mass range obtained in Fig. 11.
The Lamb shift in muonic hydrogen fixes µ and p to have the same sign, therefore, we choose µ and p to be positive, and n and e are allowed to have either sign. From our previous work [12,45], the allowed values of p and n are shown in Figs. 1 and 2 between the solid black, dotted red, and dashed blue lines. Since n can take either sign, we present n as a ratio to p . The allowed values of µ are shown in Fig. 3. We can find the allowed regions of u and d in Fig. 4 and 5, using p and n in Figs. 1 and 2.    . Exclusion plot for u (shaded region is excluded). The region in the black lines is the allowed u obtain from the region between the solid black, dotted red, and dashed blue lines in Figs. 1 and 2 from our previous work [12,45]. The vertical line indicates the allowed mass range obtained in Fig. 11.

B. A Concrete UV Model
Coupling a single scalar to u and d quarks in an effective Lagrangian at the MeV scale is not consistent with the SU(2) L × U(1) Y gauge invariance of the SM above the electroweak symmetry breaking scale. An explicit model that specifies the SU(2) L × U(1) Y quantum numbers of new particles and the Lagrangian above and below the weak scale was constructed in the light "Dark Higgs" solution to the (g − 2) µ puzzle provided in Ref. [20]. In their work, the low energy Lagrangian (that of our Eq. (1)) appears as their low-energy theory. (We need only the CP-even dark Yukawa coupling here.) Ref. [20] provided a possible UV completion of the low energy effective theory of Eq. (1). In this framework, all new particles are assumed to be charged under U(1) d with the same dark charge so only their SM charges are identified. They let X , where = e, µ, τ is a flavor index, be vector-like fermions with the quantum numbers of right-handed SM leptons R (i.e. SU(2) singlets), and masses m X few × 100 GeV. They also introduce a new Higgs scalar doublet H d and a complex scalar singlet φ. It is assumed that H d and φ have nonzero vacuum expectation values (vevs) which spontaneously break U(1) d . These ingredients can be motivated within a dark Z model [61]. Then they postulate the following SM × U(1) d invariant interactions where L and H refer to SM lepton and Higgs doublets, respectively. The above interactions respect lepton flavor conservation up to soft breaking by (small) off-diagonal masses m X , which are assumed to be the only sources of lepton flavor violation.
Ref. [20] illustrates how the above can be realized in a model with flavor symmetries that allow for a realistic neutrino mass matrix. A vacuum expectation value for H d followed by charged lepton mass matrix diagonalization could result in misaligned φ and H lepton couplings which lead to interesting consequences, as discussed in their paper. Their model also includes scalar coupling to quarks. In that case the H and H d alignment with the mass matrix is maintained and flavor changing current constraints are avoided at the tree level.

III. THE ηπ 0 φ VERTEX
To use η decay to constrain φ we need to know the ηπ 0 φ vertex: g ηπ 0 φ to compute the η decay rate and constrain the possible values of u,d . To obtain the relevant vertex, the π 0 Figure 6. Diagrams that contribute to the η → π 0 φ form factor. The single line is the dressed-quark propagator, the solid-dot is the φ-quark vertex, and the crossed-circle is the appropriate Bethe-Salpeter vertex.
and η quark-model wave functions are used the coupling g reads where the overall factor of 2 takes into account both the q and q contributions; the dimension of f is that of mass. Recall that f = g f /e. Since φ actually couples to the quarks, the ηπ 0 φ vertex should be described by a form factor that accounts for meson structure. We use the Nambu-Jona-Lasinio (NJL) model [33][34][35][36] to estimate this form factor. The piece of the three-flavor NJL Lagrangian relevant to this problem reads and λ a are the Gell-Mann matrices with λ 0 = 2 3 1. The ηπ 0 φ form factor in Fig. 6 is found to be where N c is number of colors, the π 0 and η Bethe-Salpeter vertices [62] are given by Γ π = √ Z π γ 5 λ 3 and Γ η = Z η γ 5 λ 8 1, where M, M s are the dressed quark masses. These are generated by the spontaneous symmetry breaking famously (Nobel Prize 2008) inherent in the NJL model. Performing the flavor space 1 In this calculation we ignore η − η mixing. trace in Eq. (7), and using the proper-time regularization scheme [62][63][64][65] gives where q = p − p and The integral over proper-time includes both an infrared and ultraviolet cutoff, namely, τ ir = 1/Λ 2 UV and τ uv = 1/Λ 2 I R , where the former implements aspects of quark confinement [62]. Since we are interested in the decay rate of η → π 0 γγ or η → π 0ψ f ψ f process, the allowed values of q 2 range from 0 to −(m η − m π ) 2 if η is at rest (the lower bound of q 2 is exact for π 0 γγ final state and a good approximation for π 0ψ f ψ f final state with light fermions).
To constrain f (q 2 ) we consider dressed quark masses in the range 200 MeV < M < 400 MeV, with values outside this range deemed unlikely based on previous studies of meson and nucleon properties [62,66]. For each dressed-quark mass the parameters Λ UV and G π are adjusted so that the empirical values of the pion mass and decay constant are reproduced, and Λ I R should be of the order of Λ QCD because it represents the effects of confinement in the model, and we chose Λ I R = 240 MeV. Our results are summarized in Table I, where the Bethe-Salpeter vertex normalizations, Z η and Z π , are outputs of the calculation.
The range of values for f (q 2 ) are summarized in Table I. The form factor f (q 2 ) can take values from 551 to 1274 MeV including momentum dependence and model dependence. Within our range of interest, we can drop the momentum dependence and treat f as a coupling constant, with a magnitude and uncertainty given by where the central value of f is chosen to be 825 MeV, such that the associate parameters give the best description of pion and kaon system [66]. We will use the lower bound of f to constrain the scalar boson φ. One can also consider η − η mixing [67][68][69][70][71]: The influence of η−η mixing implies the following replacement in the previous equations: where θ ranges from −17 • to −12 • which corresponds to approximately a 30% change in f (q 2 ). However, η − η mixing is also model dependent in our calculation, therefore we will include the effect of η − η mixing in the estimate of the model-dependent uncertainties already given in Eq. (11). The range of f (q 2 ) corresponds to q 2 from 0 to −(m η − m π ) 2 . Dimensionful quantities are in units of MeV with the exception of G π which is in units of GeV −2 .

A. φ decay rate
If m φ > 2m e , φ can decay to two fermions and the decay width is where α is the fine structure constant and m f is the fermion mass. If φ decays to two photons through a fermion loop where the superscript f of Γ indicates the fermion in the loop; Q f is the electric charge of the fermion in the units of e, e.g. Q u = 2/3; I is obtained in Ref. [72] and reads where The total φ decay rate is where Γ interference φ→γγ is the interference of different fermions loops. We use the constituent quark mass (m u = m d = 200 MeV ) for the decay through quark loop, because the relevant scale of the φ decay process is m φ and the quark mass should accordingly be evolved to this scale. The result is shown in Fig. 7. Since the interference contribution is expected to be smaller than the leading one, we neglect Γ interference  Figure 7. The total φ decay width (shaded region is allowed to decay). The decay width with a superscript f indicates the process decaying through f fermion loop. The thin solid blue, dashed green, and dotted red lines are φ → γγ through u quark, d quark, and muon loops with couplings from Figs. 4, 5, and 3, respectively. The thick solid cyan and dotted dashed yellow lines are φ → e + e − and φ → γγ through electron loop with the coupling e from Fig. 9.

B. η decay rate
The total decay width of η, Γ η,total , is 1.31±0.05 keV [73]. Since the new scalar particle φ decays into two fermions or two photons final state, η decay may include the process where X can be either fermion or photon. The three body final state phase space integral is [74] DΦ 3 = m 2 η 64π 3pp 3 d cos θds 13 (20) wherep,p 3 , and s 13 are dimensionless and given bỹ x 2 = m 2 π 0 /m 2 η and x 1 = x 3 = m 2 X /m 2 η ; cos θ is the polar angle of p 3 with respect to p. Therefore, the three body final state decay rate is s mi n 13 ds 13 |M | 2pp 3 (24) where S is the symmetry factor taking into account how many identical particles in the final state, s min 13 = ( √ x 1 + √ x 3 ) 2 , and

V. BEAM DUMP EXPERIMENTS
In our previous work [12,45], we considered the constraints of beam dump experiments [29,[75][76][77]. However, in making the e exclusion plots we only included φ → 2e, φ → 2µ, and φ → 2γ as proceeding through the electron loop. Since µ , u and d are much larger than e , we should include φ → 2γ through muon, u quark, and d quark loops as well. Further investigation and recalculation shows that the exclusion plots change quite a lot, see Fig. 8.
The changes of the exclusion plots are easy to explain. The coupling of u and d become bigger in the large mass region as well as the decay width of φ → 2γ through quark loops. The decay length of φ become shorter than the thickness of the shield so the exclusion stops when m φ 30 MeV. On the other hand, Γ u φ→γγ becomes dominant when m φ > 10 MeV. The constraint for e need to be smaller so that φ is harder to produce from beam dump. Therefore the lower bound of the constraint is lower. The result is also shown in Fig. 9.
We emphasize that if the intermediate state φ is a real particle, its decay products are extremely difficult to detect since φ is long-lived. Therefore we consider the virtual scalar boson decay to two fermion channel. The amplitude is (25) and summing the square of the amplitude over final states gives  Figure 9. Exclusion plot for | e | (shaded region is excluded). The thin blue, thick red, thin dashed yellow, and thick dashed green lines are from our previous work [12,45] corresponding to beam dump experiments (see, Fig. 8), electron anomalous magnetic moment (g − 2) e [78,79], Bhabha scattering [80], and the Lamb shift of hydrogen [81][82][83][84]. A1 at MAMI [85] and BABAR 2014 [86] constraints are in the upper right corner. The dotted black line is our new result from η → π 0 e + e − decay. The vertical line indicates the allowed mass range obtained in Fig. 11.
In Fig. 10, using Eq. (27) we show that the new channel to muon pair through a virtual φ is much smaller than the the SM constraint, therefore no new constraint is obtained from this channel.
The decay of the virtual φ to electron-positron pairs, φ * → e + e − would contribute to this rate.  Figure 10. Decay rate of η → π 0 φ * → π 0 µ + µ − (shaded region is allowed to decay). Since the decay rate of this channel is much smaller than the SM constraint, there is no new constraint obtained from this channel.
We have to handle this process with more care. If m φ > 2m e , the virtual φ propagator can be on-shell and we need to put in the total φ decay width We can further use the narrow width approximation (NWA) Assuming m φ > 2 m e the decay rate becomes Using Eqs. (27) and (32), we show the exclusion of e in Fig.  9.
VII. ETA DECAY TO TWO PHOTONS (η → π 0 φ * → π 0 γγ) Based on the same reason argued in Sec. VI, we only consider the intermediate state φ is virtual. The amplitude η → π 0 γγ through a virtual φ and then a fermion loop with flavor f (two diagrams) is where Q f is the fermion electric charge in units of e, e.g. Q d = −1/3; λ µ is the photon polarization vector; the term in the square bracket includes the contribution of two Feynman diagrams; I is defined in Eq. (16).
The amplitude squared (assuming there is only one fermion loop) after summing over final states is Using Eq. (24), the decay rate is found to be We used the fact that there are two photons in the final state (S = 2), s min 13 = 0, andp 3 = √ s 13 /2. We can further apply the narrow width approximation The process η → π 0 γγ is observed and the value is The decay rate of a virtual φ to two photons cannot be too big to spoil the observed value. We define that the scalar boson is excluded if the decay rate is greater than the observed value plus 3σ, i.e. Γ(η → π 0 γγ)/Γ η,total > 3.22 × 10 −4 . In Fig. 11, we show a virtual φ decay to two photons through different fermion loops. We can read from the plot that the allowed scalar boson mass is from 168 keV (from our previous work) to 3.45 MeV (from muon loop decay channel). Again, we neglect the interference terms because they are expected to be smaller than the leading term. The allowed u , p , and µ are all around 10 −3 . One might think that the existence of the scalar boson could survive the constraints if its coupling to the u quark (or muon) were zero. However, investigating the allowed regions of the parameter spaces, we find that this is not the case. The exclusion of φ coupling to u quark (muon) is equivalent to excluding the existence of the φ.

VIII. COEXISTENCE WITH THE OLD AND THE NEW
To allow the new physics, the electron-proton scattering (CODATA value), and the new experiment to coexist, we assume that the actual proton radius lies within 3σ of both new and old experiments. Such a value is between 0.8568 to 0.8620 fm. First, following our previous work, we obtain the constraint of p , n , µ , u , and d , and the results are shown in Figs. 12-16, respectively. Second, we repeat the analysis in Secs. VI and VII, and obtain the φ decay rate in different channels in Fig. 17, e exclusion plot in Fig. 18, η → π 0 µ + µ − decay rate in Fig. 19, and η → π 0 γγ decay rate in different channels in Fig. 20. Finally, we obtain the allowed m φ is from 168 keV to 2.50 MeV and u , p , and µ are all around 10 −3 . Note that n is completely ruled out with the new result and this means that the scalar coupling to the neutron is zero, i.e. u = −2 d . At first glance, one might think that the effect of new physics is approximately halved, therefore the constraint should be less strict. This intuition is not correct. Comparing with the results in Sec. VII, we see that the upper bound of m φ becomes smaller. The problem is that the effect of n should be included carefully. To correctly analyze it, examine Eq. (36). The η → π 0 γγ decay rate is proportional to ( u − d ) 2 (this factor comes from the ηπ 0 φ vertex). In the previous scenario, d can be positive (see Fig. 5), so ( u − d ) 2 can be smaller than 2 u ; in this section, d stays negative in the allowed m φ range (see Fig. 16), so ( u − d ) 2 is bigger than 2 u . The root cause of this strange behavior can be traced back to n can be non-zero in the previous scenario, but strictly zero in this scenario.

IX. CONCLUSION
Our previous work limits the existence of the new scalar boson to the mass range of about 160 keV to 60 MeV. Here we reanalyze the beam dump experiments and find that the exclusion region is quite different than that of the previous work    Δa μ Figure 14. Exclusion plot for µ (shaded region is excluded). See caption of Fig. 3 for more details. [45]. With strong interaction input to η decay from the NJL model, we present a tighter constraint on the new scalar boson φ: The mass range is now from 160 keV to 3.5 MeV, u , p , and µ are all around 10 −3 , n is from −0.4 to 0.2. We also considered the scenario where the new physics    Figure 17. The total φ decay width (shaded region is allowed to decay). See caption of Fig. 7 for more details.
coexists with the new regular hydrogen laser spectroscopy experiment and the old CODATA value. Most of the parameter space is similar, except n is zero, meaning u = −2 d .  Figure 18. Exclusion plot for e (shaded region is excluded). See caption of Fig. 9 for more details.  Figure 19. Decay rate of η → π 0 φ * → π 0 µ + µ − (shaded region is allowed to decay). See caption of Fig. 10 for more details.  . Decay rate of η → π 0 φ * → π 0 γγ (shaded region is allowed to decay). See caption of Fig. 11 for more details.
One might expect that if we change the quark masses used in loop calculation in Sec. IV A and VII, the result might change drastically. However, this is not the case. Although the allowed decay channels in Fig. 11 and 20 change accordingly when varying u and d quark masses, the resulting upper bound of m φ is always around few MeV. After thorough investigation, we find the constraint of upper bound of m φ is quite robust.
There are experiments aiming to explore η decay with higher precision, such as recently approved the JLab Eta Factory (JEF) experiment [87], and proposed and the REDTOP project [88]. There are several experiments that study the muonic puzzles: The MUSE experiment [89] plans to investigate the proton radius puzzle by measuring µ ± and e ± − p elastic scattering at low energies. The muon g − 2 experiments at Fermilab [90] and J-PARC [91] are of high interest for their own sake and for their bearing on the proton radius puzzle. The COMPASS collaboration is planning a radius measurement using their muon beam [92].
The present work severely constrains the mass range of the possible new scalar boson φ. Measurements aimed at investigating this particle for masses between 160 keV and 3.5 MeV could either discover the particle or completely rule it out. Our constraints are based on the assumption that only φ couples to SM particles through a simple Yukawa coupling. In the case of ruling out the scalar boson, this assumption becomes a constraint for model construction as an explanation for muonic puzzles. Indeed, there might be a subtle symmetry to forbid Yukawa terms, or there might be more complicated mechanisms for new physics to interact with the Standard Model.