331 Models Facing the Tensions in $\Delta F=2$ Processes with the Impact on $\varepsilon^\prime/\varepsilon$, $B_s\to\mu^+\mu^-$ and $B\to K^*\mu^+\mu^-$

Motivated by the recently improved lattice QCD results on the hadronic matrix elements entering $\Delta M_{s,d}$ in $B_{s,d}^0-\bar B_{s,d}^0$ mixings and the resulting increased tensions between $\Delta M_{s,d}$ and $\varepsilon_K$ in the Standard Model and CMFV models, we demonstrate that these tensions can be removed in 331 models based on the gauge group $SU(3)_C\times SU(3)_L\times U(1)_X$ both for $M_{Z^\prime}$ in the LHC reach and well beyond it. But the implied new physics (NP) patterns in $\Delta F=1$ observables depend sensitively on the value of $|V_{cb}|$. Concentrating the analysis on three 331 models that have been selected by us previously on the basis of their performance in electroweak precision tests and $\varepsilon^\prime/\varepsilon$ we illustrate this for $|V_{cb}|=0.042$ and $|V_{cb}|=0.040$. We find that these new lattice data still allow for positive shifts in $\varepsilon^\prime/\varepsilon$ up to $6\times 10^{-4}$ for $M_{Z^\prime}=3~TeV$ for both values of $$|V_{cb}|$ but for $M_{Z^\prime}=10~TeV$ only for $|V_{cb}|=0.040$ such shifts can be obtained. NP effects in $B_s\to\mu^+\mu^-$ and in the Wilson coefficient $C_9$ are significantly larger in all three models for the case of $|V_{cb}|=0.040$. In particular in two models the rate for $B_s\to\mu^+\mu^-$ can be reduced by NP by $20\%$ for $M_{Z^\prime}=3~TeV$ resulting in values in the ballpark of central values from CMS and LHCb. In the third model a shift in $C_9$ up to $C_9^\text{NP}=-0.5$ is possible. We also consider the simplest 331 model, analyzed recently in the literature, in which $X=Y$, the usual hypercharge. We find that in this model NP effects in flavour observables are much smaller than in the three models with $X\not=Y$, in particular NP contributions to the ratio $\varepsilon^\prime/\varepsilon$ are very strongly suppressed.

