Purely flavor-changing Z' bosons and where they might hide

A plethora of ultraviolet completions of the Standard Model have extra U(1) gauge symmetries. In general, the associated massive $Z^\prime$ gauge boson can mediate flavor-changing neutral current processes at tree level. We consider a situation where the $Z^\prime$ boson couples solely via flavor-changing interactions to quarks and leptons. In this scenario the model parameter space is, in general, quite well constrained by existing flavor bounds. However, we argue that cancellation effects shelter islands in parameter space from strong flavor constraints and that these can be probed by multipurpose collider experiments like ATLAS or CMS as well as LHCb in upcoming runs at the LHC. In still allowed regions of parameter space these scenarios may help to explain the current tension between theory and experiment of $(g-2)_\mu$ as well as a small anomaly in $\tau$ decays.


Introduction
In searching for new physics it is prudent to explore the limits of applicability of standard tests and probe for corners in parameter space where they can be evaded. We can then turn around and check if in these regions other tests become more powerful. It is in this spirit that in this paper we want to examine flavor-changing Z bosons coupling to quarks and leptons. In this case severe constraints arise from precision tests of flavor-changing neutral currents (FCNCs), in particular on mesons. Yet we will see that there are still interesting areas of parameter space that can be probed with direct production at the LHC and "non-flavored" measurements such as (g − 2). Flavor-changing Z bosons could be a remnant of a solution to the still unsolved question of the origin of the flavor structure of the Standard Model (SM). Indeed, one of the earliest approaches towards an explanation of Yukawa coupling patterns and the family structure of the SM fermions was the introduction of so-called horizontal or gauged flavor symmetries [1][2][3][4][5][6]. For example, the different copies of up-and down-type quarks, charged leptons and neutrinos can transform as a multiplet under a new horizontal SU(2) symmetry group. Likewise, one can assign charges under a new local U(1) gauge symmetry. The breaking of such symmetries generally leads to the emergence of new massive gauge bosons mediating FCNCs. In such scenarios, care has to be taken in that the magnitude of such an effect does not violate experimental constraints [7]. Nevertheless, gauged flavor models are enjoying a renewed popularity [8,9]. Especially, the case of a new U(1) gauge symmetry has been studied extensively in the past (see Refs. [10,11] for reviews).
In this paper we take a phenomenological approach of extending the SM by a neutral massive Z boson, which is possibly the remnant of a broken gauge symmetry, with the simplest possibility being a U(1) 1 . Specifically, we consider models with exclusively flavor-changing couplings, one in the quark and one in the lepton sector Purely flavor-changing interactions provide a simple but interesting test case. On the one hand they provide a maximally flavor-changing effect. On the other hand they are often more difficult to detect. For example, if the quark part of the interaction involves a b-and an s-quark, production at proton colliders like the LHC requires reliance on the sea-quarks in the protons which are less abundant 2 . Similarly at LEP simple s-channel production of Z -bosons via the lepton couplings is not possible as the initial state is not flavored.
The paper is structured as follows. In Section 2 we will discuss collider constraints on our model from reinterpreting an ATLAS search for neutral resonances in eµ, eτ and µτ final states [12]. In Section 3 we review relevant existing constraints on our model. In this context we will discuss meson mixing, meson decays into charged leptons and neutrinos, lepton decays, muonium-antimuonium oscillations, LEP searches as well as electron and muon (g − 2) measurements. In view of the (g − 2) µ anomaly, we will also discuss a possible explanation of the observed shift ∆a µ within our model together with a small anomaly in τ -decays (cf. also [13,14]). Finally, we will interpret and wrap up our results in Section 4. The summary plots of our findings of collider and existing flavor constraints on our model can be found in Section 2 in Figs. 1 to 3 and in the Appendix C in Figs. 20 to 22. An example interpretation of the (g − 2) µ and τ -decay anomalies is depicted in Fig. 18. We focus on the situation where the Z bosons are heavy M Z M Z . Unless otherwise stated we take the lepton sector couplings to be purely right-handed. Most plots, however, are also applicable to the purely left-handed case. The additional limits present in this situation are given as dotted and dash-dotted lines.

Flavor violation from a collider point of view
One main goal of this paper is to reinterpret existing neutral resonance searches at the LHC in the context of a flavor-changing Z boson. This provides us with new constraints on the induced FCNCs, complementary to the usual bounds coming from flavor and electroweak precision experiments (see Section 3). Previously, Davidson et al. [15] have investigated flavor-changing four-quark contact interactions coming from new physics from a scale M M W in an effective field theory (EFT) approach. They consider the various four-quark operators the LHC is sensitive to, also including quark flavor violating operators, which are of the type with X, Y ∈ {L, R} and the indices (m, i, n, j) denoting flavor. Then they derive a limit on their suppression scale Λ by reinterpreting existing LHC dijet analyses. We show the corresponding limits as brown areas in the figures. In our model, however, we are considering combined lepton and quark flavor violation. While our model also contains the effective operators (2.1) we have additional operators of the type This type of operator can be generated from a Z exchange in the full theory and consequently a bound on it can be turned into a constraint on the corresponding Z couplings. In the following, we therefore want to reinterpret an existing ATLAS analysis of heavy neutral particles decaying to eµ, eτ or µτ [12] in the light of our flavor-violating Z model.

Reinterpreting collider searches
The model we consider in this paper induces ∆F = 2 flavor-violating processes of the type qq → Z → . In order to constrain the relevant couplings g qq and g we first need an expression for the corresponding cross section within our model. Introducing the non-chiral reduced couplingḡ we can derive an approximate expression for the cross section scaling as This expression gives a valid estimate for the cross section at parton level. However, we cannot access this cross section at the LHC directly. As we are dealing with a hadron collider we have to take into account parton distribution function effects and hadronization. Moreover, the observable cross section will also be affected by a number of detector effects like finite resolution, mistags, acceptance etc. Our approach to incorporate all these effects is quite straight forward. We simulate the total cross section σ MC for the process pp → Z → in our model for a given combination {qq , } of flavor-violating interactions. The ratio of the simulated cross section σ MC to the ATLAS limit on the cross section σ lim allows us to derive an approximate limit on the off-diagonal quark couplingḡ qq as a function of the lepton couplingḡ according to Hereḡ qq , MC andḡ , MC denote the values of the reduced couplings used for the simulation and denotes the ratio of them. We immediately notice that Eq. (2.5) has a pole at This simply indicates that for values ofḡ close to this pole the lepton coupling is so weak that even for a very large value of the quark couplingḡ qq the signal cannot be distinguished from background, i.e. the process is unobservable at the LHC. In other words, given the observed cross section limit σ lim , we are not able anymore to set a limit on the quark couplingḡ qq for such low values ofḡ . We still have to determine the simulated cross section σ MC . First, we calculate the leading order (LO) cross section σ LO with madgraph v2.3.3 [16] for the processes under consideration. Next, we determine a mass-dependent K-factor to take into account next-to-next-to leading order (NNLO) effects. In the auxiliary material of Ref. [12] the NNLO cross sections of the Sequential Standard Model (SSM) Z are provided. In order to determine the values of the K-factor, we calculate the LO cross section for the SSM Z in pythia v8.215 [17] and compare to the provided NNLO results. The K-factor is then found from the ratio The values we have determined in our analysis are given in Appendix B in Table 1. Finally, we determine an effective mass-dependent acceptance times efficiency A × from the ratio of the number of events that survived the detector plus analysis cuts to the expected total number of events at NNLO The number of events at NNLO N NNLO = σ NNLO × dt L is simply obtained from multiplying the cross section by the integrated luminosity. The numerical values of the determined acceptance times efficiency we used in our calculations can be found in Appendix B in Table 2.
Putting everything together we obtain the full simulated cross section as (2.10)

