Another SMEFT story: Z′ facing new results on ε′/ε, ∆MK and K → πνν¯\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \pi \nu \overline{\nu} $$\end{document}

Recently the RBC-UKQCD lattice QCD collaboration presented new results for the hadronic matrix elements relevant for the ratio ε′/ε in the Standard Model (SM) albeit with significant uncertainties. With the present knowledge of the Wilson coefficients and isospin breaking effects there is still a sizable room left for new physics (NP) contributions to ε′/ε which could both enhance or suppress this ratio to agree with the data. The new SM value for the K0 − K¯0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\overline{K}}^0 $$\end{document} mass difference ∆MK from RBC-UKQCD is on the other hand by 2σ above the data hinting for NP required to suppress ∆MK. Simultaneously the most recent results for K+ → π+νν¯\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\pi}^{+}\nu \overline{\nu} $$\end{document} from NA62 and for KL → π0νν¯\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\pi}^0\nu \overline{\nu} $$\end{document} from KOTO still allow for significant NP contributions. We point out that the suppression of ∆MK by NP requires the presence of new CP-violating phases with interesting implications for K → πνν¯\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \pi \nu \overline{\nu} $$\end{document}, KS → μ+μ− and KL → π0ℓ+ℓ− decays. Considering a Z′-scenario within the SMEFT we analyze the dependence of all these observables on the size of NP still allowed by the data on ε′/ε. The hinted ∆MK anomaly together with the εK constraint implies in the presence of only left-handed (LH) or right-handed (RH) flavour-violating Z′ couplings strict correlation between K+ → π+νν¯\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\pi}^{+}\nu \overline{\nu} $$\end{document} and KL → π0νν¯\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\pi}^0\nu \overline{\nu} $$\end{document} branching ratios so that they are either simultaneously enhanced or suppressed relative to SM predictions. An anticorrelation can only be obtained in the presence of both LH and RH couplings. Interestingly, the NP QCD penguin scenario for ε′/ε is excluded by SMEFT renormalization group effects in εK so that NP effects in ε′/ε are governed by electroweak penguins. We also investigate for the first time whether the presence of a heavy Z′ with flavour violating couplings could generate through top Yukawa renormalization group effects FCNCs mediated by the SM Z-boson. The outcome turns out to be very interesting.