The Standard Model (SM) describes globally the existing data on quark-flavour violating processes rather well [1] but with the reduction of experimental errors and increased precision in non-perturbative and perturbative QCD and electroweak calculations a number of tensions at the level of 2 − 3 σ seem to emerge in various seemingly unrelated observables. While some of these tensions could turn out to be the result of statistical fluctuations, underestimate of systematical and theoretical errors, it is not excluded that eventually they all signal the presence of some kind of new physics (NP). Therefore, it is interesting to investigate what this NP could be.
In the present paper we will address some of these tensions in 331 models based on the gauge group SU (3) C × SU (3) L × U (1) X [2,3]. As these models have much smaller number of new parameters than supersymmetric models, Randall-Sundrum scenarios and Littlest Higgs models, it is not evident that they can remove all present tensions simultaneously.
Our paper has been motivated by a recent analysis in [4] which demonstrates that the new lattice QCD results from Fermilab Lattice and MILC Collaborations [5] on B 0 s,d −B 0 s,d hadronic matrix elements imply a significant tension between ε K and ∆M s,d within the SM. The authors of [5] find also inconsistences between ∆M s,d and tree-level determination of |V cb |. But the simultaneous consideration of ε K and ∆M s,d in [4] also demonstrates that the tension between these two quantities cannot be removed for any value of |V cb |. Moreover, the situation worsens for other models with constrained MFV (CMFV), indicating the presence of new flavour-and CP-violating interactions beyond CMFV framework at work.
The question then arises how 331 models face this tension and what are the implications of new lattice results on other observables for which some departures from SM predictions have been identified. In particular, taking the results in [4,5] into account we want to concentrate our analysis on ε /ε, In this context the following facts should be recalled.
• Recent analyses in [6][7][8][9] find the ratio ε /ε in the SM to be significantly below the experimental world average from NA48 [10] and KTeV [11,12] collaborations. The recent analysis in the large N approach in [13] indicates that final state interactions will not modify this picture at least on the qualitative level. The analysis in [14] shows that CMFV models cannot cure this problem. In any case models providing an enhancement of ε /ε should be favoured from present perspective.
• The branching ratio for B s → µ + µ − measured by CMS and LHCb [15] has been always visibly below rather precise prediction of the SM [16]. The most recent result from ATLAS 1 , while not accurate, appears to confirm this picture and models suppressing the rate for this decay relative to its SM prediction appear to be favoured.
• LHCb data on B d → K(K * )µ + µ − indicate some departures from SM expectation although this issue is controversial. See [18,19] and references to the rich literature therein.
Assuming again that statistical fluctuations or underestimated errors are not responsible for these effects, significant NP contributions to the Wilson coefficient C 9 or C 9 and C 10 are required.
These three items have been already addressed by us within 331 models in the past [20][21][22][23]. In particular in [23] the issue of ε /ε anomaly has been addressed, while in [21,23] the last two items above have been considered. The main result of [23] is that among 24 331 models only three (M8, M9 and M16 in the terminology of [22]) have a chance to survive if an improved fit to electroweak precision observables relative to the SM is required and the ε /ε anomaly will be confirmed in the future. Two of them (M8 and M9) allowed simultaneously a suppression of the rate for B s → µ + µ − by 20% thereby bringing the theory closer to the data without any significant impact on the Wilson coefficient C 9 . The third model (M16) provided, simultaneously to the enhancement of ε /ε, a shift up to ∆C 9 = −0.6, softening the anomalies in B → K * µ + µ − , without any significant impact on B s → µ + µ − . While, before the ATLAS data on B s → µ + µ − , M16 seemed to be slightly favoured over M8 and M9, this data and the fact that the theoretical uncertainties in B s → µ + µ − are significantly smaller than in B d → K * µ + µ − make us believe that at the end models M8 and M9 have a bigger chance to survive.
However, the constraints from ∆F = 2 transitions used in [23], prior to the lattice QCD result in [5], were significantly weaker and it is of interest to investigate what is the impact of these new lattice results on our previous analyses and whether the increased tension between ε K and ∆M s,d within the SM pointed out in [4] can be removed in these three models. In fact one should recall that the mixing and CP-violation in B 0 s,d −B 0 s,d systems play very important roles in the determination of new parameters in 331 models [20] and it is not surprizing that our previous results will be indeed modified in a visible manner. In this context let us remark that within the SM, dependently on whether ∆M s or ε K has been used as a constraint, rather different values for |V cb | have been required to fit the data within the SM [4]: This in turn resulted in rather different predictions for rare K and B s,d decays as seen in Table 4 of [4].
In our most recent analysis in [23] we have performed numerical analysis for |V cb | in the ballpark of the higher value in (2), that is 0.042. In the present paper we will also use this value in order to see the impact of new lattice data on our previous results, but in addition we will perform the analysis with 0.040 which is in the ballpark of its lower value in (2). This will tell us whether 331 models can cope with the tensions in question for both values of |V cb | and whether the implications for NP effects are modified through this change of |V cb |.
The second motivation for a new analysis of 331 models is the following one. In our analyses and also in [24][25][26] the U (1) X factor in the gauge group SU (3) C × SU (3) L × U (1) X differed from the hypercharge gauge group U (1) Y . As various 331 models are characterized by two parameters β and tanβ defined through these analyses dealt with β = 0. Here T 3,8 and X are the diagonal generators of SU (3) L and U (1) X , respectively. Y represents U (1) Y and v i are the vacuum expectation values of scalar triplets responsible for the generation of down-and up-quark masses in these models.
Recently a special variant of 331 models with β = 0 or equivalently U (1) X = U (1) Y has been considered in [27]. Moreover, these authors set tanβ = 1 as this choice with β = 0 simplifies the model by eliminating Z − Z mixing studied by us in detail in [22] for β = 0. As this is the simplest among the 331 models, the question arises whether it is consistent with the flavour data in the setup in [27] and what are the implications for quark flavour observables for arbitrary tanβ when Z − Z mixing enters the game. In particular the comparison with our studies for β = 0 in [20][21][22][23] and in the present paper is of interest. As the authors of [27] did not address this question, to our knowledge this is the first quark flavour study of this simplified 331 model.
We will see that in the absence of Z − Z mixing the choice β = 0 provides a unique 331 model in which the phenomenologically successful relation is satisfied. Here C NP 9 and C NP 10 stand for the shifts in the Wilson coefficients relevant in particular for B → K * µ + µ − and B s → µ + µ − , respectively. This is good news. The bad news is that setting β = 0 modifies the values of all couplings relative to the ones in M8, M9 and M16 models. We find then that NP contributions to ε /ε in this simple model are at most 1 × 10 −4 for M Z = 3 TeV and decrease with increasing M Z . The effects in C NP 9 and C NP 10 are at most at the level of a few percent even if Z − Z mixing is taken into account. Thus the model fails in solving three anomalies listed above. But as we will see it is able to remove the tensions between ∆M s,d and ε K .
Our paper is organized as follows. In Section 2 we address the tensions between ∆M s,d and ε K in M8, M9 and M16 models and we update our analysis of ε /ε, B s → µ + µ − and C 9 in [23] taking new ∆F = 2 constraints from [5] into account and performing the analysis at two values of |V cb | as discussed above. In Section 3 we specify the existing formulae in 331 models to the case β = 0 but for arbitrary tanβ and we derive the results mentioned above. We conclude in Section 4.