Numerical evaluation
In this section, we present as an example the exclusion limits on the couplings g R qq and g R for the combinations {qq , } = {g bs , g eµ }, {g bs , g µτ }, {g bd , g eµ }. First, we define the ratios of leftto right-handed couplings Lepton decays Muonium LEP II Figure 1: The flavor-violating couplings g R bs and g R eµ for a Z boson of M Z = 750 GeV and a coupling ratio of ρ q = g L bs /g R bs = 0.1324 (for an explanation of ∆ρ see Section 3.1.1). The red area indicates the limit from the ATLAS analysis of the process pp → eµ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The green area is the limit coming from the meson decay B 0 s → eµ. The gold and magenta areas are purely leptonic limits coming from LEP and muonium oscillation constraints. Finally, the black dotted line is the exclusion limit from B + → K +ν ν and the black dash-dotted line the limit from lepton decays. These last two limits, however, apply only if we consider g L eµ instead of g R eµ .
We present the derived bounds on the quark limits always as a limit in terms of the right-handed coupling g R and the corresponding ratio ρ q to the left-handed coupling. Hence, the limit on the quark couplings derived in the last section reads (2.12) Examples of the limit derived from reinterpreting the ATLAS analysis [12]   The flavor-violating couplings g R bs and g R µτ for a Z boson of M Z = 750 GeV and a coupling ratio of ρ q = g L bs /g R bs = 0.1324 (for an explanation of ∆ρ see Section 3.1.1). The red area indicates the limit from the ATLAS analysis of the process pp → µτ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The green area represents the excluded region from the decay B 0 s → µτ . The black dotted line is the limit from B + → K +ν ν and the black dash-dotted line the limit from lepton decays. These last two limits, however, apply only for left-handed lepton couplings. The light and dark cyan areas depict the 4 and 5 σ exclusion bands from ∆a µ . Section 3. We show plots for these three particular combinations of couplings as they illustrate all the main features and relevant limits. First, LHC limits are generally strongest for a nonzero coupling in the eµ sector and weakest for a coupling in the µτ sector. The limits from meson decay into charged leptons, too, are in general stronger in the eµ sector than in the eτ and µτ sector. In addition, the eµ sector has further strong leptonic constraints, namely the one from LEP (which is also present in the eτ sector) and the one from muonium oscillations. In contrast to all other quark combinations the bs sector is not tested by the dijet limits of Ref. [15]. This can be seen by comparing Figs. 1 and 2 with Fig. 3. The brown dijet region is present in the latter, but not in the former two. Lastly, only the µτ sector receives sizable constraints from the (muon) anomalous magnetic moment as shown in Fig. 2. In the plots we have indicated also Lepton decays Muonium LEP II Figure 3: The flavor-violating couplings g R bd and g R eµ for a Z boson of M Z = 750 GeV and a coupling ratio of ρ q = g L bd /g R bd = 0.1335 (for an explanation of ∆ρ see Section 3.1.1). The red area indicates the limit from the ATLAS analysis of the process pp → eµ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The green area is the limit coming from the meson decay B 0 → eµ. The gold and magenta areas are purely leptonic limits coming from LEP and muonium oscillation constraints. The black dotted line is the exclusion limit from B 0 → π 0ν ν and the black dash-dotted line the limit from lepton decays. These last two limits, however, apply only for left-handed lepton couplings. The brown area depicts the limit from the four-quark contact interaction [15].
limits from decays of mesons and leptons into neutrinos (black dotted and dash-dotted lines). They are only applicable if we replace the right-handed lepton coupling by a left-handed one. Furthermore, we note that the limits from meson decays are absent in the cu sector (see plots in Appendix C). The limits shown in Figs. 1 to 3 were derived for a Z mediator of M Z = 750 GeV. In general, we have derived these limits for various masses in the range of 200 GeV to 3 TeV for the flavor combinations qq ∈ {sd, bs, bd, cu} and ∈ {eµ, eτ, µτ }. In Appendix C an example exclusion plot for each flavor combination is given for a Z boson with M Z = 1 TeV.  Figure 4: (Left) Observed exclusion limits on the branching fraction times cross section in the eµ channel [12]. The limits have been extrapolated down to masses of 200 GeV by applying a constant continuation. (Right) The blue curve shows the exclusion limits on the branching fraction times cross section at √ s = 8 TeV given in Ref. [12]. The black curve shows the preliminary limit at √ s = 13 TeV given in Ref. [18]. The green curve depicts the 8 TeV limit rescaled to the 13 TeV luminosity. The order of magnitude agreement between this curve and the 13 TeV limit in this region indicate that our luminosity rescaling is sensible.
supplementary material to this paper the full set of exclusion plots can be found. From these plots one can see that with increasing mass all limits weaken. Except for the case of direct production at the LHC, limits typically scale as g ∼ M Z .
We want to point out, however, that in order to obtain limits for Z bosons with masses below 500 GeV we have extrapolated the ATLAS limits rather optimistically. As shown in the left panel of Fig. 4, we assume a constant scaling of limits down to masses of M Z = 200 GeV. Nevertheless, this seems to be justified as the ATLAS resonance search under consideration was designed for heavy mediators [12]. Thus, in principle, a dedicated analysis in the low invariant mass range should yield better limits than those given in the analysis [12] we used.
Furthermore, we have projected the limits deduced from the ATLAS search at √ s = 8 TeV and 20.3 fb −1 of data to a Run 2 (LHC Run II) scenario with √ s = 13 TeV and 100 fb −1 and a high luminosity scenario (HL-LHC) with √ s = 13 TeV and 3000 fb −1 . For this purpose we have rescaled the exclusion limits on the cross section by the respective luminosities σ (13) assuming that scaling of the limits is only due to statistics. As a cross-check for this prescription to work, we have compared luminosity rescaled limits from the ATLAS 8 TeV analysis to the preliminary limits from the ATLAS analysis at 13 TeV [18] in the eµ-channel. As shown in the right panel of Fig. 4, this method seems to give sufficiently accurate results for the purpose of a rough projection. The projection to the LHC Run II scenario is illustrated by the red dash-dotted line in Figs. 1 to 3, the projection to the HL-LHC scenario by the red dotted line.

Constraints
In this section we want to give a brief overview over the different constraints that already restrict our model. We will first consider pure quark-sector, then mixed quark-and lepton-sector and finally pure lepton-sector constraints.

