Standard Model estimate of $K^+\to\pi^+4e$ branching ratio

The branching ratio of the $K^+\to\pi^+e^+e^-e^+e^-$ ($K^+\to\pi^+4e$) decay is calculated at leading order in the Standard Model. The dominance of the neutral-pion pole determines the overall branching ratio to be $B(K^+\to\pi^+4e)=B(K^+\to\pi^+\pi^0)B(\pi^0\to4e)\approx7.0(3)\times10^{-6}$. The significance of this contribution is very much concentrated in the context of the whole available phase space, throughout most of which the one-photon-exchange topology is prevalent in turn. It is thus interesting to present branching ratios for only parts of the allowed kinematical region. We find, for instance, $B(K^+\to\pi^+4e,\,m_{4e}>150\,\text{MeV})=6.0(6)\times10^{-11}$.


I. INTRODUCTION AND SUMMARY
The flavor-changing (or, in particular, strangenesschanging) neutral-current weak transitions are absent at tree level in the Standard Model (SM). At the same time, they are manifest, among others, in radiative nonleptonic kaon decays like K + → π + + − (γ), = e, µ, and have become an interesting probe of SM quantum corrections and beyond. The related underlying radiative modes (transitions) K + → π + γ * (γ) have been studied before: They have been calculated in Chiral Perturbation Theory (ChPT) [1][2][3], enriched with electroweak perturbations [4,5], at leading order (LO) (at oneloop level) and beyond, including the dominant unitarity corrections from K → 3π [6,7]. The long-distancedominated K + → π + γ * transition turns out to be essential also for the decay mode to which this work is dedicated: K + → π + 4e. Moreover, one needs to consider the K + → π + γ * γ * transition, i.e. when both photons are off-shell, which entails a significant challenge.
Out of all possible contributions to K + → π + 4e, the neutral-pion exchange clearly dominates when the π 0 becomes on-shell.
The overall branching ratio for this 5-track decay is then saturated by the contribution of the associated narrow π 0 peak and can be directly determined as B(K + → π + 4e) = B(K + → π + π 0 )B(π 0 → 4e). (This has also been checked here.) It then becomes challenging to observe the K + → π + 4e process away from m 4e = M π 0 , i.e. when the leptons do not originate from the decay of the (nearly) on-shell π 0 . However, it is exactly the suppressed decay rate that makes such reactions attractive to study, at least from the point of view of beyond-Standard Model (BSM) physics. It is clear that to identify a new-physics-scenario contribution taking part in a decay, one should have a rough estimate of the SM rate. Only then the possible new-physics effects can be spotted and become visible as deviations from such SM predictions.
Related to K → π4e decays, however, there has not been any number presented in the literature which would be more than an order-of-magnitude estimate. In particular, for instance, it is naturally believed that it is unlikely for these branching ratios to exceed the O(10 −10 ) benchmark [8], since they are suppressed with respect to e.g. B(K + → π + e + e − ) ≈ * tomas.husek@thep.lu.se 3 × 10 −7 or B(K + → π + γγ, non-res.) ≈ 1 × 10 −6 simply due to phase-space factors and additional QED vertices by O(α 2 ).
(Indeed, in the present work, we find that the non-resonant topologies give rise to B(K + → π + 4e, non-res.) = 7.2(7) × 10 −11 .) This leads to increasing attention of theorists, and possible BSM scenarios are being explored. For example, in Ref. [8], a model is introduced in which the K → π4e decays proceed via K → π(X → XX) intermediate states, with a cascade of dark-sector particles X ( ) and underlying dynamics potentially significantly enhanced compared to the SM case. In turn, the searches then follow in suitable experiments, and more precise knowledge of SM background, ideally at the level suited for Monte Carlo (MC) implementation, becomes essential. This is the main motivation behind the present work. Besides the expressions for the matrix elements, the squares of which are the vital input for the MC event generators, we also provide here the branching ratios related to the respective contributions for reference.
At present, the K + → π + e + e − single-event sensitivity at NA62 [9] is O(10 −11 ) [10]. The acceptance for the K + → π + 4e channel is likely to be lower, leading to O(10 −10 ) or even O(10 −9 ) expected sensitivity, i.e. up to two orders of magnitude above the SM estimate presented here. Thus, unless a possible BSM signal would manifest itself inside the detector, a SM signature of a K + → π + 4e decay outside the m 4e = M π 0 kinematical region should not be observed beyond statistical fluctuations.
Further work is currently in progress: It is planned to add the neutral channel, i.e. to treat the branching ratios for K 0 → π 0 4e, and to include the expressions for matrix elements squared.
In Sec. II, we first present the contributions that appear at LO in the Standard Model and then provide a brief description of how the relevant matrix elements are obtained. Some suitable ChPT Lagrangians employed to calculate the pion-pole contribution are included. Section III then combines and compares all the terms, and provides numerical results for respective contributions to the branching ratio, obtained in terms of a simple Monte Carlo method. Some of the details are deferred to Appendix A so as not to interrupt the flow of the section. The one-fold differential decay widths as functions of the 4-lepton invariant mass m 4e , as well as the branching ratios, are presented in Fig. 3 and Table I At lowest order in the QED expansion (at O(α 2 )) and in ChPT (at O(G F ), or, equivalently, O(G 8 ) and O(G 27 )), there are 3 distinctive underlying topologies present on the meson side of the amplitude: 1. K + → π + γ * conversion, as shown in Fig. 1, which we will call the single-photon exchange (between the meson and lepton parts); 2. K + → π + γ * γ * conversion, as shown in Fig. 2, which we dub the two-photon exchange, manifesting itself as two sub-topologies: (a) K + → π + γ * conversion with an extra radiative photon coming either from the meson legs or from the vertex itself, as shown in Fig. 2a, (b) K + → π + π 0 * , π 0 * → γ * γ * , as shown in Fig. 2b.
A few remarks are in place here. Regarding the topology (2a), when one of the two photons becomes on-shell, such a contribution is an important part of the radiative corrections in K + → π + + − decays [11,12] and is naturally considered when the K + → π + γ * transition form factor is measured. The topology (2b) might become important for the K + → π + e + e − (γ) decay and the related form-factor extraction, and as a background for the measurement of the branching ratio of the rare decay π 0 → e + e − . The most nontrivial ingredient of topology (1), the K + → π + γ * conversion, has already been extensively studied in the literature, as mentioned in the introduction. Its contribution can be represented via a single form factor, for which several parametrizations are being used. In what follows, we denote this form factor simply as F (s), with s being the virtuality of the photon. It is related to the standard W + (z), introduced in Ref. [7], through For the parameters in Eq. (2), we employ a + = −0.584 (8) and b + = −0.700(35), the weighted average of available experimental inputs [13][14][15]. For the (longdistance-dominated) transition matrix element, we can then simply write, based on the gauge and Lorentz symmetries, note that P − r = k was retained for the manifestly gauge-invariant form. The Lorentz structure of the transition is further simplified when coupled to a conserved current, such as the lepton electromagnetic current: Figure 1: The one-photon-exchange topology. There is a cross diagram where the additional off-shell photon is radiated from the positron line (obtained for instance by changing the direction of the fermion-number flow). These two then serve as a gauge-invariant building block that comes in 4 permutations of external legs. At tree level, the conversion γ * → 4e can be resolved by drawing 8 diagrams in total: When the photon converts into an electron and positron, each of these can radiate another off-shell photon converting into another e + e − pair; see Fig. 1. This accounts for 2 diagrams that form a gauge-invariant building block. There are then 4 independent ways how to assign momenta to the leptons (2 for electrons and 2 for positrons). For every such interchange, one needs to take into account a relative sign for crossing the electron lines (anticommuting two fermion fields). The gauge-invariant building block in question is with and J α (p, q) ≡ū (p)γαv(q) (p+q) 2 . The lepton part of the amplitude then amounts to and the overall amplitude for the topology (1) reads where tilde denotes the matrix element with e = 1.
In principle, the analogous techniques employed for the calculation of the K + → π + γ * γ within the ChPT framework can be used to obtain the form factors accompanying the Lorentz structures of the general twophoton-conversion matrix element. In this work, only a simplified expression is utilized, based on the approximate model used originally mainly for the K + → π + γ * γ transition. The two-photon transition of topology (2a) can thus be written as This form, dependent on a single form factor (the same one as in Eq. (3)), is particularly useful for the case when it provides the contribution to radiative corrections for the K + → π + + − decay as part of the process of measuring F (s). There, one takes one of the photons onshell. In the soft-photon regime, such an approximation is justified. For the hard photons, the free parameter |κ| 1 is introduced to cover the model uncertainty. In practice, the desired physical results do not seem to be sensitive to this parameter. Further details can be found in Ref. [12]. For our purposes, we assume that the form in Eq. (9) is good enough as an order-of-magnitude guess, at least. Since it turns out to be numerically negligible compared to the topology (1) (with a relative suppression of one order of magnitude), we will not discuss it in greater detail here. Finally, we look at the topology (2b). The K + π + π 0 vertex can be obtained, in the framework of ChPT, from the Lagrangian [16,17] with F ≈ 92 MeV the pion decay constant, L µ ij = i(U † ∂ µ U ) ij and U = exp iΦ F . Above, the mesons are represented by the matrix Φ = 8 a=1 λ a φ a , with φ a being the eight meson fields. The LO matrix element for the K + → π + π 0 * (q) transition following from Lagrangian (10) is On-shell and in the isospin limit, this, of course, becomes the standard result for the charged amplitude 1 Besides the chiral SU(3) 27-plet contribution presented in Eq. (10), there is also the octet Lagrangian L ∆S=1 c., which for K + → π + π 0 leads at LO to (isospin-breaking) terms proportional to (M 2 π + − M 2 π 0 ) or the π 0 -η mixing angle (2) ≈ 1%. Such a contribution is therefore numerically small compared to the G 27 term, as G 8 /G 27 20 [17], and it is neglected here.
Since B(K + → π + π 0 ) ≈ 20.7 % [18], Γ K + M K + ≈ 1.08 × 10 −16 and the phase space amounts to nothing else than 1 , we find |G 27 | ≈ 0.53 TeV −2 . It can also be rewritten in terms of the SM parameters involved as G 27 = − 1 2 1 v 2 V ud V * us g 27 , with |V ud | ≈ 0.974, |V us | ≈ 0.224, and v ≈ 0.246 TeV [18]. Hence, we get |g 27 | ≈ 0.29, in agreement with Refs. [17,19,20]. For completeness, for the virtual η we would have The π 0 γγ vertex stems from the Wess-Zumino-Witten (WZW) term [21,22], which can be reduced for our application to the form Above, we also include the η term and one thus have (together with Eq. (13)) all the ingredients to include the virtual-η-exchange contribution, too. However, this contribution can be safely neglected compared to the π 0 case and we will not consider it from now on: It does not affect the shape of the overwhelming π 0 peak and the correction to the tail within the kinematically allowed region is rather significant but only cosmetic in the global picture when it comes to the branching-ratio estimate.
In the pion case, the Lagrangian (14) leads to the LO transition amplitude the form of which is exact for on-shell pion and photons. For on-shell pion and off-shell photons, this can be further modulated by the doubly off-shell neutral-pion electromagnetic transition form factorF(k 2 1 , k 2 2 ), witĥ F(0, 0) = 1. In particular, As an example, see, for instance, the model described in Ref. [23]. Further complications arise for the off-shell pion, but due to the pion pole enhancement, only the events with neutral-pion invariant mass in the vicinity of M π 0 will matter. Hence, we could safely proceed with the on-shell pion form factor and for the sake of simplicity, we will actually only consider the LO formula (15) in what follows. Combining Eqs. (11) and (15) with a π 0width-regulated propagator, we find the matrix element for the two-photon transition of the topology (2b): 2 Above, we used the shorthand notation ρσαβ k α 1 k β 2 = ρσ(k1)(k2) .
The two-photon-exchange topology is then the sum of the two sub-topologies from Eqs. (9) and (17): Above, the repeating process labels were suppressed. The two pairs of leptons are then coupled in the following way: Finally, we add up the one-and two-photon-exchange topologies: As was already mentioned, it is observed that the topology (2a) is rather suppressed compared to at least one of the other two, depending on the kinematical region, and this applies also to the (antisymmetric) interference term of topologies (2a) and (1) which vanishes upon symmetric integration; 3 the interference of (2a) and (2b) is exactly zero. This means that, in particular, we can eventually write with s = (P − r) 2 = (p 1 + p 2 + p 3 + p 4 ) 2 utilized also later on, and Moreover, the two amplitudes in Eq. (21) do not interfere either, so, to ≈ 10 % approximation of the branchingratio estimate, the square of the total LO matrix element is obtained as a direct sum of the matrix elements squared of topologies (1) and (2b). This simplifies the analytical structure significantly. Nevertheless, for the numerical results that follow we use the complete expression given by Eq. (20).