Preliminaries
Let us recall that in these three models new flavour-violating effects are governed by tree-level Z exchanges with a subdominant role played by tree-level Z exchanges generated through Z −Z mixing. All the formulae for flavour observables in these models can be found in [20][21][22][23] and will not be repeated here. In particular the collection of formulae for Z couplings to quarks and leptons for arbitrary β are given in (17) and (18) of [21].
New sources of flavour and CP violation in 331 models are parametrized by new mixing parameters and phasess 13 ,s 23 , δ 1 , δ 2 withs 13 ands 23 positive definite and smaller than unity and 0 ≤ δ 1,2 ≤ 2π. They can be constrained by flavour observables as demonstrated in detail in [20]. The non-diagonal Z couplings relevant for K, B d and B s meson systems can be then parametrized respectively within an excellent approximation through s 13 and δ 1 can be determined from ∆M d and CP-asymmetry S ψK S whiles 23 and δ 2 from ∆M s and CP-asymmetry S ψφ . Then the parameters in the K system are fixed. This correlation tells us that the removal of tensions between ε K and ∆M s,d is not necessarily automatic in 331 models and constitutes an important test of these models. The remaining two parameters, except for M Z mass, are as seen in (3), β and tanβ. Moreover, the fermion representations of SM quarks under the SU (3) L group matter. The three models in question are then characterized by with F 1 and F 2 standing for two fermion representations. In F 1 the first two generations of quarks belong to triplets of SU (3) L , while the third generation of quarks to an antitriplet. In