Meson mixing
One of the strongest probes of flavor violation in the quark sector is provided by meson mixing, where a meson M oscillates into its conjugate stateM. In the Standard Model these processes are loop-suppressed as they require FCNCs and thus the matrix elements are rather small. Therefore, meson mixing is very sensitive to new physics models that have flavor-changing couplings in the quark sector. This is the case for our model where meson mixing arises at tree level via Z exchange from the diagrams shown in Fig. 5.
As the mass splitting ∆M M of the conjugate meson states is directly proportional to the transition matrix element M 12 , it is the appropriate observable for testing flavor violation. As meson mixing is a low-energy effect we follow Ref. [19] and investigate the Z effects in an EFT approach, where we will integrate out the Z at the high scale µ in ∼ M Z . The resulting four-quark operators describing the low-energy phenomenology of the Z -induced FCNCs are given by [19,20] It should be noted that the operator O LR 2 is generated only through QCD-loop effects from operator mixing due to the running of operators from the high to the low scale.
After matching the operators of Eqs. (3.2) to (3.5) to the full theory we find with the off-diagonal matrix element M * 12 = M|H ∆S=2 eff |M for the mass splitting [19] 3 Lepton decays Muonium LEP II Figure 6: The flavor-violating couplings g R bs and g R eµ for a Z boson of M Z = 750 GeV with purely right-handed quark couplings (ρ q = 0). The red area indicates the limit from the ATLAS analysis of the process pp → eµ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The green area is the limit coming from the meson decay B 0 s → eµ. The gold and magenta areas are purely leptonic limits coming from LEP and muonium oscillation constraints. The black dotted line is the neutrino exclusion limit from B + → K +ν ν. The dash-dotted line originates from tests of lepton decays. However, these latter two limits only apply for left-handed lepton couplings. The gray area represents the B s −B s mixing limit.
where the P i denote the hadronic matrix elements corresponding to the operators O i . We calculate the hadronic matrix elements for K, B d and B s mesons mainly from the relations given in Refs. [19,20] and the lattice bag parameters from quenched QCD calculations given in [21,22]. For D mesons we rely on the relations given in Ref. [23]. The R i (µ) are the renormalization group evolution coefficients encoding the running of the operators O i due to NLO QCD effects. They are normalized such that R i (µ in ) = 1 at the scale µ in where the Z is integrated out. The coefficients are given by [19] R VLL From the measurement of the meson mass splitting ∆M M we can derive a limit on the coupling g R ij . In practice, we use the values provided by the UTfit collaboration [24] and more specifically the maximally allowed deviation between the measurement and the SM prediction.
In Fig. 6 we show as an example again the limits for the combination {qq , } = {bs, eµ}. We consider purely right-handed quark couplings, i.e. ρ q = 0, and include the mixing limit as the gray area. It can be seen that the mixing limit in the quark sector, in general, is so strong that it excludes the whole region of parameter space that can be probed with multipurpose detectors at the LHC.

Cancellation
We have just seen that for a general quark coupling configuration it seems hopeless to test flavor violation with ATLAS at the LHC. However, there is an important subtlety to these considerations that alters the picture just enough to serve as motivation for a search of Z induced FCNCs at the LHC.
The term in brackets in Eq. (3.6) is a quadratic form in the parameter ρ q . Thus, if the discriminant is greater than zero, we have two solutions ρ 0 for which the mass splitting due to Z exchange vanishes exactly. We will only consider the solution that has mostly right-handed couplings. The other solution is simply given by 1/ρ 0 and has essentially the same behavior.
As we are interested in comparing flavor bounds on the Z couplings to collider bounds, we additionally define an upper and a lower tolerance ∆ρ − and ∆ρ + such that the upper limit on the coupling g lim ij derived from mixing is less stringent than some reference limit g * , i.e.
With this definition we can find an interval I 0 around the exact cancellation point ρ 0 , where the limit due to mixing is subdominant compared to the reference limit g * = min(g col ij ), i.e. to the strongest bound we can set on g q ij from the ATLAS search at a given mass M Z . For example, the exclusion plot in Fig. 1 shows the limits for a Z boson with M Z = 750 GeV coupling to bs and eµ. As can be seen in the plots the limits hold for ρ 0 = 0.1324 with ∆ρ − = 0.0005 and ∆ρ + = 0.0019, i.e. the tolerance interval in this case is I 0 = [0.1319, 0.1343].
So far we have treated cancellation effects only at tree level. However, we have to make sure that these effects are persistent even at higher orders. Therefore, we investigate the impact of one-loop corrections on the cancellation solution ρ 0 and the associated tolerance interval I 0 .
The details can be found in Appendix A. We find that the cancellation solution ρ 0 and the interval I 0 both receive an overall shift δρ at one-loop level. However, this shift is negligible in the sense that δρ < ρ 0 in the region of parameter space that is probed by LHC bounds.

Meson decays
Another process in the meson sector that directly constrains our model is rare neutral meson decays M 0 → + − , where M 0 can be the K, D, B d or B s meson. These decays involve two flavor-changing vertices and therefore are highly suppressed in the SM, whereas we can generate these processes on tree-level in our Z model. From the Lagrangian in Eq. (1.1) we can immediately construct the relevant four-fermion operators [25] with corresponding Wilson coefficients (3.14) With knowledge of the relevant operators and their associated Wilson coefficients we can calculate the branching ratio BR(M 0 → + − ) for the different mesons [25,26] BR where f M 0 , M M 0 and Γ M 0 are the decay constant, the mass and the total width of the decaying meson. Furthermore, we have assumed that is the heavier of the two leptons and we have neglected the mass of the other one. Of course, these limits only exist when the mass of the relevant meson is bigger than the combined mass of the two leptons. Finally, based on Eq. (3.15) we can derive a limit on the Z coupling where we have made use of the relations in Eq. (2.11). The corresponding limits due to meson decays are depicted in Figs. 1 to 3 by the green area. The power of the decay limits comes from the fact that it constrains the product g qq g . Therefore, meson decays can probe regions in parameter space where one of the two coupling is very small while the other one is big, a region hard to probe at the LHC. However, the LHC limits are generally more stringent in the direction of parameter space where both couplings become small but are of comparable size. Especially with Run II or HL-LHC data one can expect rather big gains along that direction. Figure 7: Leading order SM contributions to the decay K + → π + νν. Figure 8: Leading order Z contribution to the decay K + → π + νν.

Neutrino limits
It seems reasonable that an extension of the Standard Model should preserve the SU(2) L gauge symmetry at high energies. In return this means for our effective model that the Z gauge boson should couple to the quark and lepton doublets Q L and L L , if left-handed couplings are present. In this scenario the Z couples to neutrinos ν i and ν j with equal strength g L i j as to charged leptons i and j . However, the coupling to neutrinos opens up a whole new class of constraints to our model. Especially, meson decays of the form M 0,± 1 → M 0,± 2ν ν can be a sensitive probe for the presence of left-handed Z couplings.
In this section we will now investigate such decays for the different neutral mesons that can have an impact on our model. In particular we consider decays of kaons and B-mesons. The corresponding measurements in the D-meson sector are not yet very restricting.

