\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varepsilon '/\varepsilon $$\end{document}ε′/ε in the Standard Model at the Dawn of the 2020s

We reanalyse the ratio \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varepsilon '/\varepsilon $$\end{document}ε′/ε in the Standard Model (SM) using most recent hadronic matrix elements from the RBC-UKQCD collaboration in combination with most important NNLO QCD corrections to electroweak penguin contributions and the isospin-breaking corrections. We illustrate the importance of the latter by using their latest estimate from chiral perturbation theory (ChPT) based on the octet approximation for lowest-lying mesons and a very recent estimate in the nonet scheme that takes into account the contribution of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\eta _0$$\end{document}η0. We find \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(\varepsilon '/\varepsilon )^{(8)}_\text {SM} = (17.4 \pm 6.1) \times 10^{-4}$$\end{document}(ε′/ε)SM(8)=(17.4±6.1)×10-4 and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(\varepsilon '/\varepsilon )^{(9)}_\text {SM} = (13.9 \pm 5.2) \times 10^{-4}$$\end{document}(ε′/ε)SM(9)=(13.9±5.2)×10-4, respectively. Despite a very good agreement with the measured value \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(\varepsilon '/\varepsilon )_\text {exp} = (16.6 \pm 2.3) \times 10^{-4}$$\end{document}(ε′/ε)exp=(16.6±2.3)×10-4, the large error in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(\varepsilon '/\varepsilon )_\text {SM}$$\end{document}(ε′/ε)SM still leaves room for significant new physics (BSM) contributions to this ratio. We update the 2018 master formula for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(\varepsilon '/\varepsilon )_\text {BSM}$$\end{document}(ε′/ε)BSM valid in any extension beyond the SM without additional light degrees of freedom. We provide new values of the penguin parameters \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$B_6^{(1/2)}(\mu )$$\end{document}B6(1/2)(μ) and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$B_8^{(3/2)}(\mu )$$\end{document}B8(3/2)(μ) at the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu $$\end{document}μ-scales used by the RBC-UKQCD collaboration and at lower scales \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(1\, \text {GeV})$$\end{document}O(1GeV) used by ChPT and Dual QCD (DQCD). We present semi-analytic formulae for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(\varepsilon '/\varepsilon )_\text {SM}$$\end{document}(ε′/ε)SM in terms of these parameters and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hat{\Omega }_\text {eff}$$\end{document}Ω^eff that summarizes isospin-breaking corrections to this ratio. We stress the importance of lattice calculations of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(\alpha _{\text {em}})$$\end{document}O(αem) contributions to the hadronic matrix elements necessary for the removal of renormalization scheme dependence at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(\alpha _{\text {em}})$$\end{document}O(αem) in the present analyses of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varepsilon '/\varepsilon $$\end{document}ε′/ε.


Introduction
The direct CP-violation in K → ππ decays, represented by the ratio ε /ε, plays a very important role in the tests of the Standard Model (SM) and more recently in constraining its possible extensions [1]. In the SM ε /ε is governed by QCD penguins (QCDP) but receives also an important contribution from the electroweak penguins (EWP), pointed out already in 1989 [2,3], that entering ε /ε with the opposite sign to QCDP suppress this ratio significantly. The partial cancellation of these two contributions in addition to the evaluation of the hadronic matrix elements of QCDP and EWP operators is the reason why even today a precise prediction for ε /ε in the SM is not available. Yet, significant progress has been made during the last years and the goal of our paper is to update the SM value of ε /ε taking into account all available informations both from lattice QCD (LQCD) and analytic approaches most relevant for the evaluation of the Wilson coefficients but presently also for the estimate of the isospinbreaking corrections to the isospin amplitudes.
The situation of ε /ε in the SM before April 20, 2020 has been summarized by us in [4]. In short there are presently three approaches to calculate hadronic matrix elements entering ε /ε: • Lattice QCD, lead by the RBC-UKQCD LQCD collaboration. Using their results from 2015 for K → ππ matrix elements [5,6] and including isospin-breaking corrections from [7,8] as done in [9,10], leads to a value for ε /ε in the ballpark of (1 − 2) × 10 −4 . Although exhibiting a large error of 5 × 10 −4 the result lies one order of magnitude below the data. Taking these analyses at face value one could talk about an ε /ε anomaly of at most 3 σ . • DQCD approach [11][12][13], which gave a support to these values and moreover provided an upper bound on ε /ε in the ballpark of 8 × 10 −4 . The main QCD dynamics suppressing ε /ε in this approach is represented by the meson evolution, which is necessary to match long-distance contributions to short-distance ones. On the other hand it has been argued in [13] that final state interactions (FSI) should have only a minor impact on ε /ε and the quoted bound does not include them. • Chiral Perturbation theory (ChPT) [14][15][16] where, using ideas from ChPT, the authors found ε /ε = (14 ± 5) × 10 −4 attributing an important role to FSI in this result. While in agreement with the measurement, the large uncertainty, that expresses the difficulties in matching long-distance and short-distance contributions in this framework, does not allow for any clear-cut conclusions. 1 In view of the fact that LQCD calculations contain both the meson evolution 2 and FSI, while the estimate of ε /ε in the other two approaches does not include one of them, we have recently proposed the optimal strategy for the evaluation of ε /ε as of 2020 [4,19] 1. Use LQCD results for hadronic matrix elements of the dominant QCDP and EWP operators Q 6 and Q 8 , respectively. They are represented by the parameters B (1/2) 6 and B (3/2) 8 defined in Sect. 2. On the other hand the hadronic matrix elements of (V − A) ⊗ (V − A) operators can be determined from the experimental data on the real parts of the K → ππ amplitudes as performed in [9,20]. In fact this procedure has been recently adopted with slight modifications by the RBC-UKQCD collaboration [21] with the goal to decrease their errors. This procedure is clearly legitimate when testing the consistency of the SM with the data, if the reduction of the uncertainties is significant. But the RBC-UKQCD analysis shows only a little gain and therefore we will not use it in the present work, but rather base the full analysis on all hadronic matrix elements from RBC-UKQCD that are absolutely free from new physics. We will only use the experimental values of Re A 0,2 in the basic formula for ε /ε because they automatically take possible NP contributions into account. 2. Include isospin-breaking corrections from ChPT [15] that are compatible due to large uncertainties with the results obtained already 33 years ago in [22]. Very recently the latter analysis has been updated [23] and we will include these new findings as well. 3. Include NNLO QCD contributions to EWP in [24] thereby reducing the unphysical scale and renormalization scheme dependences in the matching at μ W = O(m W ), with the largest part due to the top-quark mass.
The removal of the dependence on μ c at NNLO has still to be done, see also the next point. 4. Take into account NNLO QCD contributions to QCDP [25,26]. This reduces the left-over renormalization scale uncertainties present at the NLO level, in particular those due to the matching scale μ c .
Recently significant progress in the estimate of ε /ε in the SM has been made through the improved values of the K → ππ hadronic matrix elements presented by the RBC-UKQCD collaboration [21]. Not only statistical errors have been significantly decreased but also a better agreement with the experimental values of the ππ strong interaction phases δ 0,2 has been obtained. The RBC-UKQCD collaboration, using their new results for the hadronic matrix elements and known Wilson coefficients at the NLO level [20,[27][28][29][30][31] but not accounting for isospin-breaking corrections, finds [21] (ε /ε) SM = (21.7 ± 8.4) × 10 −4 , (RBC-UKQCD − 2020) (1) to be compared with the experimental world average from NA48 [32] and KTeV [33,34] collaborations, While the result in (1) is in full agreement with the experimental value in (2) the theoretical error of 39% does not allow for clear cut conclusions whether some amount of new physics contributions is present in ε /ε or not. The same is the case of the earlier updated ChPT analysis [15], which resulted in with an error of 36%, very close to the LQCD one. But it should be remarked that with the present best values of the CKM parameters as used by us the central value in (3) would be raised to 15.0 × 10 −4 . Despite large errors both results deviate significantly from the DQCD values of ε /ε in the ballpark of 5 × 10 −4 stressed in particular in [17]. While there is no question about that meson evolution necessary for a proper matching between Wilson coefficients and hadronic matrix elements at scales O(1 GeV) must play a role in the evaluation of ε /ε it appears from present RBC-UKQCD results that precisely in the case of the matrix element of the Q 6 operator its suppression is overcompensated by other QCD dynamics which was hidden due to the contamination of the excited ππ states present in their 2015 analysis. It has been removed in the latest analysis. In fact as we will see soon the value of ε /ε obtained using the optimal procedure with hadronic matrix elements from [21], agrees very well with the one advocated in [15] and given in (3). Yet, it is not evident at present that FSI, as claimed by ChPT experts, are responsible for this agreement. Possibly other dynamical QCD effects apparently not taken into account both in the ChPT and DQCD approaches are responsible for the enhancement of ε /ε relative to DQCD expectations. 3 However, a clear-cut conclusion on this issue is difficult because of rather different techniques that are used in these three approaches. The fact that the central value in (3) differs significantly from the central LQCD value in (1) is dominantly due to the omission of isospin-breaking effects in the RBC-UKQCD prediction that are included in (3).
Even if the new improved calculation of K → ππ hadronic matrix elements in [21] is an important advance towards the accurate calculation of ε /ε, the result in (1) does not represent the present SM value of ε /ε properly. Indeed, as we emphasized in [4] the hadronic matrix elements in question are only a part of the ε /ε story. The three additional advances, listed in the context of the optimal strategy, that are not taken into account in the result in (1) are also important, in particular because they all lower the value of ε /ε. As we will demonstrate below, the final result for ε /ε differs significantly from the one obtained by the RBC-UKQCD collaboration. Indeed after including isospin-breaking effects from [15] that include the effects from the octet of lowest-lying mesons and NNLO QCD corrections to EWP contributions, we find using the hadronic matrix elements of RBC-UKQCD (ε /ε) (8) On the other hand including the singlet η 0 in this estimate one arrives at [23] (ε /ε) (9) Both results agree very well with experiment and with the ChPT expectations but in view of our comments on the ChPT analysis are on a more solid footing. We expect further reduction of ε /ε by roughly (5-10)% when NNLO QCD correc- 3 See also the discussion in [23] on this point. tions to QCD penguin contributions will be taken into account [25,26]. We look forward to the final results of these authors. Our paper is organized as follows. In Sect. 2, after recalling a number of basic formulae, we determine the parameters B (1/2) 6 and B (3/2) 8 using the recent RBC-UKQCD results and compare them with the expectations from ChPT [15] and DQCD [12,13]. It turns out that while there is a good agreement on the value of B (1/2) 6 between LQCD and ChPT, the rather precise value of B (3/2) 8 from RBC-UKQCD is by a factor of 1.5 larger than the ChPT one when both are evaluated at μ = 1 GeV. On the contrary, while there is a good agreement on the value of B (3/2) 8 between LQCD and DQCD [12], the most recent value of B (1/2) 6 from RBC-UKQCD is by a factor of two larger than the values quoted in [12,13]. We close this section with an updated formula for ε /ε in terms of B . In Sect. 3 we derive the results in (4) and (5) which take into account the updated isospinbreaking effects [15,23] and also NNLO QCD corrections to EWP contributions [24]. We also perform a detailed anatomy of various contributions. In Sect. 4 we update the BSM master formula for ε /ε [35,36] in view of the new RBC-UKQCD results. A brief summary and an outlook are given in Sect. 5. Some additional information on the numerical analysis are given in appendices. This includes the values of the hadronic matrix elements from RBC-UKQCD and the Wilson coefficients at various scales. We discuss in detail the effect of isospin-breaking corrections present in the renormalization group (RG) flow on ε /ε in Appendix C.

Preliminaries
The amplitudes for K 0 → (ππ) I , with I = 0, 2 denoting strong isospin of the final state, are given as where z i (μ) and y i (μ) are the S = 1 Wilson coefficients and Q i (μ) 0,2 the hadronic matrix elements of the operators Q i , both in the MS scheme at the low-energy factorization scale μ in the N f = 3 flavour theory [20]. By convention the strong phase shifts δ 0,2 are not included in A 0,2 , and therefore the Q i (μ) 0,2 are real-valued. Further The real parts Re A 0,2 are given entirely by the z i , because the y i are strongly suppressed by τ ∼ O(10 −3 ), on the other hand the imaginary parts Im A 0,2 ∝ Im(V * ts V td ) and depend only on y i . The Wilson coefficients of the QCD penguin (QCDP) operators i = 3, . . . , 6 are usually larger compared to those of the electroweak penguin (EWP) operators i = 7, . . . , 10, as can be seen in Table 7.
The scheme and scale dependences cancel between the Wilson coefficients and the matrix elements individually in A 0 and A 2 . We will take advantage of this freedom to use different scales μ 0 and μ 2 in the evaluation of A 0 and A 2 , respectively. In particular we choose the values at which the RBC-UKQCD lattice collaboration presents their results of the I = 0 [21] and I = 2 [6] matrix elements. There are only seven linearly independent Q i (μ) 0 and three linearly independent Q i (μ) 2 in the N f = 3 flavour theory [9,20], since RBC-UKQCD work in the isospin-symmetric limit, where also QED corrections are not included yet.
We remind that the amplitudes A 0,2 and the strong phase shifts δ 0,2 are related to the decay amplitudes relevant for ε /ε as follows with the experimental values of A 0,2 for h = 1 given in Table 2, whereas RBC-UKQCD works with the convention h = √ 3/2. These relations are valid also in the presence of finite QED corrections, as long as virtual infrared-divergent contributions, and also Coulomb corrections, are properly subtracted and combined with real photon radiation [8] when determining the amplitudes and phases from data.

Basic formula for ε /ε
As in [9], our starting expression is the formula where [8,15] Here a andˆ eff summarise isospin-breaking corrections. The latter include strong isospin violation (m u = m d ), the correction to the isospin limit coming from I = 5/2 transitions and electromagnetic corrections as first summarized in [7,8] and recently updated in [15] These analyses are based on the so-called octet scheme which includes only the octet of the lowest-lying pseudoscalar mesons. The inclusion of the singlet η 0 in the nonet scheme has been known already for 33 years [22,37] to give stronger suppression of ε /ε through the η − η mixing, but only very recently this estimate has been updated and put on a more solid basis than it was possible in 1987. Witĥ eff =ˆ the role of isospin-breaking effects is enhanced relative to the ChPT estimate in (12). The inclusion of the isospin-breaking corrections requires a modification in the evaluation of the Im A 0 part in ε /ε as follows [9] such that only leading isospin-breaking corrections are included.
A strong reduction of the uncertainty of ε /ε can be achieved firstly [20] by the use of the experimental values of Re A 0,2 in the denominators of (10). Secondly, the real parts of the relations (6) and (7) allow to eliminate one Q j (μ 0 ) 0 and one Q k (μ 2 ) 2 , respectively, in favour of the measured values of Re A 0 and Re A 2 , respectively. These can then be used in the numerators Im A 0 and Im A 2 , as proposed in [9]. The particular choice of j and k is subject to optimisation. However, as already announced previously we will not use this procedure here.
The real parts of the isospin amplitudes A 0,2 in (10) are then extracted from the branching ratios on K → ππ decays in the isospin limit. In the limit a = 1 andˆ eff = 0 the formula in (10) reduces to the one used by RBC-UKQCD [21], where all isospin breaking-corrections except for EWP contributions at the NLO level have been set to zero.

Extracting B
(1/2) 6 and B parameters, that enter the formula (20), are defined as follows with [22,38] in the large-N limit. We have introduced the factor h in order to emphasize different normalizations of these matrix elements present in the literature. We find from the latest RBC-UKQCD results for I = 0 [21] matrix elements at the scales μ = 1 GeV, μ = μ c = 1.3 GeV and μ = μ 0 = 4.006 GeV and for I = 2 from [6] for μ = 1 GeV, μ c = 1.3 GeV and to be compared with the 2015 values B has been changed by a factor of more than two, but with slightly larger uncertainty, which would correspond to a 2.6 σ discrepancy. However, in view that the systematic uncertainty of the 2015 results for the I = 0 matrix elements has been underestimated [21], the uncertainty quoted for the 2015 result of B (1/2) 6 (μ c ) must not be taken at face value anymore.
The new value of B is in the ballpark of values advocated in [15], but it is unclear to us at present whether this is a numerical coincidence or due to FSI dynamics. Moreover, the large uncertainty in the value of B does not yet rule out the values of B (1/2) 6 < 1.0 as expected from the DQCD approach [12]. Similar, the decrease of both parameters with increased μ, pointed out already in [20] and seen above, is also present below 1 GeV within the DQCD allowing smooth matching between hadronic matrix elements and Table 1 Coefficients entering the semi-analytic formula (20), when the amplitudes A 0 and A 2 are evaluated at the scales μ 0 and μ 2 , respectively, for different choices of (μ 0 , μ 2 ). The central values for (ε /ε) (8,9) with these approximate formulas are given in the last two lines.
8.00 8.79 9.85 Wilson coefficients. On the other hand it turns out that while there is a good agreement on the value of B (3/2) 8 between LQCD and DQCD [12], its rather precise value from RBC-UKQCD is by a factor of 1.5 larger than the ChPT one, in the ballpark of 0.55, when both are evaluated at μ = 1 GeV.

An analytic formula for ε /ε
As is well-known and shown also in the full analysis later, ε /ε is strongly dominated by the two terms ∝ Q 6 0 ∼ B . For convenience we provide a semi-analytic result of ε /ε in terms of these two parameters. Contrary to [4,9], we evaluate A 0 and A 2 at the two different scales μ 0 and μ 2 and use now for the remaining matrix elements the RBC-UKQCD results. Then with the coefficients given in Table 1 for various choices of (μ 0 , μ 2 ). The numerical input of the various parameters entering (20) is given in Table 2 and details on the Wilson coefficients at scales μ 0,2 are collected in Appendix B. The quark masses in (15) and (16) have been calculated as well at the two scales μ 0 and μ 2 , respectively. The coefficients a j i are comparable to [4], but differ because of the updated values for the remaining I = 0 matrix elements and changed values of the down-and strange-quark masses. Note that a EWP contains the I = 0 and I = 2 contributions of the EWPs at the scales (μ 0 , μ 2 ) = (4.006, 3.0) GeV, where the QCDP matrix elements for I = 2 are zero because the lattice calculation is done in the isospin limit. In general, when using the RG equations to evolve these matrix elements to dif-ferent scales μ 0,2 , the isospin-breaking in quark charges in the RG flow lead to nonvanishing I = 2 QCDP matrix elements that would also contribute to a EWP . As explained in more detail in Appendix C, we evolve the matrix elements of the operators from the initial scales μ 0,2 to the scales μ = 1.3 GeV and 1.0 GeV only with NLO QCD RG equations instead of NLO QCD × QED, which maintains isospin relations for these matrix elements.

ε /ε in the Standard Model
The new results for the I = 0 matrix elements from RBC-UKQCD imply a modification of ε /ε in the SM relative to those values presented in 2015 in [5,9,10], taking into account additional advances listed in Sect. 1. Here we include the isospin-breaking correctionsˆ eff and NNLO QCD corrections to EWPs calculated in [24]. Both contributions lead to a considerable reduction of ε /ε, as discussed previously [4]. Note that the RBC-UKQCD collaboration [21] prefers to use the magnitude of the isospin-breaking corrections from ChPT in the octet scheme (12) exclusively as an estimate of their size, thereby introducing an additional large uncertainty in ε /ε. In contrast to previous predictions [4,9], here we use in obtaining the final result for ε /ε directly the LQCD values of matrix elements Q i (μ 0 ) 0 and Q i (μ 2 ) 2 . For the interested readers, we provided the updated values of the two most important bag factors B We find for the amplitudes (h = 1) and where NNLO QCD corrections have been included in EWP parts [4]. The statistical uncertainties due to the matrix elements (ME, stat) were determined including the available correlations for I = 0, whereas the systematic ones (ME, syst) are based on the overall 15.7 % for I = 0 and (3 − 6) % for I = 2, as estimated by RBC-UKQCD in [6,21], respectively. For comparison, these values are very close to the RBC-UKQCD predictions Re A 0 = 24.44 × 10 −8 GeV, Re A 2 = 1.22 × 10 −8 GeV, Im A 0 = −5.70 × 10 −11 GeV, Im A 2 = −6.81 × 10 −13 GeV, from Eqs. (77a, 85, 90) [21] and Eq. (64) [6], respectively. 5 The scale uncertainties are obtained by varying μ c ∈ [1.0, 3.0] GeV and μ W ∈ [50,140] GeV for the NLO expressions, shown in Fig. 1. Note that we use m t (μ W ), and hence the μ W variation includes the top-mass scheme dependence. We emphasize that the μ W uncertainty for Im A 0,2 , and ε /ε, is very conservative, because we actually include here partial NNLO QCD corrections to EWPs [24], which remove the implicit μ W dependence associated with the top-quark mass and some of the explicit μ W dependence as well, see also [4] for more details. The parametric uncertainty due to the input value for the top-quark mass in Table 6 is denoted by "m t ".
The various relative contributions of the operators to Re A 0,2 and Im A 0,2 are listed in Table 3 when using μ 0,2 . These numbers show that Re A 0,2 are dominated by the current-current operators. In Re A 0 the Q 2 dominates with almost 96%, whereas the Q 1 and Q 6 contributions of about 12% cancel each other and there are subleading 2% and 1% contributions from Q 4 and Q 5 . In Re A 2 the Q 2 of 129% and the Q 1 of 27% enter with opposite signs and there is a subleading contribution from Q 8 of − 1.5%. On the other hand the Im A 0 is dominated by QCDP operators, where the 121% contribution of Q 6 is mainly reduced by Q 4 and Q 5 . The Im A 2 is dominated by EWPs, in particular by 122% due to Q 8 , which is partially cancelled by Q 9 . The 5% corrections from Q 7 and Q 10 cancel each other.
In the SM ε /ε receives contributions from QCDP and EWP via the I = 0 matrix elements and from EWP via the I = 2 matrix elements, that exhibit quite some hierarchies as can be seen in (36) and (38), respectively. These hierarchies are strongly counteracted by those present in the Wilson coefficients y i at the two scales μ 0 = 4.006 GeV and μ 2 = 3.0 GeV, where we evaluate Im A 0 and Im A 2 , respectively. This is illustrated by the following semi-analytic results of ε /ε that include the NNLO QCD corrections to EWPs [24] Here the experimental values of Re A 0,2 have been used only in the denominator of (10). As a remark on the side, we note that at the scale μ 2 = 3 GeV the relations Q 3,4,5,6 2 = 0 hold because the RBC-UKQCD calculations so far do not 5 We use here h = 1 as opposed to RBC-UKQCD collaboration that uses h = √ 3/2.  Table 3 The contribution in % of each operator to Re A 0,2 and Im A 0,2 at μ 0,2 . include isospin-breaking corrections, neither due to quark masses nor quark charges. In consequence no such contributions appear in (25). Moreover isospin relations (37) have been used to substitute Q 1,2,10 2 → Q 9 2 . As mentioned before in Sect. 2.4, a straight-forward application of the NLO QCD × QED RG equations to the matrix elements to evolve them to some different scale will generate nonvanishing Q 3,4,5,6 2 , because the RG flow includes isospinbreaking effects from quark charges. The effect on ε /ε is discussed in further detail in Appendix C, where we provide the analogous result to (25) at μ = 1.3 GeV. The hierarchy of the Wilson coefficients signaled for instance by large coefficients in front of Q 3,4 0 is strongly counteracted by a hierarchy in the hadronic matrix elements modifying the pattern of the various contributions: where the "| i,I " indicate the origin of the contribution. This shows much clearer the relevance of Q 6 0 ∼ B shows the contributions of QCDP in I = 0 and the partial cancellation of EWP contributions from I = 0 and I = 2. Note that this statement is scale dependent, i.e. at some other scales μ 0,2 the composition changes slightly due to RG flow. The final result for a = 1.017, usingˆ (8) eff = 0.17 ± 0.09 in the octet scheme (12), with NNLO QCD in EWP and other parameters as collected in Tables 2 and 6 is (ε /ε) (8) Considering the new valueˆ (9) eff = 0.29 ± 0.07 from the nonet scheme (13) we find: (ε /ε) (9) There is a statistical error due to the matrix elements from the lattice, based on covariance matrices for I = 0, propagated with Monte Carlo methods as well as individually available statistical errors for I = 2 matrix elements. The systematic uncertainty due to various sources related to the lattice approach is entirely dominated by the 15.7 % systematic error of Q 6 0 in Im A 0 . The isospin-breaking corrections to QCDP from ChPT, summarized inˆ eff in (11), contribute a relative uncertainty of 15 %. There is an overall relative uncertainty of 5.5 % from Im λ t due to the CKM input.
In our prediction we made only use of the experimental values Re A 0,2 | exp in the denominator of (10). As proposed in [9], in addition also in the numerator one of the I = 0 and one of the I = 2 matrix elements could be eliminated in favour of Re A 0,2 | exp to improve the accuracy in the framework of the SM. Here we did not adapt this strategy, because in agreement with the RBC-UKQCD collaboration [21], we did not find evidence for a substantial improvement when employing it to the I = 0 amplitude. It must be also noted that this strategy leads to a slightly reduced value of Im A 0 compared to the result without the additional information from Re A 0 | exp .
The " I = 1/2 rule" is given by the ratio and agrees with the experimental result 22.35 ± 0.05. Our value almost coincides with the RBC-UKQCD prediction. The RBC-UKQCD lattice results show that QCD dynamics, present dominantly in current-current operators, is responsible for this large ratio thereby confirming the findings within DQCD obtained many years ago [11,48]. This is also seen in another recent LQCD analysis [49].

BSM master formula
In this section we report the updated master formula coefficients describing the new physics effects beyond the SM (BSM) in ε /ε, which were first presented in [35,36]. The BSM contribution to ε /ε is given by the weight factors P i for each Wilson coefficient C i (μ EW ) of the operators and their chirality-flipped counterparts listed in Table 4. The P i (μ EW ) contain the information of the RG evolution from the low-energy scale μ to the electroweak (EW) scale μ EW and are linearly dependent on the hadronic matrix elements of the operators at the scale μ, such that the μ-dependence cancels. The master formula takes the simple form with the N f = 5 effective Hamiltonian leading to dimensionless Wilson coefficients and P i factors.
The sum runs over all Wilson coefficients of the operators in Table 4. These operators are a complete basis for non-leptonic S = 1 transitions in the absence of any other light degrees of freedom [36]. The Wilson coefficients and their weight factors are evaluated at the particular value μ EW = 160 GeV of the EW scale. For more details we refer to [35,36].
In Table 4 we summarize the updated P i factors after taking into account the most recent I = 0 matrix elements reported by RBC-UKQCD [21]. Table 4 has been obtained by taking into account the tree-level matching [50] and oneloop running [51] below the EW scale using the public codes wilson [52] and WCxf [53]. Only the P i factors of operators in Class A are affected by this change, since they depend exclusively on matrix elements of the SM operators. In all other classes the P i 's depend on matrix elements of BSM operators or the chromomagnetic operator Q 8g and remain unchanged. The central values as well as statistical and systematic uncertainties of the I = 0, 2 matrix elements of all operators are listed in Table 5 at the common scale μ = 1.3 GeV.
The changes are moderate of not more than 30% for operators that contribute directly to K → ππ, whereas changes can be larger for those operators (with s, c, b-quarks) that enter via RG running from the EW scale down to the lowenergy scale and have smaller coefficients. The last column of Table 4 shows the suppression scale that would generate (ε /ε) BSM = 10 −3 for C i = 1/ 2 , assuming the presence of only this particular operator. For comparison, the theory uncertainty of the SM prediction (28) is about 0.6 × 10 −3 . The scale is strongly dependent on the uncertainties of the matrix elements, which did not all decrease in the latest RBC-UKQCD predictions. A comparison to the previous values [35] shows a slight increase of for the first seven operators, which contribute directly to K → ππ. In general also increases for the remaining Class A operators, with a few exceptions, pushing the NP scale also in these cases up, even though they are entering only via RG mixing. This shows that the new results for the matrix elements from RBC-UKQCD will lead to stronger bounds on CP violation beyond the SM.
Eventually we point out that the large increase of the central value of (ε /ε) SM in the SM from ∼ (1 − 2) × 10 −4 with the 2015 RBC-UKQCD results to ∼ 14 × 10 −4 with the 2020 results constitutes one order of magnitude and hence has significant impact on excluded regions of parameter spaces of BSM scenarios. The 2015 SM predictions [5,9,10] suggested a strong anomaly with a constructive (ε /ε) BSM ≈ (5 − 15) × 10 −4 to reach agreement with the experimental value (ε /ε) exp = (16.6 ± 2.3) × 10 −4 . Contrary, the (ε /ε) SM predictions based on 2020 results do not show anymore an anomaly, but allow now for both, a constructive and destructive interference, that can be still sizable in view of the large theory uncertainties as a rough 1 σ range. The complete error propagation can be obtained properly for general BSM scenarios with the master formula, which is implemented in the public code flavio [54,55]. Despite the large uncertainties, ε /ε was and remains one of the strongest constraints on CP violation in the quark-flavour sector, as has been shown for different BSM scenarios in the past. The BSM studies based on the 2015 SM predictions of ε /ε used mostly the working hypothesis of a constructive (ε /ε) BSM of similar size, see references in [4], and the obtained conclusions for 0 < (ε /ε) BSM are still mostly valid.

Summary and outlook
Our final result for ε /ε in (5) differs significantly from the one of the RBC-UKQCD collaboration but in view of large uncertainties in both results they are in agreement with each other and with experiment. But as emphasized in [23] the perfect agreement of (5) with the ChPT result in (3) is a pure numerical coincidence because the latter was obtained witĥ (8) eff in place ofˆ (9) eff and with the values of B The recent advances in LQCD allow us to hope that in the coming years we should be able to have a value of ε /ε within the SM with a comparable error to the experimental one. In order to reach this goal and thereby to obtain an assessment on the allowed room for NP contributions to ε /ε it is important to perform a number of steps: • A more precise determination of Q 6 (μ 0 ) 0 or B (1/2) 6 (μ 0 ). At least a second LQCD collaboration should calculate ε /ε, in order to confirm the large enhancement of B (1/2) 6 found by RBC-UKQCD that has not been identified in DQCD. Also the errors in other matrix elements should be decreased.
• A more precise determination ofˆ eff . In particular in LQCD calculations isospin-breaking corrections and O(α em ) corrections in hadronic matrix elements required for the removal of renormalization scheme dependence at this order should be taken into account. The present status is summarized in [56].
• A more precise determination of the short distance contributions, especially in the QCD penguin sector, which in the context of the RBC-UKQCD analysis will decrease the sensitivity to the matching scale μ c . Despite the fact that the NNLO analysis of QCD corrections to EWP con- tributions practically removed the sensitivity of ε /ε to the renormalization scheme of the top-quark mass and μ W , our analysis shows that the significant μ c uncertainty in the EWP sector still has to be removed through the matching of N f = 4 to N f = 3 effective theory at the NNLO level.
• The computation of the BSM K → ππ hadronic matrix elements of four-quark operators by lattice QCD, which are presently known only from the DQCD approach [57].
Several BSM analyses of ε /ε have been performed, which are collected in [4]. A recent example of a Z model with explicit gauge anomaly cancellation has been discussed in [58]. Furthermore leptoquark models, except the U 1 model, would not be able to explain large deviations of the SM value from the data due to constraints coming from rare K decays [59]. This underlines the importance of correlations of ε /ε with other observables in NP scenarios. The new SM value in (5) removes the difficulties of leptoquark models pointed out in [59], but these problems could return with an improved analyses of ε /ε within the SM.
Furthermore the lessons from the SMEFT analysis in [36] should be useful in this respect. Such general analyses allow to take into account constraints from other processes such as collider processes, electroweak precision tests, neutral meson mixing as well as electric dipole moments. Finally the master formula for ε /ε presented in [35] valid for any BSM scenario should facilitate the derivation of constraints on CP-violating phases beyond the SM imposed by ε /ε. In this respect we point out that also Re A 2 has a very precise SM prediction and can be predicted rather precisely also in BSM scenarios, providing thus a second observable besides ε /ε to constrain also real parts of the Wilson coefficients of non-leptonic S = 1 operators.

Data Availability Statement
This manuscript has no associated data or the data will not be deposited. [Authors' comment: No data was created during course of this project].
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/. Funded by SCOAP 3 .

A Hadronic matrix elements
Here we collect the input for the K → (ππ) I isospin I = 0, 2 hadronic matrix elements of the relevant operators in the traditional SM basis [20]. They are given for the MS scheme in order to combine them with the Wilson coefficients in Appendix B at the scale chosen by RBC-UKQCD. In addition we provide these matrix elements together with the complete set of non-leptonic S = 1 operators beyond the SM in Table 5 at the common scale μ = 1.3 GeV. These results can be used for new physics studies using the master formula of ε /ε in [35,36] that we updated in Section 4.
The new results for I = 0 matrix elements of the SM operators from the year 2020 are from the RBC-UKQCD lattice collaboration [21]. They are given at the scale μ 0 = 4.006 GeV in the N f = 2 + 1 flavour theory. As the current lattice calculation works in the isospin limit, out of the ten Q 1...10 0 there are only seven linearly independent (for h = 1): (18)(14), Q 2 0 = +0.120 (12)(19), The first and second errors are of statistical and systematic origin, respectively. In particular the statistical error comprises also a covariance matrix provided in [21] that we use for the uncertainty propagation in Table 5 and predictions of ε /ε. The systematic uncertainty due to various sources in the lattice approach was estimated to be 15.7 % (see table XXV [21]) for each matrix element without providing correlations. The I = 2 matrix elements of the SM operators are also from RBC-UKQCD [6] from the year 2015 for N f = 2 + 1. In particular we use the results from the RI-SMOM ( / q, / q) scheme (Table XVI) and convert them to the MS scheme with the scheme conversion factor Eq. (66) in [60]. The matrix elements are given at μ 2 = 3 GeV where they fulfill the isospin relations and reduce to three independent ones (for h = 1) Q 7 2 = 0.2340(52)(70), Q 8 2 = 1.072(28)(64), where we have increased the systematic uncertainty of the results of the RI-SMOM ( / q, / q) by adding in quadrature the difference of the results in the RI-SMOM ( / q, / q) and the RI-SMOM (γ , γ ) schemes as given in [6] to account for this   Table 5, see also [36]. The matrix elements of operators beyond the SM were calculated using DQCD in [57]. The single error is of parametric and systematic origin. The matrix element of the chromo-magnetic dipole operator O 8g has been calculated in [61,62] in 2017/18. Note that we use here the normalization of [35,36].

B Wilson coefficients
Here we summarize the S = 1 Wilson coefficients at the various scales used in our analysis in the NDR-MS scheme using the NLO RG evolution from [20]. The numerical input entering the Wilson coefficients is fixed to values in Table 6 [41]. We remind that the threshold corrections enter here for the first time at NLO, hence to be able to cancel some of the renormalization scheme dependences of the bottom-and charm-quark masses, one has to go to the NNLO order, as for example done in [25] in the case of QCD penguins.
The top quark mass m t (μ t ) is in the MS scheme for μ t = μ W , obtained from the pole mass 6 value given in Table 6: m t (m t ) = 163.5 GeV and m t (μ W ) = 173.2 GeV. We follow [4] and include also important NNLO matching corrections [24] that resolve the NLO renormalization scheme ambiguities for our choice μ t = μ W via the modifications of y 7,...,10 (μ) at the low-energy scale of about 1.07, 1.07, 0.89 and 0.76 leading to the NNLO' values in Table 7, which we adapt in the numerics. For further details we refer to [4]. The prime in this indicates that still small O(α W α s sin 2 θ W ) corrections are not included.

C Isospin-breaking RG effects
In this appendix we comment on isospin-breaking effects in the RG flow present in the anomalous dimension matrices that govern the scale dependence of Wilson coefficients and matrix elements of the operators. They are due to quark charges present in the definitions of the EWP operators and due to QED corrections, and are of purely perturbative origin known up to NLO in QCD × QED [28,29,31]. To compute the Wilson coefficients at the low-energy scale in the N f = 3 theory, the isospin-breaking effects in the RG flow are combined with initial Wilson coefficients, which contain further isospin-breaking effects of the SM at the electroweak scale.
In the predictions of the amplitudes A 0,2 as well as for ε /ε the low-energy Wilson coefficients are multiplied by the corresponding matrix elements, which in principle leads to the cancellation of the renormalization scheme dependence. The matrix elements contain the dynamics of QCD × QED at scales below the low-energy scale, i.e. in the nonperturbative regime of QCD. In the calculation with lattice methods so far no isospin-breaking corrections have been included, i.e. a purely isospin-symmetric QCD setup is used by RBC-UKQCD that neglects effects due to different quark masses m u = m d as well as QED effects 7 due to different quark charges. For this reason, firstly, scheme dependences between Wilson coefficients and matrix elements can only cancel for the isosymmetric parts. Secondly, the initial values of the I = 2 matrix elements fulfill the isospin relations (37) at the scale μ 2 = 3 GeV, where they have been calculated. These isospin relations become broken by the full QCD × QED RG flow when evolving the matrix elements to a different scale, as can be seen for the example in Table 5 at μ = 1.3 GeV, as for example 3/2 Q 1 2 = 3/2 Q 2 2 = Q 9 2 = Q 10 2 .
The central values of the matrix elements in Table 5 lead then to (ε /ε) (8) = 18.1 × 10 −4 and (ε /ε) (9) = 14.6 × 10 −4 , respectively. These results differ by about +0.7 × 10 −4 from the predictions (28) and (29), which are rather small variations in view of larger uncertainties in the prediction of ε /ε. They might be viewed as an estimate of isospin-breaking corrections due to quark charges. However, the lacking contributions in the nonperturbative matrix elements could cancel them in large parts, as one would expect on the basis of renormalization scheme cancellations. We note that they are unrelated to isopin-breaking quark-mass effects, which are entirely of nonperturbative origin.
In the case of our preferred choice of calculating ε /ε at the scales μ 0,2 where RBC-UKQCD has calculated I = 0, 2 matrix elements in the NDR-MS scheme, these isospinbreaking corrections are entirely contained in the Wilson coefficients. At a different scale μ one might impose the isospin relations for matrix elements by hand, but a more consistent possibility is the use of the isospin-conserving parts of the RG equation only, i.e. the NLO QCD evolution obtained by setting α em = 0. Then the I = 2 matrix elements fulfill the isospin relations exactly, whereas for I = 0 matrix elements they are slightly broken, because we have kept the small QCD mixing of EWP into QCDP operators. With this approach we obtain at μ = 1.3 GeV the predictions (ε /ε) (8) = 17.1 × 10 −4 and (ε /ε) (9) = 13.6 × 10 −4 , respectively, differing by about −0.3 × 10 −4 from the predictions (28) and (29) at μ 0,2 . In particular we have used the NLO QCD evolution of matrix elements for the semi-analytic equation (20) and coefficients in Table 1.
In conclusion these studies show that nonvanishing I = 2 matrix elements of QCDP operators due to isospin-breaking effects can have a non-negligible impact on ε /ε. The size of such effects is completely unknown at present, contrary to isospin-breaking effects from I = 0 matrix elements of QCDP operators, contained inˆ eff . Here, the obtained variations are based on the residual scheme dependence from the isopin-breaking contributions due to quark charges in perturbation theory, contained in the anomalous dimensions. Their impact on ε /ε is numerically subleading compared to other, currently much larger uncertainties.