III. RESULTS
Having the squared matrix elements at hand, the Monte Carlo simulations can easily compare relative sizes of the respective contributions to the total branching ratio. As expected, the topology (2b) is by far the most dominant thanks to the pion-pole enhancement. Its contribution to the matrix element squared and hence the overall branching ratio can be rather simply evaluated since it does not lead to interference terms with the other two topologies. However, right outside the narrow pole region the topology (1) becomes prominent, so it is essential to retain it besides (2b); see Fig. 3 for the comparison of the individual contributions to the differential decay width.
Establishing the absolute branching ratio and its contributions without Monte Carlo methods requires, in the case of a 5-body decay, e.g. knowledge of the differential phase space in terms of at least 8 kinematical variables and subsequent multi-dimensional integration. This can be a pretty tedious job. At the same time, the spectra like those in Fig. 3 can be rather easily obtained by generating events flat-distributed in the momentum space and modulating them by the matrix element squared using, for instance, the accept/reject procedure or assigning weights given by |M| 2 to the events. After binning the (accepted or weighted) events -effectively performing the phase-space integral -and when properly normalized, this represents the differential decay width in the chosen kinematical variable (in terms of Fig. 3, this would correspond to dΓ( √ s) d √ s ). After the integral over this distribution is performed, one gets the branching ratio from the definition simply by normalizing the resulting width to Γ K + . What remains is to find the proper nor- 6.1(4) × 10 −12 6.0(6) × 10 −11 7.2(7) × 10 −11 Table I: Branching ratios excluding and including the pionpole region. The contributions of respective topologies are shown separately. They were obtained in terms of Eq. (23) (restricting the sum to the given subregions) and they correspond to the relevant areas under the curves in Fig. 3. The κ term serves as a handle on the model uncertainty; we set κ = 1. The total branching ratio for the topology (2b) is taken as a product of branching ratios B(K + → π + π 0 ) = 20.67(8) % and B(π 0 → 4e) = 3.38(16) × 10 −5 [18]. The last row is a sum of the first three rows, taking the fourth row as a model uncertainty and 6 % as the statistical uncertainty. The rightmost value in the bottom line excludes the (2b) topology.
malization so that, for instance, Fig. 3 is at scale, since the information about the phase-space volume is lost during the procedure described above; see Appendix A for further details. Inserting the (rescaled) phase-space volume Φ 5 from Eq. (A8) into Eq. (A10) leads to an intuitive formula for the branching ratio: Taking the phase-space-averaged matrix element squared we finally arrive at Taking into account the statistical uncertainty of 6 % and model uncertainty of 5 %, the last number considers 10 % relative uncertainty as a conservative estimate of the final number for the branching ratio when the topology (2b) is excluded. This gets modified when only a limited kinematical region is considered while including all the contributions. We find A more exhaustive list of branching ratios is presented in Table I. One can also check the overall consistency of the procedures used in the present work in the following way. The total branching ratio should be in general, one gets a factor 1 2π and two-body phase space accompanied by a suitable differential ds for every such extra branching. Above, the (dimensionless) two-body phase space φ 2 amounts to with the Källén triangle function λ(a, b, c) = a 2 + b 2 + c 2 − 2ab − 2ac − 2bc and the shorthand notation φ 2 (s) ≡ φ 2 (s, m 2 , m 2 ) = 1 8π 1 − 4m 2 /s. Further, we set s 12 = (p 1 + p 2 ) 2 , s 34 = (p 3 + p 4 ) 2 and s = (p 1 + p 2 + p 3 + p 4 ) 2 .
The volume of the histograms sorting the events based on specific kinematical variables is just the total number of events. If the histogram shape corresponding to N events in bins of width w b generated based on the flat phase space has a maximum n max ≡ n max (N, w b ) and is normalized so that it corresponds to (see Fig. 4 The integral over the described distribution then amounts to a sum over the weights |M K→π4e | 2 times an appropriate scale factor x taking into account the number of generated events and the phase-space volume (while being independent of w b ): Alternatively, this is equivalent to normalizing to the integral of the distribution (A5), so the scale factor becomes simply Hence, in general, the branching ratio B can be obtained as which represents a rescaled average of the matrix element squared over the phase space (with events randomly and evenly distributed in the momentum space) and is used to get the values in Table I.