Kaons
When we couple the Z to the first two quark generations we can constrain the left-handed lepton couplings from the kaon decay K + → π + νν. In order to extract a constraint on g L i j we will calculate in the following the branching ratio BR(K + → π + νν). For a detailed derivation of the branching fractions see Ref. [27]. The relevant leading-order SM diagrams contributing to this decay are shown in Fig. 7. We can see that the leading SM contributions are already loop-suppressed. Thus, with the leading order Z contribution being a tree-level effect as shown in Fig. 8, one can expect that the Z can give a sizable contribution to the branching fraction of this decay. Adopting the notation of Figure 9: Decay K + → π 0 e + ν e in the SM.
Ref. [27], where (sd) V ±A ≡sγ µ (1 ± γ 5 )d, the relevant operators for the low energy interaction read with corresponding Wilson coefficients In order to calculate the branching ratio we will make use of isospin symmetry [27] to extract the hadronic matrix element for (sd) V −A from the decay K + → π 0 e + ν e shown in Fig. 9, Additionally we take the hadronic matrix elements of the left-and right-handed currents to be the same [28] as the process of interest is purely governed by QCD and therefore should be independent of the underlying chirality structure. The effective operator for the process in Fig. 9 reads where V denotes the CKM matrix. Thus, we can write for the branching ratios To turn this into a limit on g R sd we only consider the part of the branching fraction not explained by the SM where we used for our analysis Finally, we obtain the constraint on the leptonic coupling to be given by The resulting limits on the Z couplings are depicted by the black dotted line in the lower left panels in Figs. 20 to 22 in Appendix C. The neutrino limits in the kaon sector are quite strong and like the limits from decays into charged leptons constrain the product g qq g . Especially in the eτ and µτ sector the kaon-neutrino limits exclude all regions of parameter space that can be hoped to be tested at the LHC. However, as mentioned in the beginning of this section, the neutrino limits are only valid if we take the lepton couplings to be left-handed and can be fully circumvented by only considering right-handed lepton couplings.

B mesons
We can use the transition b → dνν analogously to the case of the kaons. However, in this case the process that is induced through non-zero g bd and g i j couplings is B 0 → π 0ν i ν j . This is shown in the left panel of Fig. 10. As in the case for the kaons a limit on the coupling g bd can be derived from comparing the branching fraction of the latter B 0 decay to the one for the decay B 0 → π − e + ν e (shown in the right panel of Fig. 10). With the PDG values [31] for the respective branching ratios Figure 11: Z contribution to the decay B + → K + νν.
we can then set a limit in analogy to Eq.
To constrain the transition b → sνν we can use the decay B + → K + νν. To extract a limit on the the Z coupling g bs in the presence of left-handed lepton couplings we will again need the Z contribution to the branching ratio of this decay [32]. First, we can parametrize any contribution to this process in an EFT approach by the effective operators As for the kaons, the leading order SM contributions are coming from electroweak loop diagrams, which therefore are only involving left-handed fermions. The corresponding SM Wilson coefficient has been calculated [32] and can be written as with s w = sin θ w denoting the sine of the Weinberg angle. The leading-order contribution coming from the exchange of a Z as shown in Fig. 11 yields the Wilson coefficients (3.33) Defining the differential branching fractions for the process as one finds for the ratio [32] with the model-independent quantities From the constraint that R K ≤ 4.3 [32] we can then derive a limit on the quark coupling

Lepton decays
Another important leptonic constraint which involves neutrinos is due to charged lepton decays. Again this applies only to Z bosons with couplings to left-handed leptons. As the process involves only leptons, the limits only depend on the leptonic couplings and will therefore correspond to vertical lines in the plots.
According to the diagram in Fig. 12 there is an additional decay contribution to the ordinary non flavor-violating lepton decay iL → jL ν iLνjL . This contribution interferes with the SM contribution generated by a W boson exchange. In addition, there are three new decay channels iL → jLνiL ν jL , iR → jR ν iLνjL (ν iL ν jL ) . These channels arise from diagrams as shown in Fig. 12 and the ones with neutrino flavors interchanged. Together they modify the SM decay rate [14] into a given lepton plus neutrinos according to, where with the weak coupling constant denoted as g, The first part is the contribution to the SM-like purely left-handed, non-flavor-violating channel. The second part are the non-SM-like chirality-flipped and/or flavor-violating channels. Importantly we do not distinguish the neutrino species in the measurement of the final state. This is why all the contributions are summed. As before, we will consider the case of purely left-or right-handed lepton couplings. In the case of purely right-handed couplings we do not get any contribution from the Z (as due to g L i j = 0 also x ij , y ij = 0). In the case of purely left-handed couplings the modification of the SM decay rate simplifies to Figure 12: Diagrams of the lepton decay i → j ν iνj . The same diagram also exists for right-handed charged leptons and for the lepton flavor-violating decays i → j ν jνi , where the neutrino flavors are interchanged.

µ decays
Measurements of the µ lifetime are very precise with a relative uncertainty of the order of 10 −6 [31]. This suggests very tight constraints on x µe . However, usually the decay of the µ is used to determine the Fermi constant G F . Therefore this measurement cannot be used anymore to test new physics. To do so we need an additional measurement. The β decay of nucleons is possible only via a charged current and is therefore unaffected by our Z . However, it contains the CKM matrix element V ud which is usually extracted from those decays. The situation is similar for the decay of kaons which contain the matrix element V us . Nevertheless we can extract a limit from this comparison using the following argument. As can be seen from Eq. (3.40) the Z contribution leads to an increase in the decay rate. Using the SM extraction of V us and V ud this would lead to smaller values of these CKM matrix elements. Assuming unitarity for the CKM matrix we can then constrain x µe using the CKM matrix elements determined in the standard way [31], On the right hand side we estimate the error by adding the errors for V us and V ud in quadrature and in the next step adding the small deviation from 1.
The resulting limits are shown as dash-dotted lines in the figures and again only apply if we take the lepton coupling to be purely left-instead of right-handed.
In principle one could also derive limits from the angular dependence of the decay of polarized muons used to search for right-handed currents in Ref. [33]. However, for these constraints to be effective requires the presence of both left-and right-handed couplings which we do not consider in the eµ sector.

τ decays
In the following we want to discuss the impact of our model on various τ decay modes. We will treat the leptonic and hadronic decay modes separately.

