On the role of leptonic CPV phases in cLFV observables

In extensions of the Standard Model by Majorana fermions, the presence of additional CP violating phases has been shown to play a crucial role in lepton number violating processes. In this work we show that (Dirac and Majorana) CP violating phases can also lead to important effects in charged lepton flavour violating (cLFV) transitions and decays, in some cases with a significant impact for the predicted rates of cLFV observables. We conduct a thorough exploration of these effects in several cLFV observables, and discuss the implications for future observation. We emphasise how the presence of leptonic CP violating phases might lead to modified cLFV rates, and to a possible loss of correlation between cLFV observables.


Introduction
Neutrino oscillations signal the presence of physics beyond the Standard Model (SM). Moreover, leptonic mixings and massive neutrinos open the door to the violation of charged lepton flavour. Present in several New Physics (NP) models, heavy neutral Majorana fermions (such as right-handed neutrinos) are among the most appealing extensions of the SM, often in association with a neutrino mass generation mechanism, as for example in the case of a fermionic seesaw [1]. Such sterile states (singlets under the SM gauge group) can also be present in more complete NP frameworks as SO (10) and its extensions. Depending on their mass scale and interactions with the SM fields, these new neutral fermions can offer a minimal framework to address the SM observational problems: account for neutrino oscillation data, offer an explanation to the baryon asymmetry of the Universe via leptogenesis, and possibly contribute to explain the dark matter problem.
Heavy neutral leptons (HNL), with masses ranging from the GeV to the tens of TeV, are among the most interesting minimal extensions of the SM, as they can be at the source of significant contributions to numerous observables, both at high-intensities and at colliders. Interestingly, the most minimal tree-level mechanism for neutrino mass generation -the type-I seesaw -calls upon the introduction of at least two such states to account for oscillation data. The type-I seesaw [1] and its low-scale variants, such as the Inverse Seesaw (ISS) [2][3][4], the Linear Seesaw (LSS) [5,6] and the ν-MSM [7][8][9], all call upon extending the SM via additional sterile fermions, allowing for Dirac and Majorana mass terms for the neutral lepton sector. Irrespective of the actual mechanism of neutrino mass generation under consideration, the mixings of the new states with the active left-handed neutrinos will lead to modifications in both leptonic charged and neutral currents, with a deep phenomenological impact. To study and numerically assess the impact of the heavy states, it is often convenient to consider simplified "ad-hoc" models, in which one adds n s sterile fermions to the SM field content. Such an approach allows to identify the most relevant effects and the consequences for the observables under scrutiny, and paves the way to a subsequent thorough study of complete models of neutrino mass generation via sterile fermions.
In addition to the masses of the new states and their mixings to the active neutrinos, constructions relying on heavy sterile states also open the door to new sources of CP violation: other than new Dirac CP violating (CPV) phases, should the massive states be of Majorana nature, further phases can be present. The role of these new CPV phases has been explored in analyses dedicated to CP violating observables, as is the case of electric dipole moments (EDM) of charged leptons [41][42][43][44]. New CPV phases have been also recently shown to play a crucial role in what concerns interference effects in LNV (and cLFV) semileptonic meson and tau decays [45], when more than one HNL is involved. Noticeably, while branching fractions of semileptonic meson or tau decays into same-sign and oppositesign dileptons are expected to be of the same order in the case of SM extensions by a single heavy Majorana fermion, this is no longer the case when the SM is extended by at least two HNLs, due to the possible interferences that might arise in the presence of multiple states. Depending on the CPV phases, one can have a modification (enhancement/suppression) of the rates of LNV modes and of the lepton number conserving ones. As a consequence, the non-observation of a given mode need not be interpreted in terms of reduced active-sterile couplings, but it could be instead understood in terms of interference effects due to the presence of several sterile states. (This effect is particularly amplified for processes with different charged leptons in the final state.) Likewise, an experimental signal of a lepton number conserving process and the non-observation of the corresponding lepton number violating one do not necessarily rule out that the mediators are Majorana fermions [45].
Similar studies have explored the role of a second heavy neutrino concerning the possibility of resonant CP violation [46], the effect of CP violation in high-scale seesaw scenarios in the context of renormalisation group running and leptogenesis [47,48], the impact for forward-backward asymmetries at an electron-positron collider [49], while others have compared the expected number of events associated with same-sign and opposite-sign dileptons at colliders in the framework of Left-Right symmetric models [49][50][51].
In this study we focus exclusively on the role of CPV phases (Dirac and Majorana) concerning an extensive array of cLFV observables. Although the importance of CPV phases in flavour observables is known, we provide here a dedicated thorough analysis, systematically studying the effects of CPV phases in cLFV observables. We work under the assumption of a unique source of lepton flavour violation, a generalised leptonic mixing matrix which now incorporates the active-sterile mixings. Following an analysis of the different observables and the implications for future observation, our results suggest that the presence of leptonic CP violating phases can strongly affect the predictions (either suppressing or enhancing the otherwise expected rates), and possibly lead to a loss of correlation between observables (typically present in simple SM extensions via heavy sterile fermions). As we will subsequently argue, the confrontation of unexpected cLFV patterns upon observation of certain channels could be suggestive of the presence of non-vanishing phases. As an example, one could have sizeable rates for µ → eee, and comparatively suppressed rates for µ → eγ. As we will also point out, another important point concerns accidental cancellations between distinct contributing terms to a given observable, as in the well-known case of µ − e conversion in nuclei [15]. The results of our work can be readily generalised for more complete NP models relying on the inclusion of heavy neutral fermions, provided that all complex degrees of freedom are systematically included (e.g. in Casas-Ibarra parametrisation [52]), and predictions for cLFV observables revisited.
This work is organised as follows: the cLFV observables addressed in the framework of the SM extension via sterile fermions are introduced in Section 2. In Section 3, a first insight into the role of Dirac and Majorana CPV phases is (analytically) presented, relying on a subset of representative cLFV observables. A phenomenological analysis of the impact of CPV phases on cLFV transitions is conducted in Section 4, including discussions on expected ranges and possible effects on the correlation patterns. In Section 5 we offer a comprehensive view of the model's parameter space, highlighting numerous points that could arise upon a future observation of several cLFV decays, followed by a final discussion and outlook in Section 6. Details on the model, expressions for the decay rates and other relevant amplitudes, as well as further information on the analysis, are collected in the Appendices.

