Model identification in $\mu^-\to e^-$ conversion with invisible boson emission using muonic atoms

In this article, we investigate the $\mu^-\to e^-X$ process in a muonic atom, where $X$ is a light neutral boson. By calculating the spectrum of the emitted electron for several cases, we discuss the model-discriminating power of the process. We report the strong model dependence of the spectrum near a high-energy endpoint. Our findings show that future experiments using muonic atoms are helpful to identify the properties of exotic bosons.

Another merit of muonic atoms is that the shape and the nuclear dependence of the electron spectrum are available to obtain detailed information on new physics. The model identification by measuring such characteristic observables has been discussed in another lepton-flavor-violating process, µ − e − → e − e − in a muonic atom [21][22][23]. Despite its importance, no one has studied the model dependence of observables in the µ − → e − X process.
Our goal of this article is to understand the model-discriminating power of the µ − → e − X process in a muonic atom. For a simple discussion of the model dependence, we introduce three effective models in Sec. II. Then, we formulate the rate of µ − → e − X in a nuclear Coulomb potential. In Sec. III, we show numerical results and discuss the model dependence of observables. Finally, we summarize this article in Sec. IV.

II. FORMULATION
In this section, we formulate the spectrum of an emitted electron from the µ − → e − X process in a muonic atom. Here, we assume a boson X lighter than muons. To investigate the model dependence, we consider three simple effective models, called S 0 , S 1 , and V 1 , as follows: First, we assume that X is a scalar field and the effective interaction Lagrangian to charged leptons is given as where P L/R = (1 ∓ γ 5 ) /2 is a projection operator, and g S0 L/R are dimensionless coupling constants. This type of Lagrangian was also analyzed in Ref. [12,18]. In this model, keeping an electron mass m e = 0.510999MeV, we find the rate of the exotic free muon decay µ → eX to be where r X = m X /m µ , r e = m e /m µ , and Multiplying it with the lifetime of muon GeV −2 is the Fermi coupling constant, we obtain the branching ratio for the free muon, Br (µ → eX) = τ µ Γ 0 . For reference, suppose that g S0 L = g S0 R (= g S0 ) and m X = 0. Then, using Br < 2.6 × 10 −6 [13], we obtain the constraint for the coupling constant, Second, we assume the following derivative coupling for the scalar X, where Λ S1 is an arbitrary energy scale to keep coupling constants g S1 L/R dimensionless. The rate of the free muon decay is given as Now we mention that, when both leptons are free and on mass shell, Eq. (5) is effectively equivalent to Eq. (1) due to the Dirac equation, (i/ ∂ − m) ψ = 0. Here, we have the relation of coupling constants given as g S0 L/R = 1 Λ S1 m µ g S1 R/L − m e g S1 L/R .
Applying the relation, we easily prove the equality of Eqs. (2) and (6). However, Eq. (7) no longer holds in a Coulomb potential. For the process in a muonic atom, it is worth investigating quantitative differences of observables between the two models.
Third, in addition to the scalar cases, we consider another case that X is a vector field and the effective interaction is given as where X αβ = ∂ α X β − ∂ β X α is the field strength of the X. The couplings g V1 L/R are dimensionless again due to the arbitrary scale Λ V1 . As with the previous models, the decay rate for free muon is given as Next, we formulate the rate of the µ − → e − X process in a muonic atom. We assume the independent particle model of a muonic atom and an initial muon in a 1s orbit. We define the transition amplitude M as where we take only the leading order of effective interaction. For simplicity, we omit the spin indices. Here, E X and E e are the energies of the emitted X and electron in the final state, respectively. m * µ = m µ − B 1s µN indicates the energy of the bound muon, where B 1s µN is the binding energy between the nucleus and muon in a 1s state. The M connects to the decay rate by The factor of 1/2 comes from the spin average of the initial bound muon. The transition amplitude M includes the overlap integrals of lepton wave functions that are solutions of the Dirac equation with the nuclear Coulomb potential [24,25]. In the central force system, one can represent the wave function of the bound muon as with a normalization condition The angular parts χ are two-component spinors, which is determined analytically. Furthermore, we obtain the radial part and the binding energy by solving an eigenvalue problem for the radial Dirac equations, The nuclear Coulomb potential V C in the equations is given as with a nuclear charge density ρ(r). Here, we use the reduced mass m µN = m µ m N /(m N + m µ ) with a nuclear mass m N . After obtaining the solution where E µ is minimized, we determine the binding energy of the 1s state by For the electron in the final state, it is convenient to use the multipole expansion of the state with momentum p e . The electron scattering state with the incoming boundary condition is expressed as follows: with the Clebsch-Gordan coefficients, (l κ , m, 1/2, s e |j κ , ν), and spherical harmonics, Y m lκ (p e ). Here, κ is a nonzero integer to label partial waves. For the index κ, the total angular momentum j κ and the orbital angular momentum l κ are determined by δ κ is the phase shift of a partial wave labeled by κ. To obtain the radial wave functions for a given E e and κ, we solve d dr The normalization is taken to be Using the expressions of the effective interactions, we find that the electron spectra for the three models (M = S 0 , S 1 , V 1 ) are universally represented as where E X is a function of E e determined by the energy conservation. To take into account nuclear recoil through E X , we apply the well-known prescription as follows [18,27,28]: This additional term represents the kinetic energy of the recoiled nucleus, and the term is sizable only at high E e but negligible at low E e . Thus, even though we do not completely consider the nuclear motion, we believe that this treatment yields a good approximation for any E e . After straightforward calculation, we obtain the explicit formulas for P M κ and P M κ . For M = S 0 , it is found that Here, we define the overlap integral, I where j l is the l-order spherical Bessel function. h indicates the radial wave function of the scattering electron, and H indicates that of the bound muon. This formula for S 0 is consistent with that in Ref. [18]. More complicated expressions for P S1 κ , P S1 κ , P V1 κ , and P V1 κ are given in Appendix A. If we neglect the electron mass, we find that the components of the transition probability satisfy which is valid regardless of M . Due to this symmetry, the cross term between g M L and g M R disappears after summing over κ. This observation is understandable because the interference between left-and right-handed components should vanish for the final electron if m e = 0.
The endpoint energy E mX end of the electron spectrum is kinematically determined as which is obtained by solving the relativistic relation of the energy-momentum conservation. Approximately, Eq. (27) is represented to where the third term is interpreted as the kinetic energy of the recoiled nucleus.