Leptonic mode
Due to our choice of having only a single flavor-changing coupling in the lepton sector, either eτ or µτ , a strong constraint can be obtained by comparing the branching ratios in these two channels [14] (in [34] a comparison with the SM branching ratio is used).
For definiteness, let us consider the case of a non-zero µτ -coupling first. In the following we strongly rely on the derivation done in Ref. [14]. A non-zero coupling g µτ leads to an enhancement of the partial decay rate of the process τ → µνν according to Eq. (3.38). Defining the ratio of the partial decay rates corresponding to τ → µνν and τ → eνν, we can rewrite Eq. (3.38) as Within the SM the ratio R µ/e has been very accurately calculated [35], For the experimentally determined value we follow [14,35] and quote a precise measurement 5 by the BaBar collaboration [38] yielding a value of R µ/e = 0.9796 ± 0.0039 . Recalling Eq. (3.43) we can use this observed deviation as a constraint on our model. If we allow only for a left-handed coupling g L µτ (which is necessary for the decay into neutrinos to take place), we find the limit In the presence of both left-and right-handed lepton couplings we can derive a limit as  Figure 13: The flavor-violating couplings g R bs and g R µτ for a Z boson of M Z = 1000 GeV and a coupling ratio of ρ q = g L bs /g R bs = 0.1276. The red area indicates the limit from the ATLAS analysis of the process pp → µτ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The green area represents the excluded region from the decay B 0 s → µτ . The black dotted line is the limit from B + → K +ν ν applying only for left-handed lepton couplings. The light and dark cyan areas depict the 4 and 5 σ exclusion bands from ∆a µ . The black and white hatched area is the preferred 1 σ region of the observed deviation ∆R µ/e = R µ/e /R SM µ/e − 1, also only applicable in the case of purely left-handed couplings.
Now we want to consider the case of non-zero eτ -coupling. In this scenario the arguments are essentially analogous to the µτ -case, however, using the inverted ratio It is worth noticing that in this case using the inverted ratio the relative deviation of the experimental value from the SM prediction as discussed in Eq. From Eq. (3.38) we see that structurally the Z contribution from a non-vanishing g eτ coupling will always lead to a positive shift in R e/µ . Therefore, the measured fluctuation leading to a negative shift will impose an anomalously stringent bound on g eτ (cf. Fig. 21 in the Appendix).
In both the case of non-zero eτ -and µτ -coupling we use the observed relative deviation ∆R plus the 2 σ uncertainty as exclusion bound. The corresponding limits are shown for example in Figs. 1 to 3 as the vertical black dash-dotted line. It can be seen that lepton decay limits are by far the strongest purely leptonic limits and cut far into the region of parameter space testable with multipurpose experiments at the LHC. However, it has to be noted that these limits apply only if we consider purely left-instead of right-handed lepton couplings g L i j . In the right-handed case these limits are absent.
If we assume the previously discussed 1.8 σ relative deviation in the ratios of branching fractions Eq. (3.46) not to be due to systematics or a fluctuation, we can speculate on a possible new physics origin. In order to justify such speculation, we checked with help of the accurate prediction of R SM µ/e given in Ref. [35] and the measured τ lifetime that the excess ∆R µ/e is indeed due to the observed value of Γ τ →µνν being significantly higher than its SM prediction. Indeed, in previous work [39] it has been noticed that the relative deviation even amounts to 2.4 σ. The resulting increase in the total width Γ tot is compatible with its observed value. Therefore, we additionally fitted the excess within our model at the 1 σ level. The preferred region is depicted by the black and white hatched area in Fig. 13. We can see that we can fit the observed deviation for a 1 TeV Z boson with moderate couplings g L µτ ∼ O(10 −1 ). However, this only applies in the case of purely left-handed lepton couplings also implying that the limit from semi-leptonic meson decay into neutrinos applies (dotted line). In return, this meson decay limit excludes most of the region in parameter space testable with ATLAS or CMS at the LHC.

Hadronic mode
In the context of τ decays we want to mention the hadronic decay mode τ ± → ± (π 0 K 0 /π ± K ∓ ) due to the diagrams shown in Fig. 14, where ∈ {e, µ}. This mode is present in the case of non-vanishing quark couplings in the sd sector 6 . As it is not present in the SM the detection of such a decay would be a smoking gun for a doubly flavor-changing Z . It should be noted that such a decay is only possible into pions and kaons as B and D mesons are too heavy for the τ to decay into.
The decays into charged mesons τ − → − π ± K ∓ and τ + → + π ± K ∓ have been searched for at BaBar [41] and Belle [42]. The corresponding branching fraction limits can be turned into a limit on the Z couplings g τ and g sd . The relevant operators contributing to this process are again those of Eqs. (3.12) and (3.13) with corresponding Wilson coefficients given in Eq. (3.14). In order to calculate the branching fraction due to the Z induced decay we need the hadronic matrix elements π + K − |(sd) L/R |0 , where we have introduced the shorthand (sd) L/R =s γ µ P L/R d. Similar to Section 3.3.1 we will use isospin symmetry to estimate the qd Z s q τ − Figure 14: Diagrams of a possible signature of a doubly flavor-changing Z boson with non-zero couplings in the sd and eτ /µτ sector (q denotes either a d-or a u-quark). The resulting decays are τ ± → ± (π 0 K 0 /π ± K ∓ ), where ∈ {e, µ}. Only in the sd sector such a semi-leptonic τ decay into a pion and a kaon is kinematically allowed. Especially for a non-zero quark coupling involving a b-quark such a decay is not possible. matrix element from the observed decay τ → ν τ K − π 0 . The SM operator for this decay reads In the following, we assume the electron and muon to be massless, which seems to be sensible as m e , m µ m τ . However, this implies that the outgoing leptons have definite handedness and consequently lead to distinguishable final states. Therefore, the ratio of branching fractions reads Treating the leptons as massless translates into |(¯ τ ) L |τ ν τ |(ν τ τ ) L |τ . Furthermore, we use isospin symmetry to relate the hadronic matrix elements As QCD is a non-chiral theory, the same relation holds also with the operator (sd) L replaced by (sd) R . Put together, we can thus estimate the branching fraction from Eq. (3.53) to be This relation allows us to put a limit on g τ and g sd . Using the Belle limits [42] Γ τ − →e − π + K − < 5.6 × 10 −8 , (3.56)

Muonium -Antimuonium oscillations
The  (3.59) where X, Y ∈ {L, R} and the Wilson coefficients In order to calculate the amplitude we will perform a non-relativistic field expansion and calculate the resulting effective potential from the Born approximation. Muonium is a non-relativistic Coulomb bound state and therefore the transitions M −M can be described by a non-relativistic effective potential V eff ( x). Taking into account that the two fermions can either be in a spin singlet or triplet bound state we obtain for the potentials Assuming the muonium to be in its electronic ground state we can calculate the mass splitting [44] with the Bohr radius of the muonimum a MM = 1 α m red and the reduced mass m red = me mµ me+mµ . In order to get in contact with the experiment, we need to know the transition probability P MM of an initially prepared muonium atom M to oscillate intoM . Therefore, we need to know the time evolution of the two-state system that is generally obtained by solving the Schrödinger equation [45] i d dt After diagonalizing the Hamiltonian of Eq. (3.65) one finds for the time evolution of an initially pure muonimum state If the bound state is mainly antimuonium it consists of e + and µ − . The muon will decay via µ − → e −ν e ν µ with a highly energetic electron e − hitting the detector. For muonium it would instead be a highly energetic positron. The measurement principle to detect that muonium has oscillated into antimuonium is to start with a pure muonium initial state and then look for the highly energetic electron resulting from the muon decay inside the antimuonium bound state. In essence this means counting the number of outgoing energetic e − . Their number can be determined by integrating the partial decay rate into electrons. This is given by the probability that the system is in an antimuonium state multiplied by the muon decay rate. Integrating over time we have [46], Assuming either left-or right-handed lepton couplings this can be translated into a rough limit on the off-diagonal lepton coupling of the Z with S B = 0.35 [44] a correction factor for the muonium splitting in the magnetic field coming from the (V ± A) × (V ± A) Lorentz structure of the interaction. The corresponding limits are shown for example in Figs. 1 and 3 by the magenta band. As these limits are purely leptonic they result in a vertical exclusion line. In general, it is the strongest limit for the lepton sector only. In contrast, the LHC can probe eµ couplings significantly smaller than those excluded by muonium. However, this is only true in combination with same order of magnitude quark-sector couplings whereas the muonium limits are universal. Furthermore, as muonium is a bound state of e + and µ − the eτ and µτ sectors are completely unaffected by this limit. The relevant eτ or µτ bound states in these cases would be difficult to access as the τ decays very rapidly. (Right) The blue curve is the total simulated cross section σ Z for the process e + e − → τ + τ − at √ s = 207 GeV. The green band depicts the 95% confidence interval of the measured cross section σ exp . The red area shows the excluded couplings g R eτ .