cLFV observables: decays and transitions
In this section we present the cLFV processes in the framework of a "3 + n S " ad-hoc SM extension (see Appendix A for a detailed description), providing the expressions for the decay rates and other relevant amplitudes. As further discussed in Appendix A, the generalised leptonic mixings are encoded in a (3 + n S ) × (3 + n S ) unitary mixing matrix, U, of which the upper left 3 × 3 block corresponds to the left-handed leptonic mixing matrix, the would-be PMNS,Ũ PMNS , which is no longer unitary.
The cLFV observables addressed here receive contributions at the one-loop level arising from dipoles (γ-and Z-penguins) and/or from boxes. The vertex diagrams contributing to the cLFV decays are depicted in Fig. 1 and the box diagrams (adapted for the µ − e conversion rate for the sake of illustration) are presented in Fig. 2.
In the following subsections, we provide analytical expressions for the cLFV observables in the general 3 + n s framework (from Section 3 onwards the discussion will be focused on the case n s = 2).  In the quark internal lines, j = 1, ..., 3 runs over the quark families; in the neutral fermion ones, i = 1, ..., 3 + n S .

Leptonic cLFV decays
We first provide the expressions for the branching ratios of the "pure" leptonic cLFV decays, i.e. the radiative decays β → α γ and the 3-body decays 1 β → α γ γ . For the sake of brevity, we present here the expression for β → 3 α decays (the full expression for the most general case can be found in Ref. [16]). Notice that in our numerical analysis we have considered all possible decays for β → α γ γ .
In SM extensions via n s heavy sterile fermions, the rates for the radiative and three-body decays are then given by [15] Here, M W is the W boson mass, α w = g 2 w /4π denotes the weak coupling, s w the sine of the weak mixing angle, and m β (Γ β ) the mass (total width) of the decaying charged lepton of flavour β. The form factors present in the above equations are given by [15,16] G βα 1 In the most general 3-body cLFV decay, β → α γ γ the primed final state denotes the possibility of having equal or opposite charges for both γ and in which the sums run over the neutral mass eigenstates (i, j = 1, ..., 3 + n S ). The loop functions are given in Appendix B, with the corresponding arguments defined as x i = m 2 i /M 2 W . Furthermore, U is the (3 + n s ) × (3 + n s ) lepton mixing matrix (unitary) and C ij is defined as

Neutrinoless µ − e conversion in nuclei
The coherent conversion rate 2 in the presence of a nuclei (N) can be cast as [15] CR In the above expression, Γ capt. denotes the capture rate for the nucleus N; D, V (p) and V (n) correspond to nuclear form factors whose values are given in [55], with e being the unit electric charge. The form factors present in the above equations are given by [15,16] F µe In addition to the form factors previously defined in Eqs. (3 -6), F µe u,d call upon the following box contributions, in which x q = m 2 q /M 2 W and V is the Cabibbo-Kobayashi-Maskawa (CKM) quark mixing matrix. Other cLFV transitions and decays also occur in the presence of muonic atoms. This is the case of Muonium oscillations and decays (Mu − Mu, Mu → ee), and of the Coulomb enhanced decay µe → ee. The expressions for the cLFV rates of the latter can be found in Appendix B; here we only emphasise that all depend on combinations of µe form factors already introduced for the radiative and 3-body decays (corresponding to Z-and photon-penguins, and box diagrams).

cLFV Z decays
As the topology of cLFV Z decays is closely related to contributions at the origin of several cLFV leptonic decays and transitions (Z-penguins), we thus also include Z → α¯ β decays in our study 3 .
Since current searches do not distinguish the charges of the final state leptons, in our numerical analysis we consider the CP averaged decay rate, that is

cLFV searches: current status and future prospects
As extensively argued, the observation of one (or several) cLFV processes would be a clear signal of physics beyond the SM. Currently, there is a vast world-wide array of dedicated experiments and searches, at different energy scales, aiming at discovering such transitions. In Table 1 we list current experimental bounds and future sensitivities 4 for the observables considered in this work.