Introduction
The ratio ε /ε that measures the size of direct CP violation in K L → ππ decays relative to the indirect CP violation described by ε K and the rare decays K + → π + νν and K L → π 0 νν have been already for many years together with the ∆I = 1/2 rule, K L,S → µ + µ − and K L → π 0 + − decays the stars of Kaon flavour physics [1,2]. The K L -K S mass difference ∆M K remained due to large theoretical uncertainties until recently under the shadow of these decays although it played a very important role in the past in estimating successfully the charm quark mass prior to its discovery [3]. However, recently progress in evaluating ∆M K within the SM has been made by the RBC-UKQCD collaboration [4][5][6] so that ∆M K JHEP12(2020)097 begins to play again an important role in phenomenology, not only to bound effects of NP contributions [7][8][9][10][11][12], but also to help identify what this NP could be. But as stressed in [1,2] and in particular in [13] such an identification is only possible by considering all the stars of Kaon physics simultaneously and also invoking observables from other meson systems.
The RBC-UKQCD lattice QCD collaboration presented very recently new results for the hadronic matrix elements relevant for the ratio ε /ε. Using the Wilson coefficients at the NLO level and not including isospin breaking and NNLO QCD effects they find [14] (ε /ε) SM = (21.7 ± 8.4) × 10 −4 , (1.1) where statistical, parametric and systematic uncertainties have been added in quadrature. However, as already demonstrated in [15], the inclusion of the effects in question, that are absent in (1.1) is important. Including the isospin breaking contributions, recently calculated in [16] and the NNLO QCD corrections to electroweak penguin contributions [17], the result in (1.1) is changed to [16,18]  which compared with the experimental world average from NA48 [20] and KTeV [21,22] collaborations, shows a very good agreement of the SM with the data, albeit leaving still much room for NP contributions. Presently values as low as 5 or as high as 25 in these units cannot be excluded. While this result allows for both positive and negative NP contributions to ε /ε to agree with the data, the new SM value for the K 0 −K 0 mass difference ∆M K from RBC-UKQCD [6] (∆M K ) SM = 7.7(2.1) × 10 −15 GeV, (∆M K ) exp = 3.484(6) × 10 −15 GeV , (1.4) hints at the 2σ level at the presence of NP required to suppress ∆M K relative to its SM value.
As noted already in [13] the suppression of ∆M K is only possible in the presence of new CP-violating couplings. This could appear surprising at first sight, since ∆M K is a CP-conserving quantity but simply follows from the fact that the BSM shift (∆M K ) BSM is proportional to the square of a complex g sd coupling so that The required negative contribution implies automatically NP contributions to ε /ε and also to rare decays K → πνν, K S → µ + µ − and K L → π 0 + − , provided this NP involves JHEP12(2020)097 non-vanishing flavour conserving qq couplings in the case of ε /ε and non-vanishing νν and µ + µ − couplings in the case of the rare K decays in question. But as pointed out in an important paper by Monika Blanke eleven years ago [23], in the presence of a strict correlation between NP contributions to ∆S = 1 and ∆S = 2 processes and in the absence of significant NP contributions to ε K implies two allowed narrow branches in the (K + → π + νν, K L → π 0 νν) plane to be called MB-branches in what follows, and thus a well defined pattern of the correlations between the branching ratios for these two decays. We will be more explicit about these findings in section 4. But as in [23] a possible impact of ∆M K on these correlations has not been analysed we will investigate the impact of the hinted ∆M K anomaly on the findings of [23]. Now, The most recent result for K + → π + νν from NA62 [24,25] and the 90% confidence level (CL) upper bound on K L → π 0 νν from KOTO [26] read respectively B(K + → π + νν) exp = (11.0 +4.0 −3.5 ± 0.3) × 10 −11 , B(K L → π 0 νν) exp ≤ 3.0 × 10 −9 . (1. 6) to be compared with the SM predictions [27,28] B(K + → π + νν) SM = (8. In their most recent status report [29] on K L → π 0 νν the KOTO collaboration presented data on four candidate events in the signal region, finding B(K L → π 0 νν) KOTO = 2.1 +2.0(+4.1) −1.1(−1.7) × 10 −9 , (1.8) at the 68 (95) % CL. The central value is by a factor of 65 above the central SM prediction and in fact violates the GN bound which at the 90% CL together with the present NA62 result for K + → π + νν amounts to 0.8 × 10 −9 . Theoretical analyses of this interesting data can be found in [30][31][32][33].
Evidently there is still much room for NP left in these decays. In particular, a novel pattern in which K + → π + νν is suppressed and K L → π 0 νν is enhanced by NP is allowed by the new data. As pointed out already in [23] and seen in the plots in [23,34] this pattern is only possible in the presence of both left-handed and right-handed flavour-violating Z couplings to quarks which with moderate fine-tuning allows to avoid the constraint from ε K , so that regions in the (K + → π + νν, K L → π 0 νν) plane outside the MB-Branches are possible. We will return to this issue in section 4.5, but we stress already here, following [23], that generally in NP scenarios in which NP contributions to ∆S = 1 and ∆S = 2 are not related to each other, different oases in the (K + → π + νν, K L → π 0 νν) plane outside the MB-Branches could be occupied. As evident from the plots in [34,35] the simplest example are models with minimal flavour violation (MFV). There the correlation between K + → π + νν and K L → π 0 νν results from the same real valued loop function X entering these two processes. This function is a priori unrelated to NP contributions in ∆S = 2 processes and therefore ∆S = 2 constraints are avoided. On the other hand in the absence of new complex flavour-violating phases in MFV models the suppression of ∆M K is not possible. This is reminiscent of lower bounds on ∆M s,d present in these models [36,37].

JHEP12(2020)097
It has been pointed out already in [13] that various patterns of NP in rare K decays in correlation with NP in ε /ε can naturally be realized in models with tree-level FCNCs mediated by a heavy Z with masses still in the reach of ATLAS and CMS but also for higher masses. But whereas in [13] the scenarios with enhanced K + → π + νν and K L → π 0 νν have been primarily considered, a novel pattern in which K + → π + νν is suppressed and K L → π 0 νν is enhanced by NP still allowed by the new data has not been considered there.
With the new information from RBC-UKQCD on ε /ε and ∆M K , the new analyses of ε /ε in [16,18] and the new data from NA62 and KOTO, it is of interest to ask how the Z scenarios considered in [13] and the new ones face the new developments listed above.
The goal of the present paper is to answer this question, but our paper should not be considered as the numerical update of the analysis in [13] motivated by the new input from RBC-UKQCD, NA62 and KOTO collaborations. The reason is that in contrast to [13], that included only QCD renormalization group effects, we will perform a complete SMEFT analysis, that takes in particular into account important top Yukawa effects, which modify significantly the properties of a Z responsible for the pattern of NP effects in question.
In particular we point out that the so-called QCD penguin scenario for ε /ε, considered in [13], in which at the NP scale only QCD penguin operators have non-vanishing Wilson coefficients, is excluded due to Yukawa renormalization group effects on ε K when NP contributions to ε /ε and ∆M K are considered simultaneously. We demonstrate this effect both analytically and numerically.
In models with vector-like quarks the operators ψ 2 H 2 D, listed in table 5, are generated at the matching scale, implying FCNCs mediated by the SM Z-boson. They can be enhanced through RG Yukawa top quark effects with an important impact on the phenomenology [28,38,39]. Such operators have vanishing Wilson coefficients in Z models at tree-level if the (H † iD µ H)Z µ coupling is set to zero. However, they are generated again through RG Yukawa top quark effects. To our knowledge this mechanism of generating FCNCs mediated by the Z in Z models has not been considered in the literature. Usually the FCNCs in Z models are generated through Z − Z mixing in the process of the spontaneous breakdown of the electroweak symmetry [40]. It is then of interest to investigate whether this pure RG effect is important.
Our paper is organized as follows. In section 2 we recall the strategy of [13] where the correlations between ε /ε, ∆M K , K + → π + νν and K L → π 0 νν have been analyzed in the framework of Z models taking into account the constraints from K L → µ + µ − and ε K . We refrain, with the exception of ε /ε, from listing the formulae for observables entering our analysis as they can be found in [13] and in more general papers on Z models in [41] and in [42] that deals with 331 models. On the other hand we discuss in some detail the aspects of new dynamics that enrich the analysis of [13] through the inclusion of the full machinery of the SMEFT, in particular of the renormalization group effects from top Yukawa coupling.
In section 3 as a preparation for the numerical analysis we discuss various Z scenarios and the related RG evolution patterns in the SMEFT.
In section 4 after recalling the arguments for the two branch structure in the (K + → π + νν, K L → π 0 νν) plane pointed out in [23], we present a detailed numerical analysis of JHEP12(2020)097 all observables listed above, including also K S → µ + µ − and K L → π 0 + − , in various Z scenarios.
In section 5 we analyze the generation of FCNCs mediated by the SM Z-boson. In section 6 we list the main results of our paper and present a brief outlook for the coming years. Some additional information is contained in an appendix.

Strategy
In our paper, as in [13], an important role will be played by ε /ε and ε K for which in the presence of NP contributions, to be called BSM in what follows, we have In view of uncertainties present still in the SM estimates of ε /ε, and to a lesser extent in ε K , we will fully concentrate on BSM contributions. Therefore in order to identify the pattern of BSM contributions to flavour observables implied by allowed BSM contributions to ε /ε in a transparent manner, we will proceed in a given Z scenario as follows [13]: Step 1. We assume that BSM provides a shift in ε /ε: with the range for κ ε indicating conservatively the room left for BSM contributions. This range is dictated by the recent analyses in [16,18] which implies the result quoted in (1.2). Specifically, we will consider three ranges for κ ε Only range A has been considered in [13] so that the study of ranges B and C is new with interesting consequences. This step will determine for given flavour conserving Z qq couplings the imaginary parts of flavour-violating Z couplings to quarks as functions of κ ε . But as we will see below in order to explain the ∆M K anomaly, which requires significant imaginary couplings, and simultaneously obtain ε /ε consistent with the ranges above the flavour conserving Z qq couplings must be O(10 −2 ).
We stress even stronger the usefulness of κ ε in the 2020s than it could be anticipated in [13]. The result in (1.2) governed by the hadronic matrix elements from the RBC-UKQCD collaboration has a very large error and we expect that it will still take some time before this error will be decreased down to 10-15%. In addition we need a second lattice group to confirm the 2020 RBC-UKQCD value and it is not evident that this will happen in this decade.

JHEP12(2020)097
Step 2. In order to determine the relevant real parts of the couplings involved, in the presence of the imaginary part determined from ε /ε, we will assume that BSM can also affect the parameter ε K . We will describe this effect by the parameter κ ε so that now in addition to (2.2) we will allow for a BSM shift in ε K in the range This is consistent with present analyses in [43][44][45]. But it should be stressed that this depends on whether inclusive or exclusive determinations of |V ub | and |V cb | are used and with the inclusive ones the SM value of ε K agrees well with the data. We will also investigate how our results change when a larger NP contribution to ε K corresponding to −0.5 ≤ κ ε ≤ 0.5 is admitted.
Step 3. As far as ∆M K is concerned, we will consider dominantly NP parameters which provide the suppression of the SM value in accordance with the LQCD result in (1.4). In particular this will require the imaginary Z couplings to be significantly larger than the real ones.
Step 4. In view of the uncertainty in κ ε we set several parameters to their central values.
In particular for the SM contributions to rare decays we set the CKM factors and the CKM phase δ to which are close to the central values of present estimates obtained by the UTfit [43] and CKMfitter [44] collaborations. For this choice of CKM parameters the central value of the resulting |ε SM K | is 2.32 · 10 −3 . With the experimental value of ε K in table 3 this implies κ ε = −0.09 . But we will still vary κ ε while keeping the values in (2.5), as BSM contributions in our scenarios do not depend on them but are sensitive functions of κ ε .
Step 5. Having fixed the flavour violating couplings of the Z in this manner, we will be able to calculate BSM contributions to the branching ratios for K + → π + νν, K L → π 0 νν, K L,S → µ + µ − and K L → π 0 + − and to ∆M K in terms of κ ε and κ ε . This will allow us to study directly the impact of possible NP contributions to ε /ε and ∆M K in Z scenarios on K + → π + νν and K L → π 0 νν and the remaining rare Kaon decays. In table 1 we indicate the dependence of a given observable on the real and/or imaginary Z or later Z flavour violating coupling to quarks. In our strategy imaginary parts depend only on κ ε and the choice of flavour conserving Z qq couplings, while the real parts depend on both κ ε and κ ε . The pattern of flavour violation depends in a given BSM scenario on the relative size of the real and imaginary parts of the couplings as we will see explicitly later on.
In the context of our presentation we will see that in most of our Z scenarios ε K and not K L → µ + µ − is the most important observable for the determination of the real parts of the new couplings after the ε /ε constraint has been imposed. This can be traced back to Yukawa RG effects. Additional constraint will come from ∆M K .

SMEFT at work
The interaction Lagrangian of a Z = (1, 1) 0 field and the SM fermions reads: Here q i and i denote left-handed SU(2) L doublets and u i , d i and e i are right-handed singlets.
This Z theory will then be matched at the scale M Z onto the SMEFT, generating the operators listed in table 2. In the Warsaw basis [46] the tree-level matching [47] with the couplings in (2.6) is given for purely left-handed vector operators by: For purely right-handed vector operators one finds: (2.11) JHEP12(2020)097 Finally for left-right vector operators the matching reads: (2.14) Different bases 2 for the SMEFT Wilson coefficients (corresponding to different models) can be used to perform the numerical analysis. A particular choice of basis is the downbasis, 3 in which the down-type Yukawas are diagonal and the q i fields are given above the EW scale by where V denotes the CKM matrix. Another popular basis choice is the up-basis with diagonal up-type Yukawas and (up-basis) (2.16) Changing between these two bases is achieved by rotating the corresponding parameters by CKM factors. For instance, to express the up-basis g ij q couplings in terms of the down-basis JHEP12(2020)097 ones, the following rotation needs to be performed: Since we are interested in FCNCs in the down-sector, it is more convenient to work in the down-basis, which we will adopt in the following. In a next step the SMEFT Wilson coefficients are evolved from the matching scale Λ down to the EW scale µ EW . In order to perform this RG evolution the SM parameters are first run up to the high scale Λ, such that all input parameters (Wilson coefficients and SM parameters) are evolved from the same scale down to µ EW . The procedure to obtain the SM parameters at the high scale is discussed in the next subsection.

Treatment of SM parameters
In order to solve the RGEs, assuming experimental values of the SM parameters at the EW scale, we evolve them to the input scale Λ. For this purpose we employ an iterative procedure, which was used in [50]. This procedure for solving the RGEs incorporates the correct values of the CKM parameters and the quark and lepton masses at the electroweak scale. Since in the present paper we are interested in exploring the role of Yukawa RGE effects, let us describe the iterative procedure to determine the Yukawa couplings at the input scale Λ: • We start with the Yukawa matrices in the down-basis at the EW scale: with the mass matrices given by Here the values of the quark masses can be found in table 1 of [50].
• In the first step the Yukawa matrices are evolved up to the input scale Λ while assuming constant Wilson coefficients (equal to their input values, C i (Λ)). As the chosen basis is not stable under RG running, a rotation of the fermion fields is performed to get back to the down-basis, Here the unprimed fields are in the down-basis, whereas the primed fields are in some random basis generated by the running of Yukawas just performed. The rotation JHEP12(2020)097 matrices at the input scale transform the primed mass matrices back to the down-basis obtaining the diagonal matrices M d (Λ) and M e (Λ) and the non-diagonal matrix M u (Λ) given in (2.19). The primed matrices at the input scale are given by • In the second step the Wilson coefficients are evolved down to the EW scale in the leading log (LL) approximation, using the Yukawa matrices from the previous step.
• Finally, the obtained Yukawas are evolved up to the input scale Λ using the constant Wilson coefficients obtained from the LL running.
This iterative procedure known as the shooting-method [51] is a common method to obtain numerically the initial conditions for ordinary differential equations at the same integration point, and is used to find the Yukawa matrices (and other SM parameters) at the high scale. The procedure is iterated until a precision of 0.1% is obtained.
The RGEs can then be solved with all parameters having their initial conditions at the same scale Λ. We reemphasize that the form of the Yukawa matrices is not stable under RGEs and therefore a back-rotation [52] is required to go back to the down basis at the EW scale. A crucial consequence of this is that one also needs to back-rotate [50,[53][54][55][56][57][58] the Wilson coefficients according to table 4 of [59]. We will return to one of these consequences in section 4.4.

Scenarios
For the numerical analysis we follow closely the reasoning in [13]. As we are interested in Kaon decays, we will assume different scenarios for the flavour transition d → s, to be referred to as LHS and RHS in the following. In these scenarios we allow for a flavourviolating coupling in the left-handed (LHS) or right-handed (RHS) quark sector between the second and first generation, respectively. Moreover, we choose the flavour-diagonal first generation quark couplings of both chiralities to be non-vanishing in both scenarios. Due to simultaneous presence of non-zero leptonic and first generation quark couplings the LHC direct searches in the dilepton final state put constraint M Z 5.5 TeV [60,61] for JHEP12(2020)097 the couplings strength same as the SM Z-boson [62]. However, we will consider the Z couplings to the first generation to be an order of magnitude smaller than those of the SM Z-boson which allows us to consider the lowest Z mass equal to 3 TeV. With this choice of couplings VLL, VRR as well as VLR operators given in (2.7)-(2.14) are generated in both scenarios. Furthermore, we define the LR scenario first discussed in [63,64], which is equivalent to the LHS or RHS, but without taking into account constraints from ε K and ∆M K . The justification for this procedure is given as follows. In the LR scenario containing LH as well as RH Z couplings to SM fermions, left-right ∆F = 2 operators are generated at tree-level. Their contributions to the mixing amplitudes M bs 12 , M bd 12 and M sd 12 are RG enhanced. For the d → s transition there is an additional chiral enhancement. However, by imposing a fine-tuning between the left-left and right-right contributions and the LR contributions, the constraints from ε K and ∆M K can be alleviated while giving sizable contributions to K → ππ, as has been shown in [63]. The generalization of this idea to the other meson systems has been done in [64]. We will briefly return to this scenario in section 4.5. For further details we refer to appendix A of [65].
In the LHS, the flavour change is achieved by the non-zero complex coupling g 21 q and in the RHS by a complex-valued g 21 d . In each scenario we allow for diagonal (real) couplings to first generation quarks. Furthermore, to also accommodate for the decays K L,S → µ + µ − a real non-zero value of g 22 is chosen. For K L → π 0 + − with = e, µ we also need non-zero g 11 and g 22 . All other couplings are assumed to vanish. Therefore we have at the high scale Λ the following three setups: Such scenarios are in general subject to gauge anomalies, which are assumed to be canceled by additional heavy fields at a higher scale [40,66,67]. Z models with explicit gauge anomaly cancellation were discussed recently in [65,68,69].
Using the matching relations in section 2.2 leads to the following non-zero four-fermion Wilson coefficients in the three different scenarios at the BSM scale: qd 1111 , C   Figure 1. Running of four-fermion operators into operators contributing to ∆F = 1 non-leptonic observables. Here the red, green and black lines indicate the operator mixing due to strong, weak and Yukawa couplings respectively. The self-mixing for all couplings is shown by a dashed black line.
operators contributing to non-leptonic ∆S = 1 and ∆S = 2 observables and semi-leptonic ∆S = 1 decays. We present the charts for LHS and RHS, where LH and RH refer to flavour-violating currents. The structure of these charts is as follows: • At the BSM scale those operators are listed which on the one hand receive a nonvanishing matching contribution and on the other hand imply through RG evolution contributions at the electroweak scale. The latter can come from the same operators with modified Wilson coefficients and from new operators generated through RG evolution. These new operators are placed on a lighter background than the original operators.
• As an example consider the first chart in figure 1. The goal is to generate at the electroweak scale four-quark operators contributing to non-leptonic ∆S = 1 processes which is indicated by the indices (2111) . The operators present already at the BSM scale contribute also at the electroweak scale but whereas the indices of the Wilson coefficients C    Figure 2. Running of four-fermion operators into operators contributing to ∆F = 2 observables.
Here the red, green and black lines indicate the operator mixing due to strong, weak and Yukawa couplings respectively. The self-mixing for all couplings is shown by a dashed black line. Figure 3. Running of four-fermion operators into operators contributing to ∆F = 1 semi-leptonic observables. Here the red, green and black lines indicate the operator mixing due to strong, weak and Yukawa couplings respectively. The self-mixing for all couplings is shown by a dashed black line.

JHEP12(2020)097
scale, while not contributing directly to ∆S = 1 non-leptonic observables, can do it indirectly via Wilson coefficients of non-leptonic operators through electroweak interactions.
• The same logic is used in the remaining charts. But one should note that in the RHS the flavour-violating indices are on the right-handed currents so that e.g. on the top of the lower charts in figures 1-3 the indices are now (1121) instead of (2111).
• The distinction between strong, weak and Yukawa interactions is made with the help of colours as described in the figure caption.

ε /ε
Since ε /ε is one of the key observables in our analysis we discuss here explicitly the impact of the LHS and RHS on this observable. The relevant SMEFT matching contributions for ε /ε can be found in [73]. Adopting the same short distance basis as therein, namely with colour indices i, j, chiralities A, B = L, R, and Dirac structures X = S, V, T with Γ S = 1, Γ V = γ µ , Γ T = σ µν , one finds at the high scale Λ: where we have neglected small contributions. However, as indicated by the red arrows in figure 1, the Wilson coefficients C qu 2111 and C (8) qd 2111 ( C qd 1121 ) are induced through QCD running down to the EW scale in the LHS (RHS). At LL one finds [72,74]: and similar expressions for C (8) qd 2111 and C (8) qd 1121 . Therefore, taking QCD RGE effects into account the matching at the BSM scale in (3.9)-(3.10) is modified at the EW scale as follows qd 2111 , (3.12) qd 1121 , (3.14) qd 1121 ,

Rare decays
For completeness we recall general ∆F = 1 Hamiltonian for the semi-leptonic FCNC transition of down-type quarks into leptons and neutrinos below µ EW with a, b being lepton indices, i, j down-quark indices and There are eight semi-leptonic operators relevant for d i a → d j b when considering UV completions that give rise to SMEFT above the electroweak scale [77] Q baji and two for The SM contribution to these Wilson coefficients is lepton-flavour diagonal where v = 246 GeV and a normalisation factor has been introduced for the NP contribution that proves convenient for matching with SMEFT in a given model. The non-vanishing SM contributions

JHEP12(2020)097
are given by the gauge-independent functions X 0 (x t ), Y 0 (x t ) and Z 0 (x t ) [78] that can be found in [2]. The NP contribution to the Wilson coefficients of the ∆F = 1 semi-leptonic operators in (3.20) and (3.21) at µ EW in terms of the semi-leptonic SMEFT Wilson coefficients at µ EW is given as follows [49,77,79] C baji 9,NP = C qe + C (3.24) We caution the reader that C qe , used in the literature, should here be written as C eq so that the lepton indices come first as in the remaining WCs in these equations. But as we already stated after (3.18) a, b are lepton indices and i, j down-quark ones, so that no confusion should result from this notation. Here contributions from Z-mediating ψ 2 H 2 D- Hq to C 9,10,L and O Hd to C 9 ,10 ,R , respectively, have been omitted. In rare FCNC Kaon decays scalar and pseudo-scalar Wilson coefficients are negligible but are relevant in B s → µ + µ − .

Observables
In our numerical analysis we investigate the following quantities: For the numerical analysis the input parameters in tables 3 and 4 are used. The constraint from B(K L → µ + µ − ) at the 2σ level is taken into account. The SM predictions for K + → π + νν and K L → π 0 νν are given in (1.7) and for the remaining decays one finds [28,[80][81][82]: where for the K L → π 0 + − decays the numbers in parenthesis denote the destructive interference case. The experimental status of these decays is given by [83][84][85]:    Finally, for the LHS and RHS we impose the constraint from ε K in the following way: where κ ε is defined in (2.4). But we will investigate what happens for a larger range κ ε ∈ [−0.5, 0.5].
Concerning uncertainties, the main error in ε /ε from QCD penguins in the SM contribution is incorporated in the parameter κ ε , which we vary in a large range. Uncertainties for the BSM contribution to ε /ε are dominated by the A 2 amplitude with a much smaller error estimated in [76] to be in the ballpark of 10%. Similar comments apply to the CKM parameters. Consequently, we stress that except for κ ε and R ∆M K , which are both varied in a large range in our analysis, the general pattern only weakly depends on theoretical uncertainties and we refrain from showing them explicitly in what follows.
4.2 Impact of ε K and ∆M K on K + → π + νν and K L → π 0 νν As pointed out in [23], in the presence of a strict correlation between NP contributions to ∆S = 1 and ∆S = 2 processes and assuming no significant NP contributions to ε K implies two narrow branches in the (K + → π + νν, K L → π 0 νν) plane, namely • a branch parallel to the Grossman-Nir (GN) bound [86] on which both branching ratios can either simultaneously increase or decrease relative to SM values, • a horizontal narrow branch on which there is no NP contribution to K L → π 0 νν because of the absence of flavour-violating complex couplings.

JHEP12(2020)097
An example of this pattern is seen in the left plot in figure 4 to which we will return soon. This is in particular the case of NP entering already at tree-level with only lefthanded or right-handed flavour-violating NP couplings with the prominent example of Z models in which the Z sd coupling enters both K + → π + νν and K L → π 0 νν as well as ε K .
Here we would like to point out that the hinted anomaly in ∆M K , requiring the imaginary couplings to be present, excludes in such a scenario the horizontal branch so that the full action of NP in this case happens only on the second MB-branch that is parallel to the GN bound.
It should be stressed that this is a new insight in this pattern beyond the one found in [23] where a possible impact of ∆M K has not been discussed. There, under the assumption of significant NP contributions to ε K but still considering only scenarios with lefthanded or right-handed flavour-violating Z couplings significantly broader branches have been found than when ε K from the SM agreed with the data. However, as we will demonstrate below, the removal of the ∆M K anomaly combined with renormalization group Yukawa top effects implies still a rather narrow MB-branch parallel to the GN bound.
Furthermore it should be noted, that we perform the following analysis with the intention to investigate possible correlations of ε /ε and ∆M K with rare Kaon decays, while considering significant NP contributions to the former ones. In this article we refrain from performing parameter fits to the observables in question, but rather fix the couplings of a given scenario to exhibit such correlations. A more rigorous analysis involving a fitting procedure might be appropriate once the experimental and theoretical precision are comparable, which should be the case in the upcoming years.

Electroweak penguin scenario: left-handed
We start with a LHS (i.e. g 21 q = 0), where the effect in ε /ε is achieved through electroweak penguin (EWP) operators such as Q 8 . To generate such operators we choose the quark couplings in the following way: In figure 4 (left), we plot the correlation between the ratios for the decays K + → π + νν and K L → πνν. Here the horizontal and vertical branches correspond to purely real and imaginary values respectively of the flavour violating coupling g 21 q . Simultaneous presence of both real and imaginary parts, which correspond to the small area at the meeting point of the two branches, are strongly constrained by the allowed range of κ ε (4.4). Furthermore, requiring the suppression of ∆M K excludes the horizontal branch, indicating the dominance of the imaginary part over the real part of g 21 q . This implies a strong correlation between K + → π + νν and K L → πνν on the MBbranch, so that they can be enhanced or suppressed only simultaneously as shown by the orange colour in this figure. Out of the three κ ε scenarios A, B and C, which are defined in (2.3), in scenario A, large departures from SM expectations for K L → π 0 νν are possible. Similarly, in figure 4 (right), the correlations between the ratio for the decay K L → πνν and the ones for K → π + − and K S → µ + µ − are shown. The upper range for R 0 νν corresponds JHEP12(2020)097   Figure 4. LH-EWP scenario for a Z of 3 TeV. The correlation between the ratios for the process K + → π + νν, K L → πνν defined in (4.1) is plotted (left). The blue (orange) lines are allowed by κ ε (∆M K ) constraints and the black line represents the GN bound. The correlations between the ratio for K L → π 0 νν and the ones for K → π + − and K S → µ + µ − are shown (right). roughly to the GN bound. If the values from KOTO given in (1.8) will be confirmed in the future, large departures from the SM predictions for the three rare decays are to be expected. Also the K + → π + νν branching ratio could be enhanced. Figure 4 (right) admits two solutions for each decay, corresponding to different values of κ ε . The upper branch results from positive values for Im(g 21 q ) and the lower one from negative ones, since positive (negative) values of Im(g 21 q ) enhance (reduce) the corresponding ratios. In figure 5 we show the results for the first three different ratios defined in (4.1) as functions of κ ε for a Z of 3 TeV and 10 TeV respectively. For the running below the EW scale we use the complete 1-loop QCD and QED running [87,88] and above the EW JHEP12(2020)097  Figure 6. LH-EWP scenario for a Z of 3 TeV (left panel) and 10 TeV (right panel). The predictions for the ratios of the decays K S → µ + µ − , K L → π 0 µ + µ − and K L → π 0 e + e − defined in eq. (4.1) are plotted against κ ε . The yellow, green and red shades correspond to the κ ε scenarios A, B and C as defined in (2.3).
scale the full SMEFT RGEs for the solid and only QCD for the dashed lines are used. Clearly, the running is dominated by QCD effects. For 3 TeV both K → πνν branching ratios are enhanced over their SM values, except for a small region around κ ε ≈ 0. For 10 TeV, significant BSM effects are only observed for κ ε ≥ 0.5. ∆M K is visibly suppressed for sufficiently large κ ε . The choice of very small values of g 11 u,d of O(10 −2 ) is implied, as noticed already in [13], by the desire to suppress ∆M K in the presence of NP contributions to ε /ε in the EWP sector. For g 11 u,d of O(1) considered in the latter paper, ∆M K is enhanced by BSM rather than suppressed which is disfavored by the present LQCD data. In figure 6 we show predictions for the remaining ratios given in (4.1), where we allow for additional couplings to left-handed electrons (g 11 ). We observe that for a lighter Z an enhancement for R 0 π 0 µµ and R 0 π 0 ee processes is predicted for negative values of κ ε , while for its positive values both suppression as well as enhancement are possible. On the other hand for heavier Z these decay modes are suppressed (enhanced) for positive (negative) values of κ ε . The ratio R S µµ is always enhanced. The difference between solid and dashed lines is mainly due to QED RG effects on κ ε , generated by semi-leptonic operators.
In figure 7 the correlations between κ ε and R ∆M K and between the ratios for K + → π + νν and K L → π 0 νν and R ∆M K are shown. As expected, κ ε and K L → π 0 νν are much more sensitive to variations of R ∆M K than it is the case of K + → π + νν.  Figure 8. LH-EWP scenario for a Z of 3 TeV (left panel) and 10 TeV (right panel). The predictions for the ratios of the decays K S → µ + µ − , K L → π 0 µ + µ − and K L → π 0 e + e − defined in (4.1) are plotted against R ∆M K .

QCD penguin scenario: left-and right-handed
Next we describe the effects related to the required basis rerotation at the electroweak scale, as described in the last point of section 2.3. This has important phenomenological consequences in any scenario, as for example in the QCD penguin (QCDP) scenario, in which a sizable imaginary coupling is present in scenarios A and C for κ ε . The LH-QCDP scenario is defined as follows: Starting with a set of non-zero Wilson coefficients in the down-basis at the high scale Λ we evolve them to the EW scale. Along with the Wilson coefficients we also need to evolve the SM parameters including the mass (or Yukawa) matrices as discussed in section 2.3. But the JHEP12(2020)097 running of the mass matrices is flavour dependent [71], and consequently after the evolution the mass matrices are not guaranteed to remain in the original basis that we started with. As a result, we need to rotate the mass matrices and hence the Wilson coefficients to adhere to our choice of the down-basis [50]. This issue is discussed in generality in a recent paper [52] but here we confine our discussion focusing on QCDP.
We illustrate this effect and its phenomenological consequences with a concrete example by considering the LH-QCDP scenario studied in the case of significant BSM contributions to ε /ε in [13], but now in contrast to that paper including RG SMEFT effects. Considering the LHS, at the high scale Λ the operators [O (1) qu ] 2111 and [O (1) qd ] 2111 are generated. They are then evolved down to the EW scale. But the simultaneous evolution of the mass matrices generates off-diagonal entries in the down-quark Yukawa matrix Y d at the EW scale. This is due to the fact that the running of Y d is proportional to the up-quark Yukawa matrix Y u , which is non-diagonal in the down-basis [89]. Indeed, we have To revert to the down-type basis, a rotation of the operators is necessary, as already explained in section 2.3. Applying this back-rotation to the Wilson coefficients generates C (1) qd 2121 at the EW scale in the down-basis as: where C (1) qd 2111 denotes the Wilson coefficient in the RGE basis and the rotation matrices U d L , U d R satisfy the following equation: Here the (non-diagonal) down-quark mass matrix M d at the EW scale is obtained by evolving Y d from the high scale Λ down to µ EW . In the LL approximation we have: However, the Wilson coefficient C qd 2121 is strongly constrained by ε K due to the large hadronic matrix element multiplying it. This scenario is illustrated in figure 9, where at the high scale we vary the input values of the Wilson coefficients C as shown on the x-axis. On the y-axis we show the output value of the Wilson coefficient C (1) qd 2121 at the EW scale which is generated through the back-rotation of (4.8). However, this LR operator gives a large contribution to ε K [13] κ ε 3.1 × 10 16 × Im C (1) qd 2121 (µ ew )GeV −2 . (4.11) The allowed values for the Wilson coefficients of the three mentioned operators are shown in the red region, given the constraints from ε K . This shows that in the LHS with QCDP dominance, significant BSM contributions to ε /ε imply a large contribution to ε K JHEP12(2020)097 qd ] 2121 is generated after RGE running of Y d and back-rotation to the down-basis at the EW scale. The allowed regions for the Wilson coefficients are in red for ε K , and vertical bands represent the three ε /ε scenarios.
inevitably generated by the running of Yukawas and subsequent back-rotation of the Wilson coefficients at the EW scale. Consequently, the QCDP scenario for ε /ε, considered in [13] is ruled out, since in this case significant BSM contributions to ε /ε would be required to fit the data. Similar comments apply to the RHS scenario defined by In this case only the QCDP scenario can be constructed. Due to SU(2) L gauge invariance the coefficient of the so-called Q 8 operator, which otherwise would give a leading contribution to ε /ε, vanishes.
On the other hand, in the case of the EWP dominance i.e. C (1) qu 2111 , also considered in [13], this effect is negligible. This is simply because in this case a much smaller value of C (1) qd 2111 is needed to enhance sufficiently ε /ε.

Left-right scenario
We have just seen that in the LHS there was a very strong correlation between K + → π + νν and K L → π 0 νν branching ratios on the MB-branch. As explained in [23] this strict correlation originates in the same complex phase present in NP contributions to ε K and rare Kaon decays in question provided NP contributions to ε K are small. This is in fact evident in our case because the same Z sd coupling enters both K + → π + νν and K L → π 0 νν and ε K . Now, 13) and to make sure that this contribution is small either Re(g sd ) or Im(g sd ) must be small. If Im(g sd ) is small the horizontal line in figure 4 results with NP basically only in K + → π + νν. If Re(g sd ) is small then there are NP contributions to both K L → π 0 νν and K + → π + νν correlated on the MB-branch. In our case this second solution is chosen by the desire to explain the ∆M K anomaly. However, such a correlation precludes the pattern of simultaneously enhancing K L → π 0 νν and suppressing K + → π + νν possibly still allowed by the NA62 and KOTO results. It is known from various studies that such a pattern can be obtained through the introduction of new operators and the most effective in this respect are scenarios in which both left-handed and right-handed flavour-violating NP couplings are present, breaking the correlation between K 0 −K 0 mixing and rare Kaon decays and thereby eliminating the impact of the ε K constraint on rare Kaon decays. The presence of left-right operators requires some fine-tuning of the parameters in order to satisfy the ε K constraint but such operators do not contribute to rare decays and the presence of new parameters does not affect directly these decays. Examples of such scenarios are Z models with LH and RH couplings considered in [64] and the earlier studies in the context of the general MSSM [90][91][92][93][94] and Randall-Sundrum models [95,96]. See in particular figure 6 in [95] and figure 7 in [64]. Needless to say also the correlations between NP contributions to ∆M K and rare decays are diluted, although the necessity of non-vanishing complex couplings required by the hinted ∆M K anomaly will certainly have some impact on rare Kaon decays.
The Left-Right (LR) scenario at 3 TeV is defined by which is equivalent to the LH-EWP scenario without imposing ∆F = 2 constraints [64]. In figure 10 correlations between ratios for K L → π 0 νν and K + → π + νν as in (4.1) are considered. Clearly no strong correlation is observed when both LH and RH couplings are allowed as shown in the green region. Where the different curves correspond to the different absolute values of the coupling g 21 q and the area inside the curves is allowed. Similarly, the strong correlation between K + → π + νν and κ ε observed in the LH-EWP scenario is absent in the LR scenario because R + νν also depends on the real part, which is not fixed through ε /ε.

JHEP12(2020)097
Imposing however the constraint from ε K and therefore studying a LH-EWP scenario limits the allowed parameter space drastically. Furthermore, as shown in section 4.3 out of the two branches in the R + νν -R 0 νν plane allowed by ε K , the horizontal branch shown in blue is disfavored by the requirement of suppression of ∆M K . In the red area we show the allowed region for the LH-EWP scenario with ε K ∈ [−0.5, 0.5].
Importantly, as evident from figure 10, the simultaneous enhancement of K L → π 0 νν and suppression of K + → π + νν branching ratios is only possible in the presence of both LH and RH flavour-violating couplings. Also, the observables R 0 νν and κ ε only depend on the imaginary part of the flavour violating coupling. Therefore they are strongly correlated in the LR as well as in the LHS scenario.
This agrees with the findings in [34], in which only QCD has been considered. The correlation between R + νν and R 0 νν in this setup is therefore invariant under Yukawa running effects.

Preliminaries
In this section we consider flavour violating (FV) Z couplings induced by FV Z couplings through SMEFT RG running effects. Let us consider the LL running from the BSM scale Λ to the EW scale µ EW . For the Wilson coefficients of the ψ 2 H 2 D operators defined in table 5 keeping only the top Yukawa coupling y t and neglecting the terms of O(V 2 ts ) and O(V tb V ts ) one finds [71,74] C (1) whereas O H and O Hud are not generated in this approximation. Yukawa running effects therefore generate modified Z-couplings to the SM fermions.
We can now express the usual FC quark couplings of the Z in terms of C (1,3) Hq , C Hu and C Hd . We have first JHEP12(2020)097 with ψ = u, d distinguishing between up-and down-quark couplings. These complex-valued couplings are related to the SMEFT Wilson coefficients through [38] [ where g Z = g 2 1 + g 2 2 and v = 246 GeV is the electroweak vacuum expectation value. In the Z scenario considered here the ψ 2 H 2 D operators are generated through RG effects and are smaller than in the case where these operators are already present at the high scale [38,39,47]. For the time being we assume that this is not the case here but we will comment briefly on their possible impact on our analysis below.

Impact of RG-induced Z on LH-EWP scenario
In this subsection we study an explicit example of FV Z couplings induced by FV Z couplings through SMEFT RG running effects and its effect on the ratios in (4.1). For this purpose we assume two scenarios: in the first one only direct contributions from a Z are generated at the matching scale. This corresponds to the LH-EWP setup in subsection 4.3.
In the second one we allow for additional non-zero couplings to the third generation quarks. The up-type quark coupling will then generate through (5.1) modified Z-couplings, which induce an additional effect compared to the Z -only case. We choose the various couplings at the matching scale as follows: Since the usual SM Z couplings obey the relation the operators Q 7 and Q 8 are generated through matching and QCD running, respectively. The Z contributions to ε /ε generated from a Z via RGE running are therefore of the EWP type.
The results for the above two scenarios are shown in figure 11, where for a Z of 3 TeV the same values for the couplings as in figure 5 are assumed. In addition we have for the Z + Z case. The dashed and solid lines correspond to the Z and Z + Z case respectively. The additional contributions due to the modified Z-couplings are destructive to κ ε in this setup, so that a larger value of g 21 q is needed in order to obtain the same value of κ ε in the presence of Z contributions. Therefore, for a given value of κ ε the effect in semi-leptonic decays and ∆M K is enhanced as compared to the Z -solo scenario. By changing the sign of the third-generation couplings, a constructive effect can be achieved for κ ε .
In the left chart of figure 11 R 0 νν and R + νν are enhanced whereas ∆M K is suppressed. The modified Z contributions can have large influence on K L → π 0 νν which is less pronounced for K + → π + νν for moderate values of ε /ε. The effect in ∆M K is also less pronounced since the modified Z coupling enters quadratically. For the predictions of the (semi)-leptonic decays in the right chart in figure 11 the effect of the generated FV Z coupling is significant for larger absolute values of κ ε and predicts enhancements of all considered ratios.

ε /ε and rare decays from RG-induced Z
In our previous discussion we found that in order to have significant BSM contributions to ε /ε within the EWP scenario right-handed flavour diagonal couplings to the first generation JHEP12(2020)097  Figure 12. This figure shows how the Z-contributions to ε /ε and other Kaon observables are generated from a Z with purely left-handed quark couplings through RG running.
quarks are required. However, in this subsection we show that one can also get BSM contributions to ε /ε even from purely left-handed Z couplings. This can happen through top-Yukawa RG running effects. For this purpose we assume a scenario in which at the high scale the diagonal couplings to the first generation quarks vanish and allow for a rather large third generation coupling, namely This choice ensures vanishing of the direct Z contribution to ε /ε through EWPs. In this setup the Wilson coefficient C qq 2133 is generated at the BSM scale, which in turn generates C (1) Hq 21 at the EW scale through top-Yukawa RGEs, as shown in (5.1). This leads to the flavour violating coupling of the Z-boson (5.8) which along with the usual SM Z couplings (5.12), generate the operators Q 7 and Q 8 . This effect is displayed in figure 12. The different ratios of (4.1) are shown as a function of κ ε . A strong suppression of ∆M K and correlation with ε /ε is possible. The large effect in ∆M K is simply due to the sizable value of the flavour violating coupling present at the BSM scale. Except for R 0 νν all other ratios are almost at their SM values and do not depend on κ ε . In LHS or RHS R 0 νν goes down (up) with increased (decreased) κ ε in Z-scenarios. This is because of special values of flavour diagonal Zqq couplings that equal the SM ones in this scenario. See the plots in [13,34].
In a similar fashion with different combinations of Z couplings at the NP scale the Z couplings can be modified through other ψ 2 H 2 D 2 operators given in (5.1)-(5.6).
Finally, it should be emphasized following [38,39] that Z contributions to ε K and ∆M K considered by us correspond really to dimension-eight operators, but the fact that the FV Z couplings in rare decays and Wilson coefficients of these operators are the same implies JHEP12(2020)097 correlations between ∆S = 1 and ∆S = 2 observables [41]. These correlations are strongly modified, even broken, in the presence of non-vanishing Wilson coefficients of ψ 2 H 2 D operators already at the NP scale. Indeed, through top-Yukawa RG effects dimension-six operators contributing to ε K and ∆M K are generated, implying in particular in the case of the O Hd operator strong constraints on rare Kaon decays [28,38,39].

Summary and outlook
The main goal of our paper was to confront Z scenarios with the pattern of BSM contributions hinted by recent results on ε /ε, ∆M K , K + → π + νν and K L → π 0 νν that appear to • allow significant positive or negative BSM contributions to ε /ε relative to its SM value, • suppress the mass difference ∆M K relative to the recent SM value obtained by the RBC-UKQCD collaboration, • suppress the branching ratio for K + → π + νν relative to the precise SM predictions as still allowed by the recent result from the NA62 collaboration, although significant enhancements are still possible, • enhance the branching ratio for K L → π 0 νν relative to the precise SM prediction as hinted by the recent result from the KOTO collaboration.
Taking into account the constraints from ε K and K L → µ + µ − we have calculated ∆M K and the branching ratios for K + → π + νν and K L → π 0 νν as functions of the parameter κ ε introduced in [13] for the choices of Z couplings to quarks and leptons that can reproduce the pattern of deviations from SM expectations summarized above. For these choices of couplings we have calculated the implications for K S → µ + µ − and K L → π 0 + − again as functions of κ ε . Moreover, we have investigated correlations between all these observables in various Z scenarios.
While an analysis of this sort has been already presented in [13], prior to the last three hints for the pattern of BSM contributions, and earlier analyses can be found in [41,42], this is the first analysis of this set of observables to date that took into account RG effects in the framework of the SMEFT, in particular the effects of top Yukawa couplings.
In this context we have also investigated for the first time whether the presence of a heavy Z with flavour violating couplings could generate through top Yukawa renormalization group effects FCNCs mediated by the SM Z-boson. Our results can be found in numerous plots. Here we want to list the most important lessons from our analysis. Lesson 1. While the correlation between the enhancement of ε /ε with the suppression of ∆M K has been already pointed out in the context of the QCD penguin scenario for ε /ε for flavour diagonal Z couplings to quarks of O(1) in [13], we find that the inclusion of RG top quark Yukawa effects rules out this scenario through the ε K constraint.

JHEP12(2020)097
Lesson 2. While, as noticed already in [13], the suppression of ∆M K in the presence of the enhancement of ε /ε could in the EW penguin scenario be only obtained for flavour diagonal Z couplings to quarks of O(10 −2 ), a numerical analysis of such a scenario has not been presented there. Our analysis demonstrates that the expectations from [13] are confirmed in the presence of the full RG SMEFT analysis. In particular the ε K constraints are satisfied.
Lesson 3. We point out that the suppression of K + → π + νν and simultaneous enhancement of K L → π 0 νν would give in the context of Z models some indication for the presence of right-handed flavour violating currents at work. The confirmation of these findings requires in particular a much more accurate measurement of the K + → π + νν branching ratio by NA62. Otherwise a strong correlation between K + → π + νν and K L → π 0 νν branching ratios on the MB-branch is implied by the hinted ∆M K anomaly. In this case if the large enhancement of K L → π 0 νν branching ratios signaled by the KOTO experiment is confirmed one day, also significant enhancement of the K + → π + νν branching ratio over its SM value is to be expected. As seen in figure 4, even larger departures from SM predictions should then be observed in K → π + − and K S → µ + µ − .
Lesson 4. We have demonstrated that RG effects can in the presence of Z contributions generate flavour-violating Z contributions to ε /ε and rare decays that have significant impact on the phenomenology as shown in figure 11. What we also find is that in the presence of O(1) diagonal Z top-quark couplings, the (V − A) × (V + A) EWP operators can be generated solely through the RG induced flavour-violating Z couplings. As shown in figure 12 this effect is sufficiently strong to provide significant BSM contributions to ε /ε, if required, while simultaneously suppressing ∆M K .
Lesson 5. The impact of BSM effects on rare Kaon decays depends both on the scenarios discussed and on the values of the couplings involved. With improved measurements it will be possibly to select the favorite scenarios. In this context the determination of the parameter κ ε through improved LQCD calculations will be important because, as seen in several plots, some of the rare branching ratios depend sensitively on this parameter.
We are looking forward to experimental and theoretical developments in the coming years. Our plots will allow to monitor them and help to identify the successful Z scenarios.  Table 6. Numerical values of K → ππ SM hadronic matrix elements used in our analysis.

A Hadronic matrix elements
In this appendix we report the hadronic matrix elements we use for the numerics of ε /ε, which have been updated recently by the RBC-UKQCD collaboration [14]. They are given in table 6.