LEP limits
The LEP collider run at CERN from 1989 to 2000 produced a large amount of clean e + e − collisions. In particular the analyses of the processes e + e − → µ + µ − and e + e − → τ + τ − provide constraints on the Z couplings g eµ and g eτ . In order to derive constraints on these couplings in our model, we use the total inclusive cross sections σ(e + e − → µ + µ − ) and σ(e + e − → τ + τ − ) as measured by the ALEPH collaboration [47]. Therefore, we have simulated the total inclusive cross section σ Z of these two processes with madgraph v2.3.3 [16] at √ s = 207 GeV, including the Z diagram shown in the left panel of Fig. 16. We have simulated σ Z for a number of different values of the lepton couplings and the Z mass allowing for an additional hard photon in the final state. We then set a limit on g eµ and g eτ by using a two-sided hypothesis test. For this we assume that the measured cross section σ exp (i.e the number of signal events) follows a Gaussian distribution. For the relative large total number of events N µµ = 683 and N τ τ = 402 at √ s = 207 GeV [47] this seems to be justified. For a given Z mass we can then scan the simulated total inclusive cross section σ Z for the lepton couplings g eµ and g eτ and exclude all coupling values that correspond to a cross section σ / ∈ [σ exp − 1.96 ∆σ, σ exp + 1.96 ∆σ], which corresponds to a two-sided 95% confidence interval for a Gaussian distribution.
The corresponding limits are found e.g. in Figs. 1 and 3 and are depicted by the golden regions. As was the case for the muonium, these limits are purely leptonic (and therefore correspond to vertical bands in the g qq − g plane). As LEP was an electron-positron collider, these limits do not concern the µτ sector. Generally, the LEP limits are weaker than the muonium limits and therefore are only of concern in the eτ sector (cf. Appendix C Fig. 21) where the muonium limits are not present. A particular feature of the LEP limits is the gap in the excluded region of parameter space. The origin of this gap can be understood with help of the right panel of Fig. 16. For small couplings g e the total cross section is mainly SM-like and agrees very well with the measurements (i.e. it lies within the 95% confidence interval). For moderate couplings g e 1 the interference term, which is linear in g e and has negative sign, starts to become important and eventually drives the cross section σ Z below the confidence interval. This leads to the first exclusion band. With increasing couplings g e > 1 the pure Z contribution, which is quadratic in g e , starts to dominate and drives the cross section σ Z well above the 95% confidence interval. This leads to the second exclusion band. In between those two regimes we have a transition region where σ Z lies within the 95% confidence interval -the gap in the exclusion region.

Magnetic dipole moments and the Z
A further possible constraint for the lepton sector of our model is coming from the measurements of (g − 2) of the electron and the muon. However, as (g − 2) µ exhibits a deviation between theory and experiment of about 3 σ or even more [48][49][50][51] it is tempting to speculate on a new physics origin. Hence, in this section we additionally want to explore the potential of the Z boson to play this role similar to earlier work [13,14,52,53].

Experimental status
We want to motivate our discussion by looking at one of the most precise measurements in the electroweak precision era: the determination of the gyromagnetic ratio g of the muon at the E821 experiment at the Brookhaven Alternating Gradient Synchrotron [54]. The naive SM tree-level calculation, i.e. the Dirac equation, yields a value for the gyromagnetic ratio of g = 2. Radiative corrections such as higher-order QED processes, electroweak loops or hadronic vacuum polarization lead to a shift of the gyromagnetic ratio, the so-called anomalous magnetic moment Much interest has been triggered by the findings of the E821 experiment that point towards a mismatch between theory [50,51] and experiment [48] of up to ∼ 3.6 σ or ∆a µ = a exp µ − a SM µ = (2.87 ± 0.80) × 10 −9 . (3.70) For a suitable chirality structure 7 the anomalous magnetic moment of the muon can get a positive shift due to radiative corrections from a Z loop through a nonzero µτ -coupling. Therefore such a Z can potentially reconcile the experimental value with the theory prediction. In previous work [13,55] models have been studied, where the Z boson couples to the L µ − L τ current, with Q being the overall lepton charge, L 2 = (ν µ , µ L ) and L 2 = (ν τ , τ L ). In such models an explanation of the (g − 2) µ tension is ruled out for Z bosons with mass M Z GeV by neutrino trident production ν N → ν µ + µ − N in the Coulomb field of a nucleus N [13,52,53,56]. Nevertheless, neutrino trident production is not possible in our model as it requires diagonal couplings to µ and τ of the Z on tree level. Consequently neutrino trident constraints are not Figure 17: Anomalous magnetic moment of fermion f a due to Z exchange and f b running in the loop.
applicable because they are looking at two muon/tau signatures. Hence, in this case (g − 2) µ can be explained by a suitable pure µτ coupling [14]. In the future it may be possible to look for flavor changing trident signals with eµ, eτ or µτ in the final state at the DUNE [57] and SHiP [58,59] facilities [60,61]. New physics can also be probed via the electron anomalous magnetic moment. Here the picture is different as with a deviation of only ∼ 1.3 σ between theory and experiment [62], the experimental result is in good agreement with the theoretical prediction and The uncertainty in ∆a e is expected to be reduced in the near future, enhancing its potential as a test of new physics. One subtlety that has to be taken into account is the fact that often the value of the fine structure constant α EM is deduced from the electron magnetic moment measurement.
To be sensitive to new physics in a consistent manner one should use in the calculation of the magnetic moment a value of α EM determined by another independent measurement, as for example by interferometry of rubidium atoms [63,64]. In light of these prospects we will also investigate the shifts of the anomalous magnetic moment a e of the electron due to Z loops.