Phases do matter
In what follows we offer a first insight into the role of CP violating phases regarding a subset of (representative) observables. All other observables considered in the full phenomenological analysis of BR(Z → e ± τ ∓ ) < 5.2 × 10 −6 (OPAL [77]) O(10 −10 ) (FCC-ee [70] BR(Z → µ ± τ ∓ ) < 5.4 × 10 −6 (OPAL [77]) O(10 −10 ) (FCC-ee [70] Table 1: Current experimental bounds and future sensitivities on cLFV observables considered in this work. All limits are given at 90% C.L., and the Belle II sensitivities correspond to an integrated luminosity of 50 ab −1 . Section 4 can be understood from a generalisation of the discussion here carried. We work under the hypothesis of having n s = 2 massive sterile states; the neutral spectrum thus comprises 5 states, with masses m i (with i = 1, ..., 5), including the 3 light (mostly active) neutrinos and two heavier states, with masses m 4,5 . The leptonic mixings (whose precise origin is not specified) are parametrised by a 5 × 5 unitary mixing matrix, U. As described in Appendix A, the full mixing matrix U can be cast in terms 5 of 10 real mixing angles θ ij , 6 Dirac phases δ ij and 4 Majorana phases, ϕ i .
The study carried in this section is accompanied by a brief analytical discussion, relying on a simplified approach to this "3+2 toy model". In particular, and for the purpose of deriving clear (and compact) analytical expressions, in this section we work under the following assumptions: firstly, the mixings between the active and the sterile states (i.e., θ αi with α = e, µ, τ and i = 4, 5) are assumed to be sufficiently small so that cos θ αi ≈ 1 to a very good approximation. The 3 × 2 rectangular matrix encoding the active-sterile mixings can then be parametrised as with s αi = sin θ αi , and where δ αi (ϕ i ) denote Dirac (Majorana) phases. We further assume θ α4 ≈ θ α5 , and take the heavy states to be nearly mass-degenerate 6 , m 4 ≈ m 5 Λ EW (typically of the order of a few TeV).
We emphasise that all numerical results (displayed in the present section and throughout the full manuscript) are obtained without relying on any approximation, taking into account all the contributions present in the most general setup. On the left panel we present BR(µ → eγ) (blue), BR(µ → 3e) (orange) and BR(Z → eµ) (green); on the right one finds |G βα γ | (blue), |F βα Z | (orange) and |F β3α box | (green), choosing for illustrative purposes α = e and β = µ. In both panels, solid, dashed and dotted lines respectively correspond to the following heavy fermion masses: m 4 = m 5 = 1, 5, 10 TeV.
As this first discussion is dedicated to understanding and rendering visible the role of phases, no experimental constraints will be applied (certain observables might thus reach values already in disagreement with current experimental bounds).

cLFV decay rates: sensitivity to CPV phases
In what follows, we focus on µ − e sector flavour violation, and consider the following subset of observables: BR(µ → eγ), BR(µ → 3e) and BR(Z → eµ). We then devote a brief dedicated discussion to µ − e conversion in nuclei.
The role of Dirac phases In Fig. 3 we display the dependence of the above mentioned cLFV rates (and their form factors) on the Dirac phases. We set as an illustrative (benchmark) choice the following values for the mixing angles, θ 14 = θ 15 = 10 −3 , θ 24 = θ 25 = 0.01 and θ 34 = θ 35 = 0. Moreover, all phases are set to zero except the Dirac phase δ 14 . We also consider three representative values of the heavy fermion masses m 4 = m 5 = 1, 5, 10 TeV (associated with solid, dashed and dotted lines). As can be seen in the left panel, all considered observables have a clear dependence on δ 14 (the only non-vanishing phase considered), with the associated rates exhibiting a strong cancellation (typically amounting to around four orders of magnitude) for δ 14 = π, for all considered masses of the heavy sterile states. This behaviour can be understood by considering the pattern shown by the form factors contributing to cLFV radiative and 3-body muon decays, all displaying an (analogous) suppression for δ 14 = π.
Working in the limits above referred to, in Appendix C we present analytical expressions for the form factors contributing to the purely leptonic decays, including the full dependence on all phases. Regarding the dipole contributions, and in the case in which only δ 14 = 0, one has thus implying that in the simplest case of µ → eγ decays, the corresponding branching fraction for the radiative decays is given by with x 4,5 = m 2 4 /M 2 W = m 2 5 /M 2 W , thus indeed approximately vanishing for δ 14 = π. Similar results can be obtained for the photon penguin form factor F µe γ , as well as for one of the terms in the form factor Z , see Appendix C), all contributing to the rate of µ → 3e. For F (2) Z , after carrying the sum over ρ = e, µ, τ , one finds 7 also exhibiting a suppression for δ 14 = π. In the above equation we introduced G Z (x 4,5 ) ≡ G Z (x 4,5 , x 4,5 ), which we also use in the following for loop functions that depend on 2 parameters, in the limit of degenerate masses (cf. Appendix B). Finally, the remaining contributing term can be approximately given by again revealing the same behaviour 8 , which is also present in the box diagram contributions to µ → 3e. The latter lead to This brief discussion explains the behaviour of the different form factors, as depicted in the right panel of Fig. 3. Albeit carried in a limiting case (degenerate heavy states, identical mixings, etc.), this discussion is helpful in understanding the more complex behaviours which will emerge upon the numerical analysis presented in Section 4.
Concerning cLFV Z decays, we do not explicitly discuss their phase dependence here; let us just notice that the form factors parameterising Z decays (cf. Section 2.3) can be understood as Z-penguins at non-vanishing momentum transfer, i.e. q 2 = M 2 Z , and are therefore expected to have a very similar behaviour in what concerns the impact of the CPV phases. This can be quantitatively confirmed in Fig. 3 upon comparison of the dependence of the Z-penguin form factor F µe Z and the branching fraction of the Z → eµ decay on δ 14 .
For completeness, and before moving to the possible impact of the Majorana phases, notice that the individual form factors -and hence the full rates -depend on the mass of the new fermions (already manifest from the behaviour shown in Fig. 3), as shown in both panels of Fig. 4 (which were obtained considering vanishing values of all phases).

Impact of Majorana phases on cLFV decay rates
In what follows, we set the mixing angles to the same values as before, and choose the same three values of the heavy fermion masses. All phases are set to zero except ϕ 4 . The results for the observables (the same as studied concerning the Dirac phase) and the contributing form factors are displayed in Fig. 5, as a function of the Majorana phase ϕ 4 .
As expected, there is no dependence of the radiative decays on the Majorana phase (cf. the full expression for BR( β → α γ) given in Appendix C). This is also true for all dipole and dipole-like contributions. In contrast, the three-body decays (and the cLFV Z decays) do exhibit a non-negligible dependence on the Majorana phase, as can be verified from both panels of Fig. 5. This is especially true for heavier mass regimes, in which the relative contribution of the form factors sensitive to ϕ 4 (Z-penguins and to a lesser extent box-contributions) become more important (cf. right-handed panel 7 While one can in general neglect the contribution of the light (mostly active) neutrinos to the form factors here considered, that is not the case for F (2) Z , since the associated loop function GZ (x, y) does not vanish in the limit x ∼ 0, y 1. As can be seen in Appendix C, despite being more complex, the term corresponding to the "light-heavy" contribution exhibits a similar dependence on the Dirac phases; here we only consider the dominant "heavy-heavy" contribution. 8 Due to the contributions associated with the combination Cij (sum over all flavours), mixings involving the tau sector also contribute; taking them into account would lead to a similar dependence on cos δ14/2.  and G βα γ (blue), the Z-penguin form factor F βα Z (orange) and the box form factor F β3α box (green) to the total branching ratio of decays of the form β → 3 α (red), choosing for illustrative purposes α = e and β = µ.  Figure 5: Dependence of cLFV observables and several contributing form factors on the CP violating Majorana phase ϕ 4 (with all other phases set to zero). On the left panel we present BR(µ → eγ) (blue), BR(µ → 3e) (orange) and BR(Z → eµ) (green); on the right, one has |G µe γ | (blue), |F µe Z | (orange) and |F µ3e box | (green). In both panels, solid, dashed and dotted lines respectively correspond to m 4 = m 5 = 1, 5, 10 TeV.
of Fig. 4). Indeed, in the simplified limits of the form factors (see Appendix C), one verifies that only two contributions in the form factors depend on the Majorana phase, F Z and F (1) box . In the presence of a single non-vanishing Majorana phase, their expressions are: The impact of the Majorana phase on the cLFV Z decays can be also understood in analogy from the dependence of the corresponding Z penguin form factor. This is readily visible from inspection of  Figure 6: Dependence of cLFV observables (left) and several contributing form factors contributing to the decay rates (right) on the Majorana CPV phase ϕ 4 , for non-vanishing Dirac CPV phase, δ 14 = π. We again present on the left BR(µ → eγ) (blue), BR(µ → 3e) (orange) and BR(Z → eµ) (green), while on the right one finds |G µe γ | (blue), |F µe Z | (orange) and |F µ3e box | (green). We consider fixed values of the mass of heavy sterile states: solid, dashed and dotted lines respectively correspond to m 4 = m 5 = 1, 5, 10 TeV. and of a subset of form factors, on the Majorana phase ϕ 4 (similar to what was presented in Fig. 5), but now taking δ 14 = π.
The effects arising from the presence of both phases are clearly manifest, especially when compared with the plots of Figs. 3 and 5. Recall that δ 14 = π was found to lead to a strong cancellation of the form factors (see discussion of Fig. 3); thus, the non-vanishing contributions to the observables (3-body and Z decays) are associated with the form factors exhibiting a non-trivial dependence on ϕ 4 -as mentioned before. Relying again on the simple analytical estimates (see Appendix C), one now verifies that taking δ 14 = π leads to the following modification of the form factors sensitive to the Majorana phases (and hence of the observables): explaining the modified pattern (shift of π/2 and dependence on 2ϕ 4 ) visible in Fig. 6.

Neutrinoless µ − e conversion in nuclei and CP violating phases
We begin by illustrating how the predictions for the conversion rate (in particular the dependence on the heavy fermion mass) reflect the nature of the muonic atom, i.e. the chosen target nucleus (see Eq. (8)). This is shown on the left panel of Fig. 7, where we display CR(µ − e, N) as a function of the (degenerate) heavy masses, for Aluminium, Titanium, Gold and Lead nuclei, which have been chosen for past and future experimental searches. We fix the active-sterile mixing angles as before (θ 14 = θ 15 = 10 −3 , θ 24 = θ 25 = 0.01 and θ 34 = θ 35 = 0), and set all phases to zero. We recover the behaviour originally pointed out in [15], with the distinct rates all exhibiting a sharp cancellation for a given value of the heavy fermion mass -which we denote m c 4,5 (the actual value of m c 4,5 , and the "depth" observed in the conversion rate, both depend on the considered nucleus).
Until now, and for simplicity, we have not considered mixings of the sterile states to the third generation of leptons; however, and even though we are studying cLFV in the µ − e sector, certain observables are sensitive to tau mixings through the C ij coupling (sum over all flavours) which arises from the Z-penguin contribution, cf. Eq. (5). On the right panel of Fig. 7, we consider Aluminium nuclei, and again depict CR(µ−e, Al) vs. the heavy sterile masses, for different choices of θ 3j (j = 4, 5), keeping the other mixing angles and phases as before. As expected, and despite seemingly indirect, the impact is significant, both to the rate itself, and in what concerns the value of m c 4,5 associated with the cancellation in the conversion rate. This is all the most important when θ 3j θ 1j , θ 2j (j = 4, 5), as can be inferred from the green line.
The impact of the CPV phases on the conversion rate is studied in both panels of Fig. 8: for three choices of θ 3j = 0, 0.01 and 0.1, we consider the impact of Dirac and Majorana phases on CR(µ − e, Al), displayed as a function of the masses of the heavy states. For the case of vanishing ϕ 4 , the choice of Dirac phases leads to effects which could already be expected from the discussion of the previous subsection 9 . In particular, and as can be seen on the left panel of Fig. 8, notice the very important suppression for δ 14 = π; for the latter case, and for sizeable θ 3j one also observes a significant displacement of m c 4,5 , lighter by almost an order of magnitude. This is the result of a cancellation of (numerically) very small terms.
The different plots in Fig. 9 summarise (and generalise) the previous findings of Figs. 7 and 8: for fixed values of the "standard" input parameters (i.e. the mass of the heavy states and the active-sterile mixing angles), a variation of the phases -individually or as a joint effect -can lead to significant changes in the predictions for the conversion rate. This is seen in the panels of Fig. 9, with the rates ranging from as low as 10 −18 (dark blue), to values above 10 −14 (bright yellow), or even beyond 10 −13 (white). Here, one can also observe regions of constructive interference, which are of different origin. For example, in the upper left plot of Fig. 9 we show the conversion rate as a function of δ 14 and δ 24 ; if both phases are close to π, the suppression vanishes (as it depends on cos((δ 14 − δ 24 )/2)), while the complex exponential multiplying the contributions depends on (δ 14 + δ 24 )/2, thus leading to a "signflip" for values of δ 14 δ 24 π (cf. App. C). This affects the signs of the individual contributions to the different form factors, and can lead to an overall constructive interference, as visible in the plots. In the remaining panels in which one observes an enhancement of the rate with respect to vanishing  CPV phases (conversion rate as a function of δ 14 and δ 34 , and δ 34 and ϕ 4 respectively), the source of constructive interference solely lies in the Z-penguin form factor. In this case, the interference occurs between terms that depend on the tau-sterile mixing angles θ 3j , j = 4, 5 and the remaining terms, in conjunction with the effects of other phases.

Muonium anti-Muonium oscillations
Other cLFV observables rely on combinations of the already discussed form factors, so we will not address them individually. However, a few remarks concerning Muonium oscillations 10

cLFV and CP violating phases in the
Leptonic cLFV tau decays (i.e. τ → α γ or τ → 3 α with α = e, µ) receive contributions from form factors whose structure is analogous to that of the corresponding muon decays (allowing for the different flavour composition of the final state leptons); likewise, tau leptons can be present as final states of cLFV Z decays. Since the dependence of the observables on the CPV phases is in all similar to what has been discussed for the µ − e sector, we refrain from a dedicated analytical study (these observables will be included in the numerical study of Section 4). Notice that due to the large tau mass, one can also have semileptonic cLFV tau decays, with the final state composed of a light lepton and (light) mesons. However, we will not address them in the present study.

Other possible enhancements
So far, and relying on the simplifying approximation θ α4 ≈ θ α5 , we have mostly addressed effects of destructive interference leading to a strong suppression of the cLFV observables (due to a cosine dependence of the corresponding form factors on the CPV phases, see Eqs. (22)(23)(24)(25)(26)(27)(28)). However, in the most general case different behaviours (in particular generic enhancements) can be encountered upon relaxation of θ α4 ≈ θ α5 . For example, by considering a simple sign difference in one of the flavours (θ 14 = −θ 15 ), we are led to a generic cancellation as in this case, one finds a sinus-like dependence of the observables on the phases. For instance the photon dipole form factor is now given by to be compared with Eq. (22). Thus, non-vanishing phases now generically lead to enhancements, which can be quite important. This is illustrated in Fig. 11, where we show results analogous to those displayed in Figs

Phenomenological study: interference effects
Following the simple analysis of the previous section, which allowed a clear view of the effect of the CPV Dirac and Majorana phases, we now proceed to a more complete numerical study. In what follows, we will survey in a comprehensive way the simple "3+2 toy model", carrying a phenomenological analysis of a larger set of cLFV observables, now taking into account the available experimental constraints (see Appendix A.3). The latter include limits on the active-sterile mixings, results of direct and indirect searches for the heavy states, and EW precision tests, among many others. Finally, current bounds on searches for cLFV transitions are taken into account (cf. Table 1). In order to further explore the effects of the new non-vanishing CPV phases in a realistic way, we perform a random scan of the active-sterile mixing angles 11 for m 4 = m 5 = 1, 5 TeV. Leading to the 11 Here we relax the assumption of θα4 = ±θα5 by means of adding gaussian noise with a relative 1 σ deviation of 10 % numerical results displayed in this section, we first perform a scan with the phases set to zero and (randomly) select 2000 points consistent with all experimental data at the 3 σ level (corresponding to the blue points in the plots presented in this section). For each of the selected points we then randomly vary the phases δ 14 , δ 24 , δ 34 , and ϕ 4 in the interval (0, 2π), drawing 100 samples from a uniform distribution (shown in orange). Finally, we further add to the data set additional points which correspond to having systematically varied the phases on a grid for the "special" values {0, π 4 , π 2 , 3π 4 , π} (shown in green). This procedure allows to exhaustively study the impact of having non-vanishing CPV phases, especially regarding correlations between observables. Due to computational limitations, in this section we still do not take into account ϕ 5 nor δ i5 ; in the limit of degenerate masses, and nearly degenerate angles, non-vanishing values of the latter can be understood as leading to a phase shift (e.g. ∝ cos(ϕ 4 − ϕ 5 ), see Appendix C).

Correlation of µ − e µ − e µ − e observables
We begin by considering cLFV observables in the µ − e sector which receive contributions from unique topologies (dipole, Z-penguins and boxes) at one-loop level, and address the impact of CPV phases on the expected correlations. More specifically, on Fig. 12, we consider µ → eγ, Z → eµ decays and the probability of muonium-antimuonium oscillations, ∝ G M M , displaying the results for two heavy mass regimes, 1 and 5 TeV. Focusing first on the blue sets of points 12 (corresponding to vanishing CP violating phases), one confirms that there is a strong correlation between the three considered µ − e flavour violating observables, as expected. This is particularly manifest in the upper row of Fig. 12, since both BR(µ → eγ) and G M M do not depend on θ 3j ; on the other hand, non-vanishing values for θ 3j do contribute to Z → eµ decays (through the C ij term), hence one observes a spread of the points along the central straight line. This well-known behaviour (correlated predictions for a given value of the propagator's mass) has been explored in the literature as a means to test the underlying BSM construction, see e.g. [39,40]. Once CP violating phases are (randomly) taken into account, the correlation between the observables is strongly affected -if not lost, as can be seen by the significant spread of the corresponding orange points. This is especially visible for the plots in the right (m 4 = m 5 = 5 TeV). The situation becomes even more degraded once the "special" values of the phases are considered (those leading to vanishing values of certain contributions, as discussed in Section 3). As expected, the predictions for the observables are dramatically impacted, as visible from the green points. Finally, notice that the effect of phases can lead to either an increase or reduction of the cLFV rates with respect to the corresponding ones obtained in the vanishing phase limit. In some cases (like for BR(µ → eγ), with m 4,5 = 5 TeV) the new predictions may even be in conflict with current experimental bounds.
In Fig. 13 we consider the joint behaviour of general µ − e flavour violating observables, which now receive non-vanishing contributions from various types of diagrams (dipole, penguins and boxes). In particular, on the first row we present CR(µ−e, Al) vs. BR(µ → eγ), while on the second CR(µ−e, Al) vs. BR(µ → 3e), for m 4 = m 5 = 1 and 5 TeV (left and right columns, respectively). Despite the mixings between ν τ and the sterile states also leading to a spread in the case of vanishing phases (for µ − e conversion and 3-body decays), one still observes a visible correlation between the different sets of observables (see blue points) 13 . As already verified in Fig. 12, the presence of CP violating phases -especially for certain values of the latter -leads to a strong loss of correlation, as visible from the dispersion of the orange and green points (again more important for m 4 = m 5 = 5 TeV). On the third and final row of Fig. 13, we present the prospects for the correlation between µ → 3e decays and neutrinoless muon-electron conversion, but now for Gold (Au) nuclei. Notice that there are significant differences between Al and Au, which become particularly manifest for m 4 = m 5 = 5 TeV. In this to the samples, i.e. θα5 = ±θα4 ± 10 % . 12 In several plots one observes that in some extreme cases the predictions of certain observables already lie above the experimental limit; this merely reflects having taken points which are consistent with cLFV bounds at the 3 σ level, while the denoted experimental limits correspond to 90% C.L.. 13 The spread and visible cancellations for vanishing CP-violating phases are a consequence of accidental cancellations in µ − e conversion as discussed in Section 3.2 and of accidental cancellations due to opposite-sign mixing angles (e.g. θ14 ≈ −θ15) leading to a suppression of dipole operators as discussed in Section 3.5. case the correlation between the conversion rate in Gold and the three-body decays is more prominent than for Aluminium nuclei. This is a consequence of milder cancellations of the type discussed in Section 3.2: as seen from Figs. 7 and 8, in the limit of vanishing phases (and with θ 3j = 0), the values of the heavy propagator mass for which the accidental cancellation occurs (m c 4,5 ) are considerably lower for Gold than Aluminium, and effects of CP phases and/or non-vanishing θ 3j tend to further shift m c 4,5 to lower values.
Consequently, the CP violating phases might thus have an important impact regarding a future interpretation of data: let us consider a hypothetical scenario in which collider searches strongly hint for the presence of sterile states with masses close to 1 TeV. Should BR(µ → 3e)≈ 10 −15 be measured

Prospects for other observables
In view of the diversity of cLFV observables, and of the associated (future) experimental prospects, we have so far focused our phenomenological discussion on the most promising µ − e cLFV transitions and decays. Before concluding this section, we address the impact of the CP violating phases for the µ − τ sector (i.e. for τ → 3µ and Z → µτ decays), as well as for cLFV Muonium decays. In Fig. 14 we summarise the results of a study analogous to those presented in Figs. 12 and 13, displaying the correlated behaviour of high-and low-energy µ − τ sector cLFV observables for (non- ) vanishing Dirac and Majorana CPV phases, and for two values of the degenerate heavy neutral leptons' masses. Although the general prospects for observation are comparatively less promising, one nevertheless encounters the same phase-induced distortion of the correlation between observables, which was present in the limit of vanishing phases. Being also dominated by penguin transitions, the tau-lepton decay modes τ − → µ − e + e − and τ − → e − µ + µ − (i.e. only one flavour violating vertex) do not offer any additional insight with respect to the τ → 3µ and τ → 3e counterparts, and we find similar predictions for the associated rates. On the other hand, tau-lepton decays with an additional flavour violating coupling, that is τ − → e − µ + e − and τ − → µ − e + µ − , are transitions which are purely mediated by box diagrams. Thus, these are strongly suppressed when compared to other modes, with typically very small branching ratios, BR 10 −17 , thus clearly beyond any future sensitivity reach. In Fig. 15, we display for completeness 14 the prospects for Mu → ee decays, depicting its correlation with CR(µ − e, Al). In the most optimal scenarios, one can expect BR(Mu → ee) ∼ O(10 −22 ).

Overall view and further discussion
So far, we have analysed the implications of non-vanishing CPV phases in "idealised and simplified" scenarios, assuming certain relations between the model's parameters. We followed this approach to explore and maximise the effect of phases on the predictions for the cLFV observables. Clearly, broader scenarios must be considered, and in the present section we relax several of the previous assumptions. This aims at a more comprehensive overview, and allows a better confrontation with (hypothetical) future data.

Comprehensive overview of the parameter space
We conclude our study by a final comprehensive overview of this very simple SM extension. The numerical data presented in this subsection is now obtained as follows: firstly, and in what concerns the masses of the two heavy states 15 , we no longer take them to be degenerate, but rather assume their masses to be sufficiently close to allow for interference effects 16 ; in practice, and for fixed m 4 , random values of m 5 are obtained from half-normal distributions with the scale set to a value representative of the width of the sterile states (in this case ∼ 50 GeV). Concerning the active-sterile mixing angles, these are now independently varied: more specifically, we draw samples from log-uniform distributions, further randomly varying their signs. For m 4 = 1 TeV, the ranges of the parameters to be here explored are then It is worth noticing that these ranges lead to scenarios complying with experimental bounds (see Section 2.4 and Appendix A.3). In our analysis, we thus (randomly) select 10 4 points consistent with all experimental data. For each tuple of mixing angles we then vary all CPV phases associated with the sterile states, i.e. δ α4,5 , ϕ 4,5 ∈ [0, 2π], drawing 100 values for each of the four from a uniform distribution. The upper limits on the intervals for the mixing angles are inferred from requiring agreement with the most constraining current bounds; clearly no lower limit for the mixing angles is phenomenologically relevant. However, we have limited ourselves to regimes that do not lead to cLFV predictions excessively far away from the corresponding future experimental sensitivity. Thus, it is important to stress that the resulting predictions for the observables could in principle be extended to extremely tiny values of the rates, should we have explored all the allowed ranges for the mixing angles. In summary, no conclusions concerning lower limits for the cLFV observables should be drawn from this analysis. The outcome of this comprehensive analysis is shown in Fig.16, where we display the predictions for several cLFV rates (focusing on the µ−e and τ −µ sectors), in particular in what concerns correlations between same-sector observables. As before, we present the predictions obtained in the case in which all CPV phases are set to zero (blue points), and then those corresponding to a random scan over all Dirac and Majorana phases (orange points). For a fixed scale of the heavy propagators (with the latter sufficiently close in mass), and as would be expected, one finds correlated same-sector observables. In addition to enlarging the range of the predictions for the different observables (possibly leading to a conflict with current bounds for certain µ − e observables), non-vanishing CP phases lead to a visible spread of the predicted rates. Should "special" values of the phases be imposed on the scan, the cancellation effects (and associated decrease of the rates) would have been even more striking.
Throughout the discussion we have focused on the µ − e and τ − µ sectors, in view of the most promising prospects (both theoretical and experimental). In general, the associated predictions for τ − e cLFV transitions typically lie beyond future experimental sensitivity. However, and to fully explore this very minimal SM extension, one would require a comprehensive probing of the associated cLFV predictions in all 3 flavour sectors (i.e. µ − e, τ − µ and τ − e transitions).  Table 1.
Profiting from the data collected leading to the results displayed in Fig. 16, we have tried to infer which would be the required future sensitivity for the τ − e channels so that the regimes (mixing angles and CP phases) leading to predictions for µ → eγ, µ → 3e, µ − e conversion in Al, τ → 3µ and Z → µτ , all within future experimental sensitivities, would also be within reach of τ → eγ and τ → 3e dedicated searches. Requiring that at least 68% of the previously mentioned subset be within τ − e future reach would imply the following ideal experimental sensitivities 17 : In other words, should a signal of cLFV in µ − e and τ − µ transitions be observed at the current and near-future facilities, an improvement of circa 4 orders of magnitude in the τ − e sensitivity is needed in order to obtain competitive constraints on these SM extensions via heavy neutral leptons from all flavour sectors.

Reconciling cLFV predictions with future observations
As discussed extensively in the previous (sub)sections, CPV phases can impact the predictions for the cLFV observables, enhancing or suppressing the distinct rates. To conclude the discussion, we have identified a small set of representative (benchmark) points, which reflect not only the effect on the rates, but also the impact that taking into account the CPV phases might have on the interpretation of experimental data (or negative search results).
In Table 2 we present the predictions for several cLFV observables for three configurations of active-sterile mixing angles, in the case of vanishing CPV phases (P i ), and for non-vanishing values of the phases (P i ): The variants P i have identical mixing angles, but in association with the following phase configurations: We have chosen m 4 = m 5 = 5 TeV for all three benchmark points. Table 2: Predictions for several cLFV observables in association with three configurations with vanishing CPV phases, P i (i = 1 − 3) and associated variants with non-vanishing CP violating phases, P i , see Eqs. (32,33). We have taken m 4 = m 5 = 5 TeV. The symbols ( , , •) respectively denote rates already in conflict with current experimental bounds, predictions within future sensitivity and those beyond future experimental reach.
The first point (P 1 ) represents a case for which only two cLFV observables would be within future experimental reach, µ → 3e and µ − e conversion in Aluminium; however, in the presence of CP phases (P 1 ), the predictions for the different considered observables are now all within future sensitivity. The points P 2 and P 2 correspond to a similar scenario, but for which only the two considered µ − τ observables lie within future reach in the case of vanishing phases. The third and final point (P 3 ) clearly illustrates the importance of taking into account the possibility of CP violating phases upon interpretation of experimental data. Negative search results for the different µ − e flavour violating transitions would lead to the exclusion of the associated mixing angles (for heavy masses ∼ 5 TeV); however, and should CPV phases be present, the considered active-sterile mixing regime can be readily reconciled with current bounds 18 (with µ → eγ now even lying beyond experimental reach). A similar exercise could be carried for other heavy mass regimes, leading to analogous conclusions. This demonstrates the crucial role of CPV phases in evaluating the viability of a given scenario in what regards conflict/agreement with the associated cLFV bounds.

Conclusions
In this work we have thoroughly addressed the impact of leptonic CP violating phases on the predictions of the rates of several cLFV observables, focusing on minimal SM extensions by heavy Majorana sterile fermions. Despite their minimality and simplicity, these extensions can be interpreted as representative of more complete constructions calling upon the addition of heavy neutral fermions (as is the case of several low-scale seesaw realisations). Via their mixings with the light (mostly active) states, and as a consequence of a departure from unitarity of the would-be PMNS mixing matrix, the new states open the door to contributions to numerous observables.
Here we have considered a simple case with 2 heavy neutral fermions, taking them close in mass in order to explore the potential impact of the new CPV phases on cLFV observables. These states could very well be embedded in a seesaw, and the latter even incorporated in more complete BSM frameworks. The conclusions drawn in this work are thus always valid once one considers that the source of lepton flavour violation stems from the enlarged leptonic mixing.
Building upon an analytical insight, our numerical study reveals that the CP violating phases can indeed lead to important effects in cLFV transitions and decays, with an impact for the rates (enhancement or suppression of the predictions obtained for vanishing phases). Moreover, whenever correlations between observables would be typically expected (in association with the dominance of a given topology for certain regimes of the model), one also encounters a potential loss of correlation. Furthermore, our analysis suggests that the non-observation of a given observable (usually expected to be within experimental reach in view of the measurement of another one) should not be a conclusive reason to disfavour a given regime.
The conclusions of this work can be generalised for a given BSM construction, provided that all complex degrees of freedom are consistently taken into account, and predictions for cLFV observables re-evaluated in view of the potential presence of new CP violating phases.
In the near future, should potentially new (unexpected) cLFV patterns emerge upon observation of certain processes, this could be interpreted as possibly hinting towards the presence of non-vanishing CP violating phases (under the working hypothesis of SM extensions via heavy neutral fermions).

A A minimal "3+2 toy model" for cLFV studies and constraints
Here we collect the most relevant information concerning this minimal extension of the SM. As mentioned before, we will consider an effective toy model, in which 2 sterile (Majorana) states are added to the SM, a so-called "3+2 toy model". The neutral spectrum comprises 5 states, with masses m i (with i = 1, ..., 5): this includes the 3 light (mostly active) neutrinos and two heavier states (with masses m 4,5 ).

A.1 Modified interaction Lagrangian
The leptonic mixings are parametrised by a 5 × 5 unitary mixing matrix, U; its upper left 3 × 3 block corresponds to the left-handed leptonic mixing matrix, the would-be PMNS,Ũ PMNS . The ensuing non-unitarity ofŨ PMNS will lead to modified charged and neutral currents, which can be cast in the physical basis as U αj¯ α γ µ P L ν j + H.c. , , with n S = 2, and in which we recall that In the above, the indices α, ρ = 1, . . . , 3 denote the flavour of the charged leptons, while i, j = 1, . . . , 3 + n S correspond to the physical (massive) neutrino states; P L,R = (1 ∓ γ 5 )/2, g w denotes the weak coupling constant, and cos 2 θ w = M 2 W /M 2 Z . The coefficients C V and C A parametrise the SM vector and axial-vector currents for the interaction of neutrinos with charged leptons, respectively given by C V = 1 2 + 2 sin 2 θ w and C A = 1 2 .

A.2 Parametrisation of the leptonic mixings
Following for instance [41], U can be parametrised through five subsequent rotations R ij (with i = j), and a diagonal matrix including the four physical Majorana phases, ϕ i The above rotations are of the form (illustrated by R 45 ): As already noticed (and clear from the L W ± term in Eq. (34)), the mixing in charged current interactions is parametrised via a rectangular 3 × (3 + n S ) (i.e. 3 × 5) mixing matrix, of which the 3 × 3 sub-block encodes the mixing between the left-handed leptons,Ũ PMNS . The deviations ofŨ PMNS from unitarity [2,3] can be conveniently parametrised in terms of a matrix η [85] In our numerical analysis, we fix the mixing parameters of the mostly active neutrinos to the central values of the NuFIT 5.0 results [86] (without atmospheric data), which we quote for convenience in Table 3.
The predictions for cLFV transitions are mostly independent of the ordering of the spectrum and of the lightest neutrino mass. For concreteness (and especially in what regards constraints from neutrinoless double beta decay) throughout our numerical analysis we assume the light spectrum to follow a normal ordering, and vary the lightest neutrino mass in the range m 0 ∈ [10 −5 , 10 −3 ] eV.  Table 3: Global fit results obtained by NuFIT 5.0 [86] for neutrino mixing data, not including experimental data from oscillation experiments measuring atmospheric neutrinos. In our numerical analysis we assume normal ordering of the light neutrino spectrum and fix the neutrino mixing parameters to their central values.

A.3 Constraints on SM extensions via sterile fermions
In addition to the cLFV observables, we constrain the mixing angles θ α4,5 using experimental data on various precision observables that do not depend on the new CPV phases. These can be divided into low-and high-energy observables. At high energies, one has lepton flavour universality ratios of W decays (see e.g. [87]) with 1,2 ∈ {e, µ, τ }. Of paramount importance is also the invisible decay width of the Z boson, Γ(Z → inv), which is typically reduced in scenarios with heavy sterile states that have sizeable mixings with the mostly active neutrinos. At low energy, we take into account (semi-)leptonic decays of the τ -lepton and leptonic decays of light mesons, out of which various ratios can be constructed; these are sensitive to a modified W ν vertex, and therefore indirectly sensitive to effects of the sterile states. These are (see for instance [87,88]) where R P ≡ Γ(P + → eν) Γ(P + → µν) (P = π, K) (42) Finally, we take into account the upper bounds on the entries of η as derived in Ref. [89], indirectly taking into account further constraints from modifications of the Fermi constant G F , the weak mixing angle sin 2 θ w , the mass of the W -boson among others.
For HNL with masses around the TeV scale (as these here considered), constraints from direct searches at colliders and limits from cosmological considerations such as big bang nucleosynthesis (BBN), are typically not competitive and therefore not taken into account. However, another important constraint of theoretical nature can be derived by imposing that decays of the HNL comply with perturbative unitarity [90][91][92][93][94][95], which gives a direct bound on their decay width . Since the dominant contribution arises from a W -boson exchange, one can obtain a bound on the sterile masses and their couplings to active states which can be written as Additional sterile states are known to lead to modifications of the predictions for the effective mass to which the amplitude of neutrinoless double beta decay is proportional to, m ee . In the presence of n S heavy states, the contributions to m ee can be written as [96,97] m ee in which p 2 corresponds to the virtual momentum, with p 2 −(100 MeV) 2 . In our analysis, we take into account the KamLAND-ZEN [98] upper limit on m eff ee (61 ÷ 165) MeV which is obtained using the isotope 136 Xe.

B Additional information on cLFV transitions and decays B.1 Loop functions
Below we summarise the loop functions for the cLFV transitions discussed in our study, as well as their relevant limits 19 , as taken from Refs. [15,16]. The photon dipole and anapole functions and relevant asymptotic limits, are given by The loop functions of the Z-penguins are given by a two-point function and two three-point functions which are symmetric under interchange of the arguments.
The (symmetric) box-loop-functions and their limits are given by
Muonium oscillations and decay Muonium (Mu) is a hydrogen-like Coulomb bound state formed by an anti-muon which slows down inside matter and captures an electron, e − µ + [99]. The bound state is free of hadronic interactions, and its electromagnetic binding is well described by the SM electroweak interactions. The spontaneous conversion of a Muonium atom to its anti-atom (Mu = e + µ − ) [100], as well as its decay to a pair of electrons, have been identified as promising cLFV observables.
In the presence of (V − A) × (V − A) interactions, Muonium anti-Muonium oscillations can be described by the following effective four-fermion interaction, in terms of an effective coupling G M M , In extensions of the SM with sterile neutrinos Mu-Mu conversion occurs at the loop level, being exclusively mediated by a set of 4 independent box diagrams; while a first set is common to both Dirac and Majorana neutrinos a second one is only present if neutrinos are Majorana particles, effectively amounting to two Majorana mass insertions (see Refs. [78,79]). Working in the unitary gauge, the computation of the different box diagrams allows to write the effective coupling of Eq. (54) as [78,79]: where F Xbox (x i , x j ) and G box (x i , x j ) are the relevant loop functions 20 given in Appendix B.1, with , i = 1, ..., 3 + n s (further details can be found in [22]).
In the presence of new physics, Muonium can also undergo the cLFV decay Mu → e + e − . In the SM extended by n s heavy neutral leptons, the cLFV Muonium decay rate is given by [101] in which Γ µ = G 2 F m 5 µ /(192π 3 ) denotes the muon decay width and α em = e 2 /4π the electromagnetic coupling constant, with |M tot | the full amplitude (summed (averaged) over final (initial) spins) [101], where the corresponding form factors have been previously introduced (see Eqs. (3 -6), setting β = µ and α = e).

C Phase dependence of cLFV observables -analytical expressions
In order to analytically study the phase dependence of the cLFV form factors, we work under the approximation sin θ α4 ≈ sin θ α5 1. Let us further assume that the masses of the heavy states are close to each other and of the order of a few TeV, m 4 ≈ m 5 Λ EW .

C.1 Photon penguins
The photon penguin form factors exhibit a similar structure. From inspection of the relevant loop functions (cf. Appendix B.1), it can be seen that contributions from the light neutrino mass eigenstates are negligible. Furthermore, since the heavy states are assumed to be close in mass, the loop functions are approximately equal and can be thus factored out of the sum. The form factor is then given by and inserting now the entries of the mixing matrix of Eq. (36) (in the limit in which ϕ i = 0) yields  20 We note a sign difference between the function F box in Ref. [79] and the function F Xbox in our convention.
with s ij = sin θ ij and x 4,5 = m 2 4 /M 2 W = m 2 5 /M 2 W , and where we have defined with the properties ∆ αβ i = −∆ βα i , ∆ αρ i + ∆ ρβ i = ∆ αβ i and ∆ αα i = 0. Further assuming that the mixing between the active neutrinos and ν 4 and ν 5 is approximately equal, i.e. sin θ α4 ≈ sin θ α5 , we can simplify G βα γ to G βα γ ≈ s α4 s β4 e − i 2 (∆ αβ 4 +∆ αβ such that the branching fraction for the radiative decays is given by Similar results can be obtained for F βα γ .

C.2 Z penguins
The form factor generated by Z penguin diagrams can be split into three different parts. It consists of bubble diagrams and triangle diagrams with or without Majorana mass insertions, and is given in the "3 + 2 toy model" by The first term which is proportional to the function F Z can be rewritten in the same way as the photon penguins, and is thus given by in analogy to G βα γ , see Eq. (59). The second and the third terms are more involved due the presence of C ij . The second term can be written as where we introduced G Z = G Z (x 4,5 ) ≡ G Z (x 4,5 , x 4,5 ), which is also used in the following for loop functions that depend on 2 parameters, in the limit of degenerate masses (cf. App. B.1). Due to the appearance of C * ij , in the third term of Eq. (64), the Majorana phases will also be present. This corresponds to a Majorana mass insertion in the corresponding triangle diagram. The last term can be cast as In the case of vanishing Dirac phases this can be further simplified to Moreover, and unique to the Z-penguin, there are non-negligible contributions to the form factor stemming from light and heavy virtual neutrinos in the loop. We write the corresponding limit of the loop function as G Z (0, x 4,5 ) = G Z . As usual, using the above approximations, the corresponding part of the form factor can be cast as Making use of the unitarity of the full 5 × 5 mixing matrix, i.e.
one finally has F (2)

C.3 Box diagrams
The form factor generated by box diagrams is given by where the first term corresponds to a diagram with a possible Majorana mass insertion, thus depending on the Majorana phases. The first term F (1) box can then be written as F (1) box ≈ (U 2 α4 + U 2 α5 )(U * β4 U * α4 + U * β5 U * α5 ) G box (x 4,5 ) = s 2 α4 e −2i(δ α4 −ϕ 4 ) + s 2 α5 e −2i(δ α5 −ϕ 5 ) s α4 s β4 e i(δ α4 +δ β4 −2ϕ 4 ) + s α5 s β5 e i(δ α5 +δ β5 −2ϕ 5 which again can be further simplified in the case of vanishing Dirac phases to The second term is independent of the Majorana phases and can be written as The box diagrams contributing to neutrinoless muon-electron conversion show a similar behaviour as that of photon-and Z-penguin diagrams with one neutrino in the loop.