Numerical Analysis
The difficulty in doing the numerical analysis are tensions between inclusive and exclusive determinations of the CKM elements |V cb | and |V ub |. The exclusive determinations have been summarized in [28] and are given as follows They are based on [5,[29][30][31][32]. The inclusive ones are summarized well in [33,34] |V cb | incl = (42.21 ± 0.78) · 10 −3 , We note that after the recent Belle data on B → D ν l [31], the exclusive and inclusive values of |V cb | are closer to each other than in the past. On the other hand in the case of |V ub | there is a very significant difference. Furthermore, after recent precise determinations of hadronic matrix elements entering ∆M s,d in B 0 s,d −B 0 s,d mixing by Fermilab Lattice and MILC Collaborations [5] there are significant tensions between tree-level determinations of |V cb | and |V ub | and ∆M s,d within the SM [5] and also the tensions between ε K and ∆M s,d [4] in this model. Moreover, as found in the latter paper, the value of the angle γ in the unitarity triangle extracted from the ratio ∆M d /∆M s and the CP-asymmetry S ψK S is with γ = (63.0 ± 2.1) • visibly smaller than it tree-level determination [35] γ In the present paper, as in [23], we will set first the CKM parameters to This choice is in the ballpark of exclusive determination of |V ub | in (10) and the inclusive one for |V cb | in (11). Moreover, it is in the ballpark of tree-level determination of γ. In view of new parameters in 331 models the value of γ does not follow from the ratio ∆M s /∆M d and S ψK S like in CMFV models and it is better to take γ from tree-level determinations as it is to first approximation not polluted by NP. Having the same CKM input as in our previous analysis will allow us to see the impact of new lattice data on phenomenology.
The choice in (13) is also motivated by the fact that NP contributions to ε K in 331 models are rather small for M Z of a few TeV and SM should perform well in this case. Indeed for this choice of CKM parameters we find and |ε K | in the SM only 4% below the data. Due to the presence of long distance effects in ∆M K also this value is compatible with the data. While the CKM parameters do not enter the shift in ε /ε and ε K , their choice matters in the predictions for NP contributions to ∆F = 2 observables in B 0 d,s −B 0 d,s systems and the rare B s,d decays. This is not only because of their interferences with SM contributions. The departure of SM predictions for ε K and ∆M s,d from the data depends on the CKM parameters, in particular on the value of |V cb |, and this has an impact on the allowed ranges of new parameters extracted from ∆F = 2 observables and consequently on final values of ε /ε, B(B s → µ + µ − ) and the shift in C 9 . We will illustrate this below by choosing also |V cb | = 0.040 which corresponds to its exclusive determination in (10). See (23).
Next, as in [20,23], we perform a simplified analysis of ∆M d,s , S ψK S and S ψφ in order to identify oases in the space of four parameters (5) for which these four observables are consistent with experiment. To this end we use the formulae for ∆F = 2 observables in [20,22] and set input parameters listed in Table 3 of our recent analysis in [23] at their central values. The only modifications in this input are the recently calculated parameters [5] while being consistent with the data for ε K , S ψK S and S ψφ . As we keep the input parameters at their central values, in order to take partially hadronic and experimental uncertainties into account we require the 331 models to reproduce the data for ∆M s,d within ±5% and the data on S ψK S and S ψφ within experimental 2σ ranges.
Specifically, our search is governed by the following allowed ranges: The ranges for ∆M s,d are smaller than used in [23] because of the reduced errors in (15).
We also impose the constraint on the ratio ∆M s /∆M d using [5] In the spirit of our simplified analysis we will keep this ratio at its central value in (20) but in order to take into account the uncertainty in ξ we will require that ∆M s /∆M d agrees with the data within ±5%. Specifically we will require that is satisfied.
In the case of ε K and ∆M K we will just proceed as in [23] imposing the ranges Having determined the ranges for the parameters (5) we can calculate all the remaining flavour observables of interest.
In Fig. 1 we show ∆M s,d vs. ε K in M8 and M16. Red dots represent central SM values and black dots the central experimental values. The experimental errors are negligible and the parametric and theoretical errors are represented by the allowed departure from them as explained above. These results do not depend on the fermion representation up to tiny effects from Z − Z mixing and consequently are practically the same for M8 and M16. As far as M9 is concerned all results presented in our paper are very similar to the ones in M8 and will not be shown. We observe that the tensions between ∆M s,d vs. ε K present in the SM can be easily removed in 331 models for M Z = 3 TeV.
In Fig. 2 we show ∆(ε /ε) versus ε K for M8 and M16 at M Z = 3 TeV. Taking the uncertainties due to charm contribution and CKM parameters in ε K into account the maximal  shifts in ε /ε in both models amount to 6 × 10 −4 , very similar to what we found in [23]. But NP effects in B s → µ + µ − and C 9 are smaller relative to the ones found in the latter paper by a factor of two. This is seen in Fig. 3, where we show correlations of ∆(ε /ε) with B s → µ + µ − (left panels) and with C NP 9 (right panels) for M8 and M16 and M Z = 3 TeV. In M8 the rate for B s → µ + µ − can be suppressed by 10% bringing the theory closer to the data in (24) [15]. Moreover, this happens for the largest shift in ε /ε. But the shift in C 9 is very small. In M16 the pattern is opposite with only a very small NP effects in B s → µ + µ − and a shift of −0.3 in C 9 which brings the theory closer to the data.

|V cb | and |V ub | Dependence
It is well known that ε K and ∆M s,d in the SM are sensitive functions of |V cb |. Moreover, ε K and S ψK S depend sensitively on |V ub |. Setting |V cb | and |V ub | to the values in (13) we have necessarily constrained the allowed range of NP parameters that are consistent with the data on ∆F = 2 observables. Changing |V cb | and |V ub | will necessarily modify this range and will modify NP contributions to flavour observables even if they do not depend directly on |V cb | and |V ub |. A sophisticated analysis which would include the uncertainties in both CKM elements from tree-level decays would wash out NP effects and would not teach us much about the impact of |V cb | and |V ub | on our results.
Therefore, we prefer to show how our results presented above are modified for a different value of |V cb | that we choose to be lower so that instead of (13) we now use We keep |V ub | and γ unchanged as this will allow us to see the role of |V cb | better. The dependence on γ in the observables in question is weak. The dependence of ∆M s,d on |V ub | is totally negligible. The inclusive value of |V ub | would compensate the decrease of |V cb | in ε K but would simultaneously have an impact on S ψK S shifting it in the ballpark of 0.80 within the SM. While in the SM this problem cannot be cured because of the absence of new CP-violating phases, in 331 models the presence of the phase δ 1 in (6) allows to satisfy the constraint on S ψK S in (19). We will demonstrate it below. This assures us that the tension between ε K and S ψK S for exclusive value of |V cb | present in the SM can be avoided within the 331 models. But as we would like to investigate the impact of the change in |V cb | on our results, we keep first |V ub | at its exclusive value. Moreover there is some kind of consensus in the community that in the case of |V ub | one can trust more exclusive determinations of this parameter than the inclusive ones. This is based on the fact that the exclusive determinations use formfactors from lattice QCD, which on the one hand are already rather precise and on the other hand do not require the assumptions like hadron duality necessary for the inclusive determination of |V ub |.
The results of this exercise are shown in Figs. 4-6. Comparing Fig. 4 with Fig. 1 it is evident that 331 models perform better for our nominal choice of CKM parameters in (13) than for a lower value of |V cb |. This is seen in particular in the case of ε K for which the maximal values of ε K are by 10% below the data. But, as discussed above, increasing |V ub | towards its inclusive value and taking the uncertainties in the QCD corrections to the charm contribution in ε K into account one can bring the theory much closer to data without violating the constraint on S ψK S in (19). We demonstrate this in Fig. 7 where we use |V ub | = 0.0042. Indeed ε K is now in a perfect agreement with the data. The slight increase of the maximal value of ∆(ε /ε) in this case will be analyzed in more details below. Therefore we can claim that 331 models also in this case remove the tensions in question, which is not possible within the SM. Interestingly, as we will see below for significantly higher values of M Z the removal of tensions for this value of |V cb | will be much easier. We refer to Section 4 in [23] for the explanation of this behaviour.
Next Figs. 5 and 6 should be compared with Figs. 2 and 3, respectively. We observe: • The correlation between ε /ε and ε K in Fig. 5 has a very different shape than in Fig. 2 but a shift of ε /ε of (5 − 6) × 10 −4 is possible. In fact this plot is similar to a corresponding plot in [23] obtained with CKM parameters in (13) but older hadronic matrix elements. This similarity is easy to understand. The increase of non-perturbative parameters in (15) has been roughly compensated by the decrease of |V cb |.
• The size of NP effects in B s → µ + µ − and C 9 is now larger than for our nominal value of |V cb | and similar to the ones found in [23]: suppression of the rate for B s → µ + µ − by 20% in the case of M8 and a shift of C 9 by −0.5 in M16 are possible. But what is interesting is that the decreased value of |V cb | lowers also the SM result for the B s → µ + µ − rate so that with the NP shift central values from CMS and LHCb [15] B can be reached.

This value should be compared with central SM values
for |V cb | = 0.042 and |V cb | = 0.040, respectively. Thus within 331 models, on the whole, the results for ∆F = 1 for |V cb | = 0.040 appear more interesting than for |V cb | = 0.042. As we will see below this is in particular the case for larger values of M Z . Finally, we show in Fig. 8 the maximal value of ∆(ε /ε) for |V cb | = 0.040 as a function of |V ub |. We observe that this value rises approximately linearly with increasing |V ub | and for M Z = 3 TeV and |V ub | = 0.0044 that is consistent with the inclusive determinations could reach values as high as 7.7 × 10 −4 . This possibility should be kept in mind even if such high values of |V ub | seem rather unlikely as stated above. For |V cb | = 0.042 the effects of changing |V ub | turn out to be smaller.     We will next investigate what happens when higher values of M Z , outside the reach of the LHC together with CKM parameters in (13), are considered. As an example we set M Z = 10 TeV. In Fig. 9 we demonstrate that also in this case the tension between ∆M d and ε K can be easily removed. In the case of ∆M s 331 models perform much better than the SM represented by the red point so that the inclusion of the uncertainty in F Bs B Bs in (15), can bring easily 331 models to agree with data which is not possible within the SM. The question then arises what happens with NP effects in other observables for such high values M Z . On the basis of our discussion in Section 4 in [23] we expect the effects in B s → µ + µ − and C 9 to be smaller than for M Z = 3 TeV, which can be confirmed as seen in Tables 1 and 2. On the other hand ε /ε was found in [23] to be significantly enhanced for M Z = 10 TeV as can be seen in Fig. 6 of that paper. Moreover through renormalization group effects it could be even enhanced slightly more than for M Z = 3 TeV. However, as seen in Fig. 10, with new lattice constraints, this is no longer the case and the maximal allowed shifts in ε /ε are below 1 × 10 −4 , far too small to remove ε /ε anomaly.
In order to understand this drastic change we recall the general formula for ε /ε for arbitrary M Z in 331 models in the absence of Z − Z mixing which is irrelevant in M8, M9 and M16 [23] In this manner the M Z suppression in (26) is compensated and the fate of ε /ε depends on the allowed values of sin(δ 2 − δ 1 ) that follow not only from ∆M s,d constraints but also from S ψK S and S ψφ constraints. Our analysis shows that whereas in our previous analysis values of sin(δ 2 − δ 1 ) = 1 were allowed this is no longer the case after new lattice results and maximal values of sin(δ 2 − δ 1 ) are significantly below unity. While this suppression appears to be roughly compensated by the increase of the products 13s23 3 for M Z = 3 TeV, this is no longer the case for M Z = 10 TeV and ε /ε is strongly suppressed. This feature is clearly seen in Fig. 11, where the old ranges are the ones from [23] and the new ones found here.

|V cb | = 0.040
We next consider the CKM input in (23). The results for M Z = 10 TeV are shown in Figs. 12-14. We observe: • The tensions between ∆M s,d and ε K can be much easier removed than for M Z = 3 TeV because of the increased NP effects in ε K . Comparing Fig. 12 with Fig. 8 we also observe that the agreement with data is better for |V cb | = 0.040.
• The upward shift in ε /ε up to (6−7)×10 −4 in now possible so that ε /ε with |V cb | = 0.040 can probe much higher mass scales than it is possible for |V cb | = 0.042 because of other constraints.
• The plots in Fig. 14 when compared with those in Fig. 11 explain why the NP effects in ε /ε for |V cb | = 0.040 have a different structure than for |V cb | = 0.042. sin(δ 2 − δ 1 ) can for |V cb | = 0.040 reach unity even for M Z = 10 TeV, while this is not possible for |V cb | = 0.042.
As seen in Tables 1 and 2 NP effects in B s → µ + µ − and C 9 are suppressed for M Z = 10 TeV but not as much as for |V cb | = 0.042.
All our results for M8 for different values of |V cb | and M Z are summarized in Table 1. Very similar results are obtained for M9. The corresponding results for M16 are summarized in Table 2. These tables show again how important is the precise determination of |V cb | in tree-level decays.
Finally the red curve in Fig. 8 demonstrates that ∆(ε /ε) for |V cb | = 0.040 and M Z = 10 TeV can for large values of |V ub | reach values in the ballpark of 8.8 × 10 −4 . This increase relative to M Z = 3 TeV is related to renormalization group effects as discussed in detail in [23].

Preliminaries
We will next look at the simplest 331 model recently proposed in [27] in which β = 0. We will denote it by M0. Even if this model fails to remove most of the anomalies in question, its simplicity invites us to have a closer look at its flavour structure. We will list Z and Z couplings in this model and present formulae for C NP 9 and C NP 10 as well as ε /ε. The expressions for ∆F = 2 processes and for B s → µ + µ − as functions of the couplings listed below can be found in [20][21][22] and we will not repeat them here. One only has to set β = 0 in that formulae. In this manner, in contrast to [27], we take Z − Z mixing in all observables automatically into account.

Z Couplings
Setting β = 0 in (17) of [21] we find for quark couplings with v ij given in (6). The diagonal couplings given here are valid for the first two generations of quarks neglecting very small additional contributions [20]. For the third generation there is an additional term which can be found in (63) of [20]. It is irrelevant for our analysis of FCNCs but plays a role in electroweak precision tests [22]. For β = 0 the diagonal b quark couplings differ from d and s couplings only by sign: Setting β = 0 in (18) of [21] we find for lepton couplings These definitions also apply to other leptons and quarks. All these couplings are evaluated for µ = M Z with g 2 = 0.633 for M Z = 3 TeV.

Z Couplings
The flavour non-diagonal couplings to quarks are generated from Z couplings through where using the general formula (10) in [22] we find for β = 0 Here with the scalar triplets ρ and η responsible for the masses of up-quarks and down-quarks, respectively. Thus for tanβ = 1 the parameter a = 0 and the Z − Z mixing vanish in agreement with [27]. On the other hand in the large tanβ limit we find a = −1 and in the low tanβ limit one has a = 1. The flavour diagonal Z couplings are the SM ones and collected in [21]. We evaluate them with g 2 = 0.652 and sin 2 θ W = 0.23116 as valid at µ = M Z .

C NP 9
and C NP

10
The corrections from NP to the Wilson coefficients C 9 and C 10 that weight the semileptonic operators in the effective hamiltonian relevant for b → sµ + µ − transitions are given as follows As seen in these equations C NP 9 involves leptonic vector coupling of Z while C NP 10 the axialvector one. C NP 9 is crucial for B d → K * µ + µ − , C NP 10 for B s → µ + µ and both coefficients are relevant for B d → Kµ + µ − . Here with G F being the Fermi constant. The terms R V µµ and R A µµ are generated by Z − Z mixing and are given as follows For β = 0 we find then that do not depend on M Z except for logarithmic M Z dependence of g 2 . The numerical factors above correspond to M Z = 3 TeV. We observe then that in the absence of Z − Z mixing (a = 0), independently of the parameters v ij , the following phenomenologically successful relation holds. This should be contrasted with models M8, M9 and M16 for which we found [23] C The result in (43) differs also from which is close to one of the favourite solutions in which NP resides dominantly in the coefficient C 9 . Thus already on the basis of B physics observables we should be able to distinguish between the models M0, (M8,M9) and M16. However in the presence of Z − Z mixing the relation (43) does not hold. While this effect is small in C NP 9 , it can be large in C NP 10 , in particular for a = 1, when C NP 10 becomes very small and the suppression of the rate for B s → µ + µ − is absent. More interesting is then the case of a ≈ −1, corresponding to large tanβ, as then the simultaneous suppressions of C 9 through C NP 9 and of B s → µ + µ − rate through C NP 10 are stronger. We find then C NP 9 ≈ −0.5 C NP 10 , (a ≈ −1), that on a qualitative level is still a better description of the data than the results in (44) and (45). But the crucial question is whether the values of both coefficients are sufficiently large when all constraints are taken into account. Before answering this question let us make a closer look at ε /ε in this model.

Preliminaries
The analyses of ε /ε in 331 models with β = 0 have been presented by us in [22,23]. We want to generalize them to the case β = 0. Generally in 331 models we have with the ∆(ε /ε) resulting from tree-level Z and Z exchanges. Now, as demonstrated by us in [22], the shift ∆(ε /ε) is governed in 331 models by the electroweak (V − A) × (V + A) penguin operator with only small contributions from other operators. This result applies to both Z and Z contributions with the latter ones significantly more important as demonstrated in [22]. Here we would like to point out that this pattern is not valid for β = 0. Indeed as seen in (41) of [22] the important coefficient C 7 (M Z ) generated by tree-level Z exchange for β = 0 vanishes for β = 0 and consequently, Q 8 operator cannot be generated from Q 7 operator by renormalization group effects. Contributions of other operators are very small so that Z contributions to ε /ε can be neglected. This is directly related to the fact, as seen in (28e), that the diagonal right-handed couplings of Z to quarks vanish for β = 0. But for Z such couplings are present implying that tree-level Z exchanges can provide a shift in ε /ε.

Z Contribution
The inclusion of this contribution is straightforward as the only thing to be done is to calculate the shifts from NP in the functions X, Y and Z that enter the SM model contribution to ε /ε. One finds then [22] where g 2 = 0.652 and λ t = V td V * ts . Replacing then the SM functions X 0 (x t ), Y 0 (x t ) and Z 0 (x t ) by in the phenomenological formula formula (90) for ε /ε in [8] allows to take automatically the first two contributions in (47) in 331 models into account. Inserting sin ξ in (33) and the Z coupling into (49) we obtain for β = 0 Evidently for a = 0, as done in [27], NP contributions to ε /ε from Z exchanges vanish and as the ones from Z can be neglected, ε /ε is full governed by the SM contribution. This appears presently problematic in view of the findings in [7][8][9]38] that a significant upwards shift ε /ε of at least 5 × 10 −4 is required to bring the theory to agree with the data from from NA48 [10] and KTeV [11,12] collaborations. The question then arises whether including Z − Z mixing we can obtain the required positive shift in ε /ε. But, as seen in (46), in order to preserve at least partly the pattern in (43) we are interested in − 1 ≤ a < 0, tanβ > 1 .
In order to answer this question we insert (51) into (90) in Appendix B of [8] to obtain first From Table 5 in [8] we find then for the central value of α s (M Z ): where (55) Now the results from the RBC-UKQCD collaboration imply the following values for B Using then the results in (56) we find and finally The question then arises whether for a in the range (52) one can get sufficient shift in ε /ε while satisfying other constraints. In particular the ones from ∆F = 2 transitions, where as seen in the previous section the SM experiences tensions in its predictions for ∆M s,d and ε K . In Fig. 15 we show ∆M s,d vs. ε K in M0. Red dots represent central SM values and black dots the central experimental values. We observe that the tensions between ∆M s,d vs. ε K present in the SM can be removed in the M0 model. But as seen in Fig. 16 the shift in ε /ε can be at most 1.1 × 10 −4 which is far too small to be able to remove ε /ε anomaly. Moreover, this maximal shift can only be obtained for the maximal Z − Z mixing. We have checked using the expressions in [22] that then the fit to EWPO is significantly worse than the one  in the SM, whereas the three models analysed in the previous section perform better in these tests than the SM [22].
We do not show the results for C 9 and B s → µ + µ − as NP effects are significantly smaller than in M8 and M16. Thus even if M0 can remove the tensions between ∆M s,d and ε K , it fails badly in the case of other anomalies and therefore cannot compete with models M8, M9 and M16, unless all anomalies disappear one day.