Z contribution to (g − 2)
We will now briefly recall the calculation of a general Z contribution to the anomalous magnetic moment a in our model. The fermions that are coupled via flavor-changing interactions to the Z generically receive a contribution to their anomalous magnetic moment from the diagram of Fig. 17. More specifically, we consider the Z interaction for two generic fermions f a and f b and C A = (g R ab − g L ab )/2 the calculation 8 yields Assuming only right-handed couplings (C V = C A = g R ab /2), we can use this relation to turn the observed shift ∆a in the electron/muon magnetic moment into a limit on g R ab . In order to get a better understanding of the Z contribution to the anomalous magnetic moment we derive an approximate formula. As we are mostly interested in Z bosons in the multi-GeV range we assume M Z m a , m b . Therefore, we expand Eq. (3.74) for small ratios x a and x b . Keeping only the leading powers yields the approximate formula As mentioned above we use the shift in the electron magnetic moment ∆a e to constrain the off-diagonal couplings g eµ and g eτ . From Eq. (3.75) we can see that in the case of purely right-handed lepton couplings (ρ = 0) the contribution to the electron magnetic moment a e is suppressed compared to the muon magnetic moment a µ by a factor (3.76) Comparing Eq. (3.72) and Eq. (3.70) we see that the precision of ∆a e is only four orders of magnitudes smaller than the one of ∆a µ . Hence, the constraints from ∆a e are much weaker than those from ∆a µ and play only a role for very light Z bosons 9 . On the other hand, for light Z bosons we obtain rather strong limits from either LEP or muonium-antimuonium oscillation (cf. Figs. 20 and 21). Therefore, the constraint from a e proves to be practically irrelevant for our model. For the µτ sector the situation is more complicated. For purely right-handed lepton couplings the Z contribution goes in the wrong direction compared to the measurement (cf. [14,53]). As the current deviation between SM prediction and the measured value is greater than 3 σ any contribution of such a Z is ruled out at this level. Therefore we show exclusions at the 4 and 5 σ level e.g. in Figs. 2, 13 and 22 as light and dark cyan bands. These limits are the only purely leptonic constraints in the µτ sector.

A hint for (g − 2) µ
In view of the tension between the measured value and the theory prediction of a µ we will now discuss the implications of a Z with non-zero g µτ couplings on the muon magnetic moment. In order to explain the observed positive shift ∆a µ of Eq. (3.70), we also need a positive Z contribution to a µ . This is the case when the term in square brackets in Eq. (3.75) is positive 10 . This is a quadratic form in the coupling ratio ρ and is positive only in between its two roots i.e. for 0.02 ρ 50.75. In this section we therefore consider a vector coupling scenario (ρ = 1) and an optimized scenario where the limits from τ decay are weakest while the Z contribution to a µ is still positive (ρ = 0.053).

Vector coupling scenario
The left panel of Fig. 18 shows the g R bs − g R µτ plane in a vector couplings scenario (ρ = 1) for a Z with a mass of 1 TeV. First, we notice the absence of a limit from the leptonic meson decay B 0 s → µτ . The absence of such a limit is a peculiar feature of the vector coupling scenario in the lepton sector. This can be understood with help of Eq. (3.16), which features a term |1 − ρ | in the denominator and consequently diverges at ρ = 1. Second, the limit from the meson decay B + → K +ν ν now becomes unavoidable due to the non-zero left-handed lepton coupling g L µτ . This limit (shown in yellow) is much stronger than current 8 TeV ATLAS limits (shown in red) and possibly even stronger than limits from an LHC Run II scenario (red dash-dotted line). Even a future HL-LHC run could only slightly improve this limit along the direction of both small quark and lepton sector couplings. As mentioned before we now get a positive Z contribution to a µ . Instead of using the observed deviation ∆a µ as a limit we can fit the excess. The purple, blue and green bands show the preferred 1, 2 and 3 σ regions of ∆a µ . It is worth noticing that for a 1 TeV Z the excess can naturally be accommodated with O(1) lepton couplings g µτ . However, one has to be careful whether the limits from τ decay rule out such a (g − 2) µ explanation. On the 2 σ level this is indeed the case [14]. As done in Fig. 13 we fit the observed deviation ∆R µ/e at the 1 σ level. This fit is shown by the black and white hatched area. The observed ∆a µ deviation is still compatible with the τ decay excess within 3 σ.

Optimized coupling scenario
The right panel of Fig. 18 shows the g R bs − g R µτ plane for a lepton coupling ratio of ρ = 0.053 for a Z with a mass of 200 GeV. This scenario is optimized such that for a positive Z contribution to a µ the limit of τ decays is weakest. Previously, Altmannshofer et al. have shown explicitly in Ref. [14] that for ρ = 0.1 an explanation of ∆a µ is not ruled out by τ decay limits for Z masses greater than a few GeV. Comparing to the vector coupling scenario we can see that we get a limit from the leptonic meson decay B 0 s → µτ . In addition, the relative strength of the limit from the meson decay B + → K +ν ν to the ATLAS limits is much weaker. This is due to the    Fig. 18 we can see that relevant parts of the parameter space of this scenario can be probed possibly already with LHC Run II data and definitely with a HL-LHC run.
Future more precise measurements of the branching fraction of the τ decays τ → µνν and τ → eνν, e.g. at the Belle-II experiment [66] as well as the factor of four improvement in the precision of (g − 2) µ in the upcoming E989 experiment at Fermilab [67] can test this interpretation. In addition to the purely leptonic tests the presence of sd type quark couplings could present an opportunity to test this model in unusual τ -decays, as discussed in Sect. 3.4.2. For example, assuming the maximally allowed quark coupling of g R sd ≈ 3 × 10 −3 for the 200 GeV Z discussed in this section yields a branching fraction of Γ τ − →µ − π + K − ≈ 9.0 × 10 −9 . This could directly be searched for at Belle-II, which aims at a sensitivity of 1 × 10 −9 in branching fraction for 50 ab −1 of data [66].

Summary
In this paper we have investigated simple test models of Z bosons with exclusively flavorchanging interactions, one in the lepton and one in the quark sector. For such models usually one would expect that precision tests of flavor-changing neutral currents are far superior to direct production at the LHC. The latter could then be taken as just a nice confirmation of what we already know. For a generic chirality structure of the couplings this is indeed the case and LHC limits are eclipsed by limits on meson mixing as one can see for example from Fig. 6. However, the latter limits depend on the relative strength of right-and left-handed couplings and there exist small regions where they can be evaded. Here, the chirality independent LHC limits take over and become the best probe of new physics. A similar situation arises with limits of mesons and leptons decaying into neutrinos (cf. e.g. [31,32,68,69]). These limits are applicable for left-handed couplings, but can be evaded for purely right-handed ones.
Due to the coupling to leptons our Z boson also gives a contribution to (g − 2). For couplings to µτ and a suitable chirality structure this allows for an explanation of the deviation in (g − 2) µ from the SM expectation (cf. also [14]) as well as a small excess in the decay of τ leptons into muons and neutrinos. In the future measurements of a Z decaying into µ and τ at ATLAS or CMS, or of B decays at LHCb can probe into this parameter space. In particular a dedicated ATLAS or CMS search at kinematics suitable for a relatively low mass resonance could be helpful. As this can only test part of the interesting region it is worthwhile to look for complementary probes. Here the study of τ -decays, as it can be done, e.g. at Belle-II [66], provides for interesting opportunities. For purely leptonic couplings in particular precision tests of lepton universality in these decays seem promising. In addition, flavor-violating trident production at high intensity experiments like DUNE [57] or SHiP [58,59] may allow to test this region [60,61]. Furthermore, additional couplings to relatively light quarks could provide for striking signals in unusual τ -decays into µ+hadrons.