III. NUMERICAL RESULTS
To obtain the radial wave functions of charged leptons and the binding energy of a muonic atom, we solve the differential equations, Eq. (14) for the initial muon and Eq. (19) for the final electron. In solving the differential equations, we use the fourth-order Runge-Kutta method. The correctness of our calculation code is numerically checked by comparing it with the analytic result for a point-charge density.
For reference, we focus on two kinds of nuclei as a target material. One is aluminum, 27 Al, which will be used in the coming COMET and Mu2e experiments. The other is gold, 197 Au, which was used in the SINDRUM II experiment [29]. For both nuclei, we assume the two-parameter-Fermi distribution as the nuclear charge density, given as where Z is the proton number of the target nucleus and e is the magnitude of the elementary charge. The parameters of the distribution, r 0 and a, are given in Table I, and ρ 0 is a normalization factor.By solving Eq. (14), we obtain the values of the binding energy B 1s µN shown in Table I. Substituting the binding energy into Eq. (27), we find the endpoint energy E mX end for an arbitrary m X . The values of E mX end for m X = 0, 25MeV, 50MeV are shown in Table II. The left panel (a) of Fig. 1 shows the electron spectra for the aluminum nucleus. The spectra are normalized by the rate for a free muon, whose expression for each model is given in Sec. II. Here, we plot only the spectrum of S 0 model because the differences between the models are too small to recognize in this energy scale. Each curve in Fig. 1 corresponds to m X , where the electron energy is universally normalized by the endpoint energy for massless X, E 0 end . As well as the endpoint energy, the position of the spectrum peak is shifted to lower as m X is larger. The peak position is approximately given as E e ≃ (m 2 µ − m 2 X )/(2m µ ), which is the expected signal energy if the momentum of the initial muon is assumed to be zero. We also note that the spectrum for 197 Au, shown in the right panel (b) of Fig. 1, has a larger width than 27 Al. This is because the momentum uncertainty of the initial muon is larger as the nucleus has a stronger Coulomb field.
Suppose that simultaneous searches for the µ − → e − X process with the µ − → e − conversion in an experiment which is optimized to detect high-energy electrons. Then, it is useful to focus on the spectrum near the high-energy endpoint. Hereafter we set m X = 0 as a reference because the high-energy endpoint of µ − → e − X is close to the signal energy of the µ − → e − conversion.
We plot the spectra for 27 Al in the range of 0.99 ≤ E e /E 0 end ≤ 1 in Fig. 2. In this figure, one can recognize the difference between the models of the boson X. In particular, the high-energy tail of the V 1 model, indicated by the dotted (green) curve, is larger than the others. This observation suggests that the analysis of the endpoint spectrum is more sensitive to the V 1 model than the others.
We should comment on the spectrum for the S 0 model, shown by the solid (red) curve in Fig. 2. One may find that the spectrum for the S 0 model is unnaturally suppressed near E e /E 0 end ≃ 0.998, which is clearly seen in (b) of Fig. 2. This happens due to the following three facts: First, the spectrum is dominated by the contribution of κ = −1 in Eq. (21). Second, P S0 −1 vanishes when E e = E µ [27]. Third, E µ is slightly smaller than E 0 end due to the finite nuclear mass. Organizing them, we notice that the main contribution of the spectrum vanishes at E e = E µ E 0 end , which is close to but smaller than E 0 end . This interesting property characterizes the S 0 model. In practice, after the confirmation of X, we need much more careful measurement to identify the spectrum shape.
Also, Fig. 3 shows the spectrum for 197 Au in the range of 0.99 ≤ E e /E 0MeV end ≤ 1. We find that the high-energy tail is much larger than 27 Al. As with 27 Al, the tail of the V 1 model is the largest of the three models. We cannot recognize the suppression of the spectrum near the endpoint for the S 0 model in 27 Al case, because the nuclear mass m N is so heavy that E µ is sufficiently close to E 0 end . Finally, we discuss which nucleus is preferable for the µ − → e − X search. Suppose that the new physics search using muonic atoms is performed by measuring the number of electrons with an energy close to the signal energy of µ − → e − conversion, which is equal to E 0 end . We define a net branching ratio as whereτ is the lifetime of a muonic atom, listed in Ref. [31]. This value corresponds to the number of electrons with E e ≥ xE 0 end (x < 1) coming from µ − → e − X, normalized by the created number of muonic atoms. For further convenience, we define so that The Z dependence of R0.9(Z) defined in Eq. (31). Sampled points are shown by crosses. For the simplicity of calculation, we use the uniform distribution with the nuclear radius of 1.2A 1/3 fm as the nuclear charge density. We take the mass number A of the most abundant isotope for each Z [33]. Setting x = 0.9, we find that Z dependence of R 0.9 (Z) is shown in Fig. 4. One can see that the typical value of R 0.9 (Z) is O(10 −9 − 10 −8 ). As larger nuclei, the lifetime of muonic atoms is shorter, but the high-energy tail of the electron spectrum gets larger. Due to the cancellation of the two effects [18], the Z dependence of R is not so strong above Z ≈ 30. Considering the current experimental constraint of Br (µ + → e + X), we find that the current upper limit of the net branching ratio is Br x (Z) < O(10 −15 − 10 −14 ). Since the goal of the created number of muons in the planned µ − → e − conversion searches [19,20] is O 10 18 , it would be possible to reach the constraint by the near-future muon sources.

IV. SUMMARY
We have investigated the µ − → e − X process in muonic atoms as an interesting candidate to constrain the property of light neutral bosons. Assuming three simple effective models of the unknown boson, we have discussed the model dependence of the electron spectrum. As a result, we found that the spectrum near the endpoint strongly depends on the property of the boson X. We also showed that the nuclear dependence of the net branching ratio is moderate.
A remaining theoretical problem is to include radiative corrections in the calculation for the spectrum near the highenergy endpoint, which is shown to be important for ordinary decay of muon in orbit [32]. Although we need further studies for the realistic sensitivity of experiments, we believe that careful measurements for the electron spectrum in a muon decay are useful to find unknown invisible bosons and to identify their property.