Summary
Motivated by the recently improved results from the Fermilab Lattice and MILC Collaborations on the hadronic matrix elements entering ∆M s,d in B 0 s,d −B 0 s,d mixing [5] and the resulting increased tensions between ∆M s,d and ε K in the SM and generally CMFV models [4], we have performed a new analysis of 331 models. In order to illustrate the sensitivity of the results to the modification of hadronic parameters we have first used the CKM input of our previous analysis in [23] that is given in (13). In addition, in order to illustrate the sensitivity of the results to the value of |V cb | we have also performed the analysis with the CKM input in (23), where |V cb | is lower than in (13). We also investigated |V ub | dependence.
The most important results of our analysis, summarized in Tables 1 and 2 are as follows: • The tensions between ∆M s,d and ε K can be removed in the three 331 models with β = 0 (M8, M9, M16) considered by us and this for both CKM inputs. This turns out to be also possible in the model with β = 0 (M0) in the case of the input (13) but it is much harder when |V cb | is smaller as in (23).
• Models M8, M9 and M16 can provide a positive shift in ε /ε up to 6 × 10 −4 for M Z = 3 TeV for both choices of |V cb | and |V ub | = 0.0036. But in contrast to our previous analysis this shift decreases fast with increasing M Z in the case of |V cb | = 0.042 but its maximal values are practically unchanged for M Z = 10 TeV when |V cb | = 0.040 is used. We also find that for |V cb | = 0.040 and the inclusive values of |V ub | the maximal shifts in ε /ε are increased to 7.7 × 10 −4 and 8.8 × 10 −4 for M Z = 3 TeV and M Z = 10 TeV, respectively. In the model M0, in which NP contribution to ε /ε is governed by Z − Z mixing, NP effects are very small even for M Z = 3 TeV.
• In M8 and M9 the rate for B s → µ + µ − can be reduced by at most 10% and 20% for M Z = 3 TeV and |V cb | = 0.042 and |V cb | = 0.040, respectively. This can bring the theory within 1 σ range of the combined result from CMS and LHCb and for |V cb | = 0.040 one can even reach the present central experimental value of this rate (24). The maximal shifts in C 9 are C NP 9 = −0.1 and C NP 9 = −0.2 for these two |V cb | values, respectively. This is only a moderate shift and these models do not really help in the case of B d → K * µ + µ − anomalies.
• In M16 the situation is opposite. The rate for B s → µ + µ − can be reduced for M Z = 3 TeV for the two |V cb | values by at most 3% and 10%, respectively but with the corresponding values C NP 9 = −0.3 and −0.5 the anomaly in B d → K * µ + µ − can be partially reduced.
Therefore our analysis demonstrates that in the presence of the anomalies discussed by us the U (1) X factor in the gauge group of 331 models cannot be U (1) Y .
• For higher values of M Z the effects in B s → µ + µ − and B d → K * µ + µ − are much smaller. We recall that NP effects in rare K decays and B → K(K * )νν remain small in all 331 models even for M Z of few TeV.
Even if models M8, M9, M16 still compete with each other and M0 does not appear to be phenomenologically viable from present perspective, our feeling is that eventually only models M8 and M9 have a chance to survive future tests if the anomalies discussed by us will be confirmed in the future. The point is that with present theoretical uncertainties in B d → K * µ + µ − NP effects, even in M16, will be hardly seen in this decay. The decay B s → µ + µ − is much cleaner and in the flavour precision era 15 − 20% effects from NP, which are only possible in M8 and M9, could in principle be distinguished from SM predictions but this would require very large reduction in the experimental error on its rate.
Thus the main virtue of 331 models as opposed to SM and CMFV models is the ability to remove the tensions between ∆M s,d and ε K and simultaneously provide a significant upward shift in ε /ε but only for lower values of |V cb | can this property remain for M Z beyond the LHC reach. The possibility of a significant suppression of the rate for B s → µ + µ − in M8 and M9 for |V cb | = 0.040 is also a welcome feature. In particular, as it is correlated with the maximal shift in ε /ε.
While the NP pattern in flavour physics identified by us in 331 models is interesting, we should hope that eventually NP contributions to flavour observables are larger than found in these models and are also significant in rare K decays which are theoretically very clean and in B → K(K * )νν which are cleaner than B → K(K * )µ + µ − decays. Most importantly the comparison of our results in [23], prior to the lattice results in [5], with the ones obtained using this new input demonstrates clearly how the shifts and increased accuracy in non-perturbative parameters can have important impact on the size of NP effects. Similar comment can be made in connection with |V cb |.