A Higher order effects in cancellation
In Section 3.1 we have seen how meson mixing arises from four-quark operators that are generated in our model at tree level. Furthermore, we have argued that there are solutions of the coupling ratio ρ q for which the mixing exactly cancel. In the following we will investigate higher order effects contributing to mixing and its impact on the cancellation. Figure 19: One-loop mixing diagram giving rise to higher order corrections of the four-quark operators Eqs. (3.2) to (3.4). The same diagrams exist also with a Higgs boson instead of the additional Z. A third correction comes from the right diagram with the Z exchanged by a second Z .

A.1 NLO effects
Diagrams like the ones in Fig. 19 will give rise to corrections of the Wilson coefficients C i of the four-quark operators in Eqs. (3.2) to (3.4). We define the corresponding Wilson coefficients as In the following we study the various contributions in more detail. Therefore, we consider diagrams of the type of Fig. 19, where the additional boson in the loop can either be the SM h, Z or a second Z .

Higgs loop contribution
First, we will consider a SM Higgs h as the additional boson in the loop. The important feature of the Higgs is that it flips the chirality of the fermion at the vertex. Therefore, the Higgs introduces operator mixing in the sense that the Higgs correction to the Wilson coefficient of e.g. operator O VLL 1 will be proportional to O LR 1 . A short calculation of the contributions from the one-loop diagrams in Fig. 19 with the Z replaced by the Higgs leads to the estimate where we have used m q = 1 GeV and M Z = 1 TeV to get the last relation. Hence, operator mixing due to NLO Higgs exchange is an effect roughly of the order of 10 −7 and therefore much too small to be of any concern as will become clear in the following.

Z loop contribution
Next, we consider the diagrams as depicted in Fig. 19, where a SM Z boson plays the role of the additional boson running in the loop. Structurally the coupling of the Z and Z are the same, so we do not introduce any mixing amongst operators of the kind δC i ∝ C j . Analogously, an order of magnitude calculation yields for the correction to the Wilson coefficients where we have assumed for the coupling of the quark to the Z a conservative value of g Z q = 0.1 and M Z = 1 TeV to get the last relation. First, we note that the correction coming from Z contributions is much bigger than the one for the Higgs. Second, we note that the correction of the Wilson coefficient is proportional to the Wilson coefficient itself due to the absence of operator mixing. This implies that the correction will be universal to all Wilson coefficients and therefore only shift the cancellation solution ρ 0 and thereby not destroying it.

Z loop contribution
Finally, we consider the case of a pure Z induced loop diagram. One has to note that only the right diagram of Fig. 19 will contribute to generic meson mixing. The diagram on the left only contributes for mesons consisting of a quark-antiquark pair of same flavor. As before we can perform an estimate of the correction to the Wilson coefficients yielding This contribution is in the same ballpark as the Z contribution, but it introduces higher powers of the Z coupling. This changes the structure of Eq. (3.6) fundamentally and therefore can possibly destroy the cancellation effect.

A.2 Numerical stability of cancellation
We will now examine whether the correction due to a Z loop can spoil the cancellation solution.
As the correction is of the order of δC i /C i ∼ 10 −2 we will assume that the exact cancellation solution can be approximated by a perturbation series ρ = ρ 0 + δρ + higher orders . (A.7) As we have just seen the Wilson coefficients at one-loop level due to Z corrections schematically read Hence the full cancellation equation at one-loop level becomes If we then define the tree-level and one-loop terms as we find for the correction to the cancellation solution from perturbation theory We checked these Z corrections for the B d , B s , D and K mesons. In all cases the corrections are reasonably small for reasonable values of the quark coupling g R qq 1. Especially the ratio δρ/ρ 0 < 1 for all coupling combinations {qq , }. We can calculate by the same method the correction on the tolerance δ(∆ρ). We find that this is generally much smaller than the tolerance itself δ(∆ρ)/∆ρ 1 and therefore negligible. Hence, we obtain a mere shift of the cancellation solution ρ 0 and its tolerance interval I 0 . Therefore, the cancellation solution ρ 0 is stable against higher order corrections and persists beyond tree-level. However it should be noted that in the µτ sector and for high masses in the eτ sector the ratio of the shift to the tolerance δρ/∆ρ can be greater than 1. This is not a problem, as the cancellation solution still persists. It merely means that the shifted cancellation ρ can lie outside of its original tolerance interval I 0 . This is an artefact of the extremely small tolerance interval in those channels.

B Monte Carlo simulation details
In this section we want to summarize important parameters we used for the determination of the simulated cross section σ MC for the process pp → Z → .
In Table 1   Figure 20: Flavor-violating couplings in the eµ sector for a Z boson of M Z = 1 TeV. The red areas depict the limits from the ATLAS analysis of the process pp → eµ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The brown areas show four-quark contact interaction limits from LHC dijet analyses [15]. The green areas are the limits coming from meson decays into charged leptons. The gold and magenta areas are purely leptonic limits coming from LEP and muonium oscillation constraints. The black dotted lines are meson decay limits into neutrinos and the black dash-dotted line the limits from lepton decays. These last two limits, however, apply only for left-handed lepton couplings g L eµ instead of g R eµ . The meson decay limits into neutrinos are absent in the cu sector completely.   [15]. The green areas are the limits coming from meson decays into charged leptons. In purple we show the bounds from the rare decay τ − → e − π + K − only applicable in the sd sector. The gold areas are purely leptonic limits coming from LEP constraints. The black dotted lines are meson decay limits into neutrinos and the black dash-dotted line the limits from lepton decays. These last two limits, however, apply only for left-handed lepton couplings g L eτ instead of g R eτ . The meson decay limits into neutrinos are absent in the cu sector completely. Lepton decays (g−2) µ @ 4σ (g−2) µ @ 5σ (a) Flavor-violating couplings g R bs and g R µτ . Lepton decays (g−2) µ @ 4σ (g−2) µ @ 5σ (b) Flavor-violating couplings g R bd and g R µτ . τ ± →µ ± π ± K ∓ K + →π +ν ν Lepton decays (g−2) µ @ 4σ (g−2) µ @ 5σ (c) Flavor-violating couplings g R sd and g R µτ .
10 -4 10 -3 10 -2 10 -1 10 0 10 1 (d) Flavor-violating couplings g R cu and g R µτ . Figure 22: Flavor-violating couplings in the µτ sector for a Z boson of M Z = 1 TeV. The red areas depict the limits from the ATLAS analysis of the process pp → µτ at √ s = 8 TeV. The red dash-dotted and dashed lines are projections to the LHC Run II and HL-LHC. The brown areas show four-quark contact interaction limits from LHC dijet analyses [15]. The green areas are the limits coming from meson decays into charged leptons. In purple we show the bounds from the rare decay τ − → µ − π + K − only applicable in the sd sector. The light and dark cyan areas depict the 4 and 5 σ exclusion bands from ∆a µ . The black dotted lines are meson decay limits into neutrinos and the black dash-dotted line the limits from lepton decays. These last two limits, however, apply only for left-handed lepton couplings g L µτ instead of g R µτ . The meson decay limits into neutrinos are absent in the cu sector completely.