Limiting FCNC induced by a CP symmetry of order 4

CP4 3HDM is a three-Higgs-doublet model based on the $CP$ symmetry of order 4 (CP4). Imposing CP4 leads to remarkable connections between the scalar and Yukawa sectors and unavoidably generates tree-level flavor-changing neutral couplings (FCNC). It remains unclear whether FCNC can be sufficiently suppressed in the CP4 3HDM. In this paper, we systematically explore this issue. We first develop an efficient scanning procedure which takes the quark masses and mixing as input and expresses the FCNC matrices in terms of physical quark observables and quark rotation parameters. This procedure allows us to explore the FCNC effects for all the Yukawa sectors possible within the CP4 3HDM. We find that, out of the eight possible CP4 Yukawa sectors, only two scenarios are compatible with the $K$, $B$, $B_s$ and, in particular, $D$-meson oscillation constraints. The results of this work serve as clear guidelines for future phenomenological scans of the model.


Contents
1 Introduction symmetry options for model building [27][28][29], including several non-abelian finite groups, as well as novel options for CP symmetries and CP violation. In fact, the multi-Higgs activity in late 1970's and early 1980's was mainly driven by a trial-and-error search for symmetry group and group representation assignments which would explain the emerging patterns of the quark masses and mixing parameters. Although no ideal candidate with an exact symmetry group was identified due to the reasons eventually explained in [30][31][32], 3HDMs with softly broken symmetries may do the job and remain today an actively explored direction [5,33].
In the 3HDM, the magnitude of the quark-Higgs FCNC can be well controlled, too. Already in the original proposal by Weinberg [4], the tree-level FCNCs were absent due to the global Z 2 × Z 2 symmetry imposed. Generic features of multi-Higgs-doublet models equipped with natural flavor conservation were analyzed in [34], and many other works on specific 3HDMs models followed, including the so-called democratic 3HDM [35]. Recent examples of multi-Higgs-doublet models which are free from FCNC can be found in [36,37]. Non-vanishing tree-level FCNC can also be tolerated if brought under control by additional symmetries, for example, in the spirit of the BGL model [38,39] or via alignment of Yukawa matrices [40].

CP4 3HDM: the present status
Historically, multi-Higgs-doublet models emerged from the search for new opportunities for CP violation, [1,4,41]. But they can also accommodate new forms of CP symmetry, which are physically distinguishable from the traditional CP . It was noticed long ago that the action of discrete symmetries, such as CP , on quantum fields is not uniquely defined [41][42][43][44]. If a model contains several fields with identical gauge quantum numbers, one can consider a general CP transformation (GCP) which not only maps the fields to their conjugates but also mixes them. For example, a GCP acting on complex scalar fields φ i , i = 1, . . . , N , can be written as [45,46] φ i (r, t) CP − − → CP φ i (r, t) (CP) −1 = X ij φ * j (−r, t), X ij ∈ U (N ) . (1.1) The conventional definition of CP with X ij = δ ij is only one of many possible choices and is, in fact, basis-dependent. If a model does not respect the conventional CP but is invariant under a GCP with a suitable matrix X, then the model is explicitly CP -conserving [41]. The presence of the matrix X has consequences. Applying the GCP twice leads to a family transformation, φ i → (XX * ) ij φ j , which may be non-trivial. One may find that only when the GCP transformation is applied k times that one arrives at the identity transformation; thus, a CP symmetry can be of order k. The conventional CP is of order 2; the next option is a CP symmetry of order 4, denoted CP4. Since the order of a transformation is basis-invariant, a model based on a CP4 represents a physically distinct CP -invariant model which cannot be achieved with the conventional CP .
Within the 2HDM, imposing CP4 on the scalar sector always leads to the usual CP [21]. In order to implement CP4 and avoid any conventional CP , one needs to pass to three Higgs doublets. This model, dubbed CP4 3HDM, was constructed in [47] building upon results of [48] and was found to possess remarkable features which sometimes defy intuition [47,49,50]. If CP4 remains unbroken at the minimum of the Higgs potential, it can protect the scalar dark matter candidates against decay [51] and may be used to generate radiative neutrino masses [52]. Multi-Higgs-doublet models based on even higher order CP symmetries were constructed in [53].
CP4 symmetry can also be extended to the Yukawa sector leading to very particular patterns of the Yukawa matrices [54]. The CP4 transformation, by construction, mixes generations; therefore, in order to avoid mass degenerate quarks, the explicit CP4 symmetry of the model must be spontaneously broken. Then, the Yukawa sector contains enough free parameters to accommodate the experimentally measured quark masses and mixing, as well as the appropriate amount of CP violation. The numerical scan of the scalar and Yukawa parameter spaces performed in [54] yielded CP4 3HDM examples which satisfied theoretical constraints, the electroweak precision tests, and did not violate the kaon and B-meson oscillation properties.
Tree-level Higgs-mediated FCNCs are unavoidable in the CP4 3HDM. The scalar alignment assumption used in [54] guarantees that the SM-like Higgs h SM does not change quark flavor. But the other neutral scalars certainly do. The scan performed in [54] produced parameter space points which satisfied the kaon and B-meson oscillation parameters. However it remained unclear how exactly it happened: was it due to cancellation among different Higgs contributions or was it a consequence of intrinsically small FCNCs? Also, D-meson oscillation constraints were not included in [54].
A few years later, Ref. [55] revealed that the charged Higgs boson couplings to top quarks ruled out the vast majority of parameter space points believed viable in [54]. Almost all examples found in [54] contained one or two charged Higgs bosons lighter than the top quark, which opened up non-standard top decay channels. Ref. [55] found that almost all these points conflict with the LHC Run 2 data on searches for these top decays or with the total top width. Only a couple of points survived due to a very peculiar structure of their Yukawa matrices. Once again, it remained unclear how these features emerged, what parameters controlled them, and whether one could generate other benchmark models with custom-tailored charged Higgs coupling matrices.
Normally, one would want to repeat the parameter space scan including the new constraints. However, an additional downside of [54] was that only a very small fraction of random scan points passed all the flavor constraints. This was due to the intrinsically inefficient scanning method, in which one begins with a random point, finds quark masses and mixing parameters way off their experimental values, and then iteratively searches for a suitable combination of parameters which would agree with quark masses and mixing values.

The goals of the present work
After Ref. [54], CP4 3HDM emerged as a viable model with remarkable cross-talk between the scalar and Yukawa sectors, which accomplishes unexpectedly much for a 3HDM equipped with a single symmetry. At the same time, [54] left many questions unanswered. What is the typical magnitude of the FCNC in each class of Yukawa models? How small the FCNCs can in principle become in each Yukawa sector, without compromising quark masses and mixing? What parameters control their smallness? Are there robust structural predictions for the amount of FCNC or for correlations among various entries? Is there any efficient cancellation among different FCNC contributions to the meson oscillation parameters?
In order to answer all these questions, a new analysis of the CP4 3HDM phenomenology is definitely needed, which must draw lessons from the above results. This new analysis should be based on a more efficient scanning procedure in the Yukawa parameter space, in which one takes the physical quark sector parameters as input and then, using the remaining freedom, arrives at the FCNC matrices. We call this approach the inversion procedure. Using it, we will be able to investigate the typical values of the off-diagonal FCNC entries as well as the minimal values one could in principle achieve.
These are the main goals of the present work. We will explicitly construct the inversion procedure for each type of the CP4 invariant Yukawa sector. Using it, we will explore the magnitude of various FCNC couplings and check whether they are compatible with the meson oscillation constraints for neutral kaons, B, B s , and D-meson systems. To stay focused, we deal with quarks only, although leptons can be included in a similar fashion.
We stress that, in this work, we do not aim at a full parameter space scan of the model including the scalar sector. Before such a scan can be undertaken, one first must develop a new procedure and clarify the above issues on FCNC pattern. The objective of this work is to understand the structural features of the Yukawa sector. The results of this work, both the new procedure and the insights on the FCNCs, will guide the full phenomenological scan of the CP4 3HDM, which is delegated to a future publication.
The structure of the paper is the following. In the next Section, we will remind the reader of the CP4 3HDM scalar and Yukawa sectors and provide general expressions for the matrices N d and N u which describe coupling of neutral scalars with physical quarks. Section 3 describes the inversion procedure, namely, how one can recover the matrices N d and N u starting from the physical quark masses and mixing parameters. In Section 4, we apply this procedure for the three non-trivial cases of the CP4 invariant Yukawa sectors. Then, in Section 5, we describe a toy model which helps understand the minimal amount of FCNC effects and the parameters which control them. Equipped with all these results, we describe in Section 6 numerical results for FCNC in all CP4 3HDM Yukawa sectors. In Section 7, we discuss the grand picture which emerges from this analysis and present guidelines on how to perform, in future, a full parameter space scan which has a chance to satisfy all flavor constraints. The Appendix contains technical details on the inversion procedure.
2 3HDM with CP -symmetry of order 4

The scalar sector of CP4 3HDM
The 3HDMs make use of three Higgs doublets φ i , i = 1, 2, 3 with identical quantum numbers. CP4 acts on Higgs doublets by conjugation accompanied with a rotation in the doublet space. Following [47,54], we use the following form of the CP4: Applying this transformation twice leads to a non-trivial transformation in the space of doublets: 3 . In order to get the identity transformation, we must apply CP4 four times, hence the order-4 transformation. It is known that any CP-type transformation of order 4 acting in the space of three complex fields can always be presented in the form (2.1) by a suitable basis change [44]. The most general renormalizable 3HDM potential respecting this symmetry was presented in [47]. It can be written as with all parameters real, and with real λ 5 and complex λ 8 , λ 9 . Out of the two remaining complex parameters, only one can be made real by a suitable rephasing of φ 2 and φ 3 . However, in this work we prefer to keep λ 8 , λ 9 complex because we will rely on the residual rephasing freedom to make the vacuum expectation values (vevs) real. Also, two additional terms with coefficients λ 6 and λ 7 invariant under the same CP4 could be inserted in (2.3) but, without loss of generality, they can be eliminated with a suitable CP4 preserving rotation, see details in [54]. Minimization of this potential and the resulting scalar bosons mass matrices were studied in [54]. In order to avoid pathologies in the quark sector, CP4 must be spontaneously broken. The three doublets acquire vevs, which can in principle be complex, but the rephasing freedom allows us to render the three vevs real: With v = 246 GeV fixed, the position of the minimum is described by two angles β and ψ (in the expression above, we used the shorthand notation for sines and cosines of these angles). Spontaneous breaking of CP4 leads to four physically equivalent minima linked by the CP4 transformation. Since we are allowed to pick up any of them, we can safely assume that the angle ψ lies in the first quadrant.
In the standard minimization procedure, one begins with the scalar potential V and locates the minimum. For a phenomenological analysis, it is more convenient to choose angles β and ψ are free parameters. In this way, there will be two relations among the parameters of the potential, which could also be found in [54]. For example, knowing ψ allows us to relate the imaginary parts of λ 8 and λ 9 : s 2ψ Im(λ 8 ) + c 2ψ Im(λ 9 ) = 0.
Next, it is convenient to rotate the three doublets to a Higgs basis as (2.5) In this Higgs basis, the vev is located only in H 1 :  [56,57]. However if unsuppressed FCNC couplings leak into the SM-like Higgs interactions, they may immediately run in conflict with the meson oscillation constraints and with the Higgs-top FCNC measurements. What amount of scalar misalignment is still tolerable is a compex issue as it depends on the interplay between the FCNC matrices, the exact neutral scalars mixing matrix, and the SM-like Higgs interactions. We definitely plan to investigate these issues in future. What we do in the present work -exploring generic FCNC couplings allowed within the CP4 3HDM -can be seen as the starting point to this study. Thus, for now, we adopt scalar alignment, and this is all the information from the scalar sector we need in this work.

General FCNC matrices in 3HDM
To set up the notation, let us begin with the general expressions for the quark Yukawa sector in the 3HDM: Here, we use the notation of [58] extended to the 3HDM. The three generations of quarks are implicitly assumed everywhere, their indices suppressed for brevity. The superscript 0 for the quark fields indicates that these are the starting quark fields; when we pass to the physical quark fields by diagonalizing the quark mass matrices, we will remove this superscript. Using the real vev basis (2.4), we write the quark mass matrices as They are, in general, non-diagonal and complex. The interaction of the neutral (complex) scalars with the quarks can be described both in the initial basis and in the Higgs basis we chose; the relation between the two bases is given by For the up-quark sector, we obtain where As usual, these mass matrices are diagonalized by unitary transformations of the quark fields, 14) The same quark rotation matrices also act on the matrices N : and, similarly, for N d3 , N u3 . The four Yukawa matrices N d2 , N d3 , N u2 , N u3 are the key objects we study in this work. They describe the coupling patterns of the neutral complex fields H 0 2 and H 0 3 with the three generations of physical quarks. Their off-diagonal elements indicate the strength of FCNC. Within the 2HDM, we would only get N d2 and N u2 , see e.g. [58]. The additional matrices N d3 and N u3 arise in the 3HDM, and their patterns can be very different from N d2 and N u2 .
Notice that the complex fields H 0 2 and H 0 3 contain four real neutral degrees of freedom. All these components mix, and to get the physical Higgses, we would need to diagonalize the neutral scalar mass matrix. In the present work, we do not aim at a full phenomenological scan of the model; we only want understand how FCNCs can in principle be limited within the CP4 3HDM. For that purpose, it is sufficient and more transparent to work directly with N d2 , N d3 and N u2 , N u3 without invoking mixing among neutral bosons.

The Yukawa sector of the CP4 3HDM
CP4 symmetry can be extended to the Yukawa sector of 3HDM. This extension is not unique, but there is a limited number of structurally different options. This problem was solved in [54] and yielded four distinct cases, labeled A, B 1 , B 2 , and B 3 , separately in the up and down-quark sectors. For example, in the down-quark sector, case A represents the trivial solution which is completely free from FCNCs. Cases B 1 , B 2 , B 3 involve all three Yukawa matrices with a non-trivial structure: • Case B 1 : • Case B 2 : In all cases, all the parameters apart from g 33 can be complex. A similar set of cases is found in the up-quark sector.
Since the up and down sectors involve the same left-handed doublets, one can only combine cases A or B 2 in the up sector with A or B 2 in the down sectors, and cases B 1 or B 3 in the up sector with B 1 or B 3 in the down sector. This leads to eight possible pairings. However one of them, (A, A), should be disregarded: as shown in the Appendix, it is unable to generate the CP -violating phase in the CKM matrix. Thus, we are left with seven viable combinations for the CP4 invariant Yukawa sectors: Since case (A, A) is excluded, we arrive at the important conclusion that FCNCs are unavoidable in CP4 3HDM.
3 Controlling FCNC in the CP4 3HDM: the strategy

The inversion procedure
Fitting the quark masses and the parameters of the CKM mixing matrix is, in general, a non-trivial task for multi-Higgs-doublet models with generic Yukawa sectors. One begins with several Yukawa matrices Γ i and ∆ i , which usually have many free parameters, multiply them by vevs v i and sum them to produce the mass matrices M 0 d and M 0 u , see Eqs. (2.8). After the bidiagonalization procedure (2.13), we obtain physical parameters m q , V CKM , as well as the quark coupling matrices N d and N u , see Fig. 1. Figure 1. In a generic multi-Higgs model, one begins with Γ i , ∆ i , computes the mass matrices and diagonalizes them. This procedure, indicated by thin arrows, is usually irreversible. In certain models, one can perform inversion (thick light blue arrows), which allows one to pass directly from the quark properties to the FCNC matrices.
In general, the passage from and vevs, one cannot recover individual Γ i . Technically, it comes from the fact that Γ i are not linearly independent. Indeed, in the general case, for each set of v i , it is possible to modify Γ i by some matrices δΓ i which sum up to the zero 3 × 3 matrix: i v i δΓ i = 0 3 . Thus, the Yukawa sectors based on Γ i and on Γ i + δΓ i lead to the same mass matrix M 0 d . Another problem is that, if matrices Γ i , ∆ i are not generic, it is not guaranteed that they can reproduce the known quark masses and mixing parameters at all, especially when the vev alignment is also constrained by the scalar sector. A nice illustration of this situation for the A 4 and S 4 symmetric 3HDMs can be found in [59,60].
In the light of these difficulties, one is often forced to scan the multi-dimensional parameter space of the model in a way which is intrinsically inefficient. If one knows vevs and randomly selects Γ i , ∆ i , one obtains M 0 d and M 0 u , which lead to quark masses and mixing very different from their experimental values. One needs to repeat the scan many times, trying to iteratively approach the measured values.
However, in certain classes of models one can invert the above procedure: that is, knowing M 0 d and vevs, one can uniquely reconstruct individual Γ i . This situation takes place, for example, if each Γ i lives in a different subspace of the general 3 × 3 complex matrix space. In particular, in models with abelian symmetry groups with different scalars transforming with different charges, the non-zero entries of the Yukawa matrices Γ i or ∆ i are non-overlapping and thus satisfy the above conditions. The possibility of inversion significantly facilitates the phenomenological study of the model. Instead of a random scan over Γ i , ∆ i , one takes the physical quark masses and mixing parameters as input, parametrizes V dL , V dR , and V uR in a suitable way, and directly obtains a parameter space point which automatically agrees with the experimental quark properties. No parameter sets are wasted anymore. Moreover, one can express the physical quark coupling matrices N d and N u via quark masses and mixing as well as V dL , V dR , and V uR . In this way, one may derive certain predictions for FCNC or use some of their entries as additional input parameters instead of V dL , V dR , V uR . Such a procedure will bring more control over FCNC within a given class of models.
The inversion procedure is available not only in models with natural flavor conservation, such as the Type-I or Type-II 2HDMs, but also in the Branco-Grimus-Lavoura (BGL) two-Higgs-doublet model [12], which allows for small FCNCs controlled by the third row of V CKM . A similar inversion procedure was recently constructed in the U (1) × Z 2 -symmetric 3HDM [39], where the quark sector closely resembled the BGL model.

The inversion procedure in CP4 3HDM
The first scan of the CP4 3HDM parameter space reported in [54] was done in the traditional way, by randomly choosing Γ i , ∆ i and iteratively approaching the physical quark sector parameters. In the CP4 3HDM, the matrices Γ 2 and Γ 3 (or ∆ 2 and ∆ 3 ) live in the same subspace, so that, at the first glance, the inversion procedure is not expected to work. However the entries of the two matrices are related by the CP symmetry of order 4. Thus, as already mentioned in [54], the inversion procedure does apply to the CP4 3HDM, although it was not used in that scan.
In the present work, we develop this procedure for each case of the Yukawa sector (2.17)- (2.19). Starting from the physical quark parameters m q , V CKM and parametrizing the quark rotation matrices, we are able to calculate N d2 , N d3 , N u2 , N u3 for the physical quark couplings with the second and third Higgs doublets H 0 2 and H 0 3 in the Higgs basis. Although these scalars are not yet the mass eigenstates, the matrices provide a clear picture of the FCNC magnitude and patterns to be expected in the model.
As we will show below, certain features of the matrices N d2 and N u2 are universal for B 1 , B 2 , B 3 and closely resemble the BGL-like models or the U (1) × Z 2 -symmetric 3HDM [39]. In particular, we will see that their FCNC couplings are controlled by the third row of the V CKM , which makes them naturally small. However there are also important differences with respect to the original BGL models and the U (1) × Z 2 -symmetric 3HDM. One is that the symmetry transformation of our model is unavoidably mixes quark generations, leading to FCNC patterns which do not always follow the familiar behavior.
Another difference with respect to the BGL-type models is that the Yukawa matrices of CP4 3HDM are less constrained. As a result, the quark rotation matrices V dL , V dR , V uL , V uR , do not, in general, exhibit the block-diagonal structure characteristic of the BGL models.
However, since the CP4 transformation acts on three generations by the irreducible representation decomposition 2+1, it automatically singles out one quark generation, which we assume to be the heaviest one. As a result, we find it natural to explore matrices V dL , V dR , V uL , V uR in the vicinity of the block-diagonal form. With the inversion procedure implemented, we will study the generic magnitude and the patterns of the FCNC matrices in the CP4 3HDM and check what is the smallest FCNC one can achieve for the physical values of m q and V CKM .

Target values for the FCNC magnitude
How small do the off-diagonal elements of N d and N u need to be in order to have a chance to satisfy the meson oscillation constraints? The exact bounds depend not only on the masses of the new bosons but also on other details of the scalar sector, such as their mixing angles, partial cancellations among different scalars, as well as self-cancellation between the scalar and pseudoscalar Yukawa couplings for each neutral Higgs boson [13,16]. For example, a scalar S and a pseudoscalar A of close masses originating from the same neutral complex field (S + iA)/ √ 2 coupled to quarks tend to compensate each other's contribution to meson oscillation amplitude due to the relative i 2 = −1 factor [13]. In this work, we are not going to invoke these delicate mechanisms. Instead, we address the question: is it possible to achieve, within CP4 3HDM, sufficiently small FCNC couplings which would satisfy all the neutral meson oscillation constraints for a 1 TeV Higgs boson without relying on additional cancellation?
If we answer it in the affirmative, then there are good chances that a full phenomenological scan of the scalar and Yukawa sectors of the CP4 3HDM, with all the cancellations included, will identify viable points with reasonably heavy Higgs bosons. Conversely, if finding such examples turns out impossible, it means that the future phenomenological scan must either heavily rely of strong cancellations or push all new Higgs bosons to multi-TeV range. In this case, the chances to find a phenomenologically acceptable version of the explicitly CP4-invariant 3HDM will be bleak.
Let us now define the target parameters. Following [16], we first rewrite the coupling matrix of a generic real scalar S of unspecified CP properties with down quarks as Both A and B are hermitean matrices. For example, K 0 -K 0 oscillations place constraints on |a ds | = |A 12 | and |b ds | = |B 12 |, and so on. Notice that, for a very asymmetric matrix 21 |, we obtain |a ds | ≈ |b ds |. A similar construction for the up-quark sector allows us to constrain the (uc) elements with the aid of D 0 -D 0 oscillations.
We consider the off-diagonal FCNC elements acceptable for a 1 TeV scalar if they satisfy the following upper limits borrowed from [16]: Notice that in all cases, the upper bounds on |b| are comparable with the Cheng-Sher Ansatz √ m i m j /v, which would give 0.9 × 10 −4 , 6 × 10 −4 , 25 × 10 −4 , 2 × 10 −4 , respectively. For a lower mass of the scalar S, the upper limits on these couplings will decrease proportionally. We stress once more that these constraints are not to be taken as strict pass-or-fail bounds. They just indicate the typical FCNC values which, for a 1 TeV scalar, could be compatible with meson oscillations without resorting to additional cancellation mechanism. However models violating these constraints by a large factor will most likely fail the meson oscillation test even with possible cancellation effects included.
Below, we will check one by one all seven Yukawa sectors given in Eqs. (2.20) or (2.21). If a sector manages to produce "viable" points, that is, parameter sets which pass all the above constraints, we consider this case promising, generate a sample of these viable points and will include it in a future phenomenological study. If a sector is unable to produce points passing simultaneously all the constraints, we say the this sector is "ruled out". Our goal is to relate them, that is, to express N 0 d2 and N 0 d3 via M 0 d , which will then allow us to express N d2 and N d3 in terms of quark masses and quark rotation matrices.
We begin with N 0 d2 , Eq. (2.10), which can be expressed as Here, with 1 3 being the unit 3 × 3 matrix and P 3 = diag(0, 0, 1) being the projector on axis 3. After the mass matrix bidiagonalization, we get This leads to a simple expression for the FCNC matrix N d2 , which makes it clear that the off-diagonal elements of N d2 are fully controlled by the third row of the matrix V dL : A similar analysis holds for the up-quark sector, leading to These expressions are familiar from the BGL model and the U (1) × Z 2 -symmetric 3HDM studied in [39]. This is not surprising: with our definition of the Higgs basis, the second doublet H 2 of our model matches the second doublet of the 2HDM in the Higgs basis. For N d3 , which has no counterpart in the BGL model, we obtain where M ij stand for the elements of the mass matrix M 0 d . The properties of the CP4 transformation lead to the unusual feature of this expression, namely, that the vev alignment angle ψ disappears at the price of using The subscripts in these matrices indicate that P 4 is a projector onto the subspace (1, 2) together with an order-4 rotation and R 2 is just a reflection within this subspace. After the quark field rotations, we get Here, we used the fact that the diagonal matrix M d is real. We also defined The matrix P (dL) 4 can be recast in a more revealing form. First we write (P 4 ) ij = − ij3 , where epsilon is the fully antisymmetric invariant tensor satisfying for any unitary matrix V ∈ U (3). In particular, it holds for V = V dL . Now, understanding V as V dL , we can write P 4 as (4.11) is also constructed with the aid of the same third row of V dL . Finally, In summary, using the quark masses and V CKM as input and assuming any form of

Case B 2
The Yukawa matrices for case B 2 , Eq. (2.18), are similar to case B 1 , with the role of rows and columns interchanged. The analysis proceeds in a similar way, with V dL and V dR being exchanged.
with the same R 0 3 as in Eq. (4.2). After the quark field rotations, we obtain A similar expression holds for the up-quark sector: The differences with respect to case B 1 are: m d i instead of m d j and V dR instead of V dL (with the corresponding differences for the up sector). Switching from the left to the righthanded fields is important: V dR and V uR can be chosen independently from each other, while V dL and V uL used in case B 1 were not independent. The expression for N 0 d3 resembles Eq. (4.6): with the same P 4 and R 2 as in (4.7). After the quark field rotations, we get Here, as before, we defined For the up-quark sector, we have The structure of Yukawa matrices in case B 3 is different, Eq. (2.19). Let us begin with N 0 d3 : which can be expressed as using the same matrix P 4 as before and the same projector on the third axis P 3 = diag(0, 0, 1). A similar expression holds for the up sector. After the quark field rotations, we get .

(4.25)
As before, for any matrix A, we understand A L as V † L A V * L and A R = V T R A V R ; for example, (P (uR) 3 For N d2 , we can construct two different presentations thanks to the properties of the matrix Γ 1 : The latter form can be compactly written as which, with the aid of the same matrices as in (4.24), leads to A similar expression holds for the up sector: . (4.29)

A toy model
In order to see how small the FCNCs can in principle be, let us begin with a toy version of the CP4 3HDM, in which the true CKM matrix is replaced by a simplified blockdiagonal form parametrized with the single Cabibbo angle θ C . Let us also assume that all four quark rotation matrices V dL , V dR , V uL , V uR , which are input parameters within the inversion procedure, can also be chosen in the same block diagonal form: The parameters θ, α, ζ, δ, γ in each matrix are independent, subject only to V CKM = V † uL V dL . Under these assumptions, the expressions for the matrices N are dramatically simplified.
Let us begin with case B 1 . The matrices N d2 , N u2 become diagonal: which is well known from the original BGL model. As for N d3 , N u3 , we first notice that P where all parameters correspond to the matrix V dR . As a result, we obtain where all unlabeled angles refer to V dR . A similar expression holds for the up sector: with the unlabeled angles corresponding to V uR . We observe that, in this toy model, none of the neutral Higgses generates any FCNC involving the third generation quarks b and t, thus easily satisfying the constraints from B/B s -meson oscillations. FCNCs only appear between the first two generations and lead to |a ds | ≈ |b ds | ≈ m s 2v sin 2θ dR sin β = 2 × 10 −4 sin 2θ dR sin β , Comparing them with Eqs. (3.3a), (3.3d), we conclude that we need θ dR < 0.2 and θ uR < 0.03 to satisfy the meson oscillation constraints. Since these two angles are independent free parameters, it is technically possible to satisfy the constraints and even to eliminate FCNCs altogether by setting θ dR = θ uR = 0. Then, the only non-trivial feature of the matrices N 's is the peculiar pattern of their diagonal values: for example, H 0 3 couples tō uu proportionally to m c , not m u .
Within the toy model, case B 2 leads to exactly the same diagonal matrices N d2 and N u2 as in Eq. (5.2). The matrices N d3 and N u3 have the same form as in (5.4) and (5.5) up to transposition; the convention now is that the unlabeled angles correspond to V dL , V uL . However, these two matrices are no longer independent. In particular, it is impossible to select both θ dL = θ uL = 0 since they are related with via a sizable Cabibbo angle |θ dL + θ uL | ≥ θ C ≈ 0.22. Incidentally, the upper values on θ dL and θ uL computed for a 1 TeV scalar are barely compatible with θ C . Another way out is to select case A at least for one of the two sectors, which would eliminate FCNCs in that sector completely and would allow to keep the other angle small.
Interestingly, case B 3 turns out to be incompatible with the toy model. Indeed, limiting the quark rotation matrices to the block-diagonal form forces us to either set Γ 2 = Γ 3 = 0 or to assume sin β = 0. In either case, one would observe that M 0 d ∝ Γ 1 , which would force m d = m s , in conflict with experiment.

Implications for realistic models
The lesson we draw from the toy model is that it is possible to suppress the FCNC couplings involving b quark even for the realistic V CKM . To achieve this, we should keep all the quark rotation matrices as close as possible to the block-diagonal form (5.1).
To get some insight, let us choose case (B 1 , B 1 ) and suppose that V uR , V dR , and V uL are block-diagonal and the compute V dL = V uL V CKM . Then N u2 keeps the diagonal form as in Eq. (5.2), while N u3 is still given by the same expression (5.5) as in the toy model. In the down-quark sector, the third row of V dL is proportional to the third row of the CKM matrix: V dL,3i = V uL,33 V CKM,3i . For N d2 we can use the general expression (4.4) but replace V dL,3i with V CKM,3i . In this way one obtains for N d2 exactly the same FCNC Ansatz as in the BGL model, which is driven in our case by CP4. As for N d3 , we return to the general expression (4.12). This matrix is not of the block-diagonal form anymore, but, when evaluating its third column, we encounter m d k (R can be expressed via the third row of V CKM in the following peculiar way: Here, the symbol indicates that we omitted the inessential phase factor. It is interesting to note that, up to this phase factor, this column is independent of the quark rotation angles and, due to |V td | |V ts |, the largest FCNCs involving b-quarks come in the form of b → d, not b → s transition. Finally, all other elements of the matrix (N d3 ) i3 can include either m d or m s but not m b . Indeed, the term with m b corresponding to k = 3 in Eq. (4.12) forces j = 3, which we just considered.

Numerical results
Equipped with the general expressions and guided by the insights from the toy model, we are ready to address the main question formulated in Section 3.3: is it possible to find examples of the CP4 3HDM Yukawa sectors which satisfy all the meson oscillation constraints given in Eqs. (3.3)?
The first step of our numerical study is to implement the inversion procedure described in Section 3. We start with the experimentally known quark masses and the CKM matrix, for which we take the following values: These are not the latest measurements; we just borrowed the valued used in [54] to facilitate the comparison. Next, we choose the quark rotation matrices in such a way that, in the original basis before the quark mass matrices diagonalization, we recover M 0 d and M 0 u exactly of the type that is required for each Yukawa sector. For cases A, B 1 , and B 2 , this can be done analytically, while for case B 3 we resort to a numerical procedure. Details of this step for each of the four cases are described in the Appendix. The angles and phases of the quark rotation matrices represent the space in which we perform numerical scan. This scan can done in different ways. What we label below as a "full scan" corresponds to a random selection, with the uniform probability distribution, of all the rotation angles and phases within their full ranges. In a "restricted scan", we try to choose the quark rotation matrices as close as possible to the block-diagonal form (5.1). We do this by choosing the angles θ 13 and θ 23 from the interval [−θ max , θ max ], where θ max is some pre-defined value, such as π/100. In all plots, unless specified otherwise, the default value of the vev ratio parameter tan β = 1. The other vev alignment angle, ψ, does not appear in matrices N 's. We begin our numerical exploration with a detailed analysis of case (B 1 , B 1 ). We present in Fig. 2 the impact of the kaon oscillation constraints on the FCNC couplings. The scatter plot shows the values of |a ds | vs. |b ds | obtained in a full scan (left) and a restricted scan with θ max = π/100 (right). On both scatter plots, we show separately the values coming from N d2 (blue) and N d3 (red).
As we see, even in a full scan, we obtain many points which pass the kaon oscillation constraints Eq. (3.3a), indicated by the dashed box, for both N d2 and N d3 . Restricting the scan to nearly block-diagonal quark matrices further suppresses FCNCs (in the right plot, the blue points corresponding to N d2 shrink at zero). This is in accordance with the toy model expectations, where N d2 becomes a diagonal matrix, Eq. (5.2), while N d3 approaches the block-diagonal form (5.4). Due to (N d3 ) 21 (N d3 ) 12 , we obtain |a ds | ≈ |b ds |, which explains the straight line segment shape of the plot. The upper limit of the last plot (N d3 for small θ max ) is explained by the estimate (5.6), taking into account that sin β = 1/ √ 2. In short, the kaon oscillations constraints can be easily satisfied within case (B 1 , B 1 ).  In Fig. 3, we explore the FCNC effects in the up-quark matrices and compare them with the D-meson oscillation constraints Eq. (3.3d). The overall picture here is the same as for kaons, with the exception that only a few points fall inside the box of the allowed |a uc |, |b uc | values. This is not surprising: the estimate (5.7) shows that we need a rather small θ uR to pass the D-meson oscillation constraints. Unfortunately, with the procedure we use for building case B 1 , we cannot fully control the value of this angle. However we checked that the points inside the box indeed correspond to |θ uR | < 0.025. Thus, the D-meson oscillations constraints can also be satisfied.  In Fig. 4 we show that the kaon and D-meson constraints can also be satisfied simultaneously. Here, we only show the most challenging case: |b uc | coming from N u3 and |b ds | coming from N d3 . To increase the density of points, we run the scan with 10 5 points instead of 5000 points used in the previous plots. We see that such points do exist and correspond to small θ uR and θ dR . We also checked that increasing tan β allows for a mild additional suppression of the FCNC couplings and further increases the number of points simultaneously passing the two sets of constraints. Next, we look into the B-meson oscillation constraints. Fig. 5 shows the same sequence of plots for values |a db |, |b db | together with the B-meson constraints (3.3b). We also observed a similar picture (not shown) for the B s mesons, constraints being less tight. As expected, small θ max brings the matrices close to the block-diagonal form and, as a result, suppresses the b-quark FCNCs. Thus, B-mesons do not represent an obstacle once we keep the quark rotation matrices close to the block-diagonal form. In fact, we could find points which satisfy the B-meson oscillation constraints for scalar masses as low as 200 GeV. The above numerical results lead to the following overall situation for case (B 1 , B 1 ).
• The off-diagonal elements of the Higgs-quark coupling matrices N d2 , N d3 , N u2 , N u3 can be controlled through appropriate parameters of the quark rotation matrices.
• Different parameters control different FCNC couplings. B/B s -meson oscillation constraints for a 1 TeV scalar are easily satisfied by choosing the quark rotation matrices sufficiently close to the block-diagonal form. The kaon and D-meson constraints can be satisfied by the restricting the mixing angle θ 12 in the matrices V dR and V uR .
• D-meson oscillations place the strongest constraint for a generic scan. They must be included in a viable phenomenological analysis. Since the previous analysis [54] did not include them, all the (B 1 , B 1 ) points considered there as viable would, mostly likely, be ruled out by the D-meson constraints. Finally, we could also select a subset of points which pass the constraints from all four neutral meson systems. In Fig. 6 we show how these "viable" points are distributed on the D-meson vs kaon plane (left) and on the B s vs B-meson plane (right). This run corresponds to tan β = 7; increasing the value of tan β allowed us to get more viable points with respect to the previous scans. Curiously, the plots show that couplings of H 0 2 and H 0 3 to quarks are shaped by different meson oscillation constraints. Indeed, N d2 and N u2 (blue points) satisfy the kaon and D-meson constraints by large margin, and their main limitations come from B physics. For N d3 and N u3 (red points), we observe the opposite trend: B physics constraints play a minor role, and the strongest restrictions come from kaons and D-mesons. This is a manifestation of the structurally different forms of the corresponding matrices. Also, from the distribution of these points inside the boxes we see that some points are compatible with the Higgs bosons with masses of a few hundred GeV. In a future work, we plan to couple them with a scalar sector scan and build viable benchmark models based on (B 1 , B 1 ) Yukawa sector of the CP4 3HDM. Case B 3 , defined by the Yukawa matrices (2.19), is a peculiar one. As already mentioned, this case does not possess a limit that could correspond to the toy model considered in Section 5. Therefore, diagonalizing a realistic case B 3 mass matrix must involve all three rotation angles. As a result, the (ds) and (uc) entries of the matrices N will receive contributions proportional to the third generation quark masses. This unavoidable threegeneration mixing is especially important in the up-quark sector. The large top quark mass gives sizable contributions to (N u2 ) 12 and (N u2 ) 21 and leads to an unacceptably large D-meson oscillation amplitude.  Consider, for example the Yukawa model (B 1 , B 3 ), in which the down-quark sector displays distributions similar to what we just studied. In Fig. 7, left, we show the usual plot of |a uc | and |b uc | obtained from N u2 and N u3 in a full scan of this model. The points cover a sizable area extending up to values of 0.01. The box with the D-meson oscillation constraints (3.3d) is barely visible on this plot, but we checked that there exist points with either N u2 or N u3 falling inside the box. However there are no examples in which both N u2 and N u3 pass the constraints. This is illustrated in Fig. 7, right, where we compare D-mixing FCNC contributions from N u2 and N u3 in terms of the following combined quantities: This plot was computed for tan β = 1; for other values of tan β, it has the same shape but is stretched along one of the axes. We observe a sharp lower bound on a quadratic combination for c 2 and c 3 , which depends on tan β and which is always much bigger than the experimental constraints. In short, we could not find any point for the model (B 1 , B 3 ) which would come sufficiently close to satisfying the D-meson oscillation constraints. We conclude that CP4 3HDM with Yukawa sector of type (B 1 , B 3 ) cannot produce viable points. The scan reported in [54] claimed to find suitable parameter space points in case (B 1 , B 3 ) but those points would be ruled out by the D-meson oscillation constraints.
We also ruled out case (B 3 , B 3 ) on the same grounds; predictions for the D-meson oscillations are a way off from the experimental constraints.
For the model (B 3 , B 1 ), too, finding points which would pass all the meson oscillation constraints turned out impossible. The plots in Fig. 8 illustrate the problems we had. The left plot shows that, in stark constraint with case (B 1 , B 1 ) given in Fig. (6), it is extremely difficult to find a point which would pass constraints from kaons and D-mesons simultaneously, especially for H 0 3 . The right plot of Fig. 8 shows comparison of B and B s constraints; it is again N d3 which falls outside the box. With all meson constraints combined, we find that every point violates at least one constraint by a factor of few.   Since case B 2 admits the toy model limit, we expect to have many parameter space points with up-quark rotation matrices V uL , V uR close to the block-diagonal form. Within the procedure we use for case B 2 , described in Appendix, Section A.2, we cannot choose arbitrary rotation angles. However when performing a scan, we observed that |a uc | ≈ |b uc | for both N u2 and N u3 and that many points indeed satisfy the D-meson constraints. This is illustrated by Fig. 9, where we show the values of |b uc | calculated from N u2 and N u3 . We could find points with the FCNC contributions significantly smaller than the D-meson constraints Eq. (3.3d). Thus, case (A, B 2 ) offers many viable points. In a future analysis, we plan to update the procedure which could parametrically suppress the (uc) couplings in N u2 and N u3 .
In case (B 2 , A), we have the opposite situation: the D-meson constraints become irrelevant as there are no FCNCs in the up sector, while there arise problems for kaons and B-physics. We found that the kaon oscillation constraints do not represent a severe problem; many points passing them for N d2 and N d3 could be found. However the combined analysis of B and B s -meson represents a serious obstacle. In Fig. 10, we show a simultaneous check of the B and B s constraints on a full scan. We see that the FCNC couplings from N d2 are rather suppressed and can pass the constraints but N d3 fails by far. This check alone rules out case (B 2 , A).  Finally, case (B 2 , B 2 ) leads to problems both in the up and down-quark sectors. Just as in the previous case, the constraints coming from kaon oscillations can be satisfied in a restricted scan with θ max = π/100. However, the D-meson constraints are nearly impossible to satisfy for N u3 even within a restricted scan with θ max = π/100, see Fig. 11, left. Moreover, when one tries to suppress D-meson mixing by restricting the ranges of rotation angles, the matrix N d3 runs into a conflict with constraints for B or B s -mesons, as shown in Fig. 11, right. This, too, is similar to case (B 2 , A). We conclude that case (B 2 , B 2 ) is ruled out by the combination of the D-meson and B/B s -meson oscillation constraints.
It is interesting to highlight once again an important difference between cases (B 1 , B 1 ) and (B 2 , B 2 ), which we already discussed in Section 5. For (B 1 , B 1 ), working in the vicinity of the block-diagonal form of quark rotation matrices, we could simultaneously suppress the kaon and D-meson FCNC contributions. This was possible because those elements were governed by the right-handed quark rotation matrices V dR and V uR , which are independent from each other. As a result, we could generate many viable points in case (B 1 , B 1 ). For (B 2 , B 2 ), these FCNC elements are controlled by the left-handed matrices V dL and V uL , which are not independent. This leads to the clash between down-sector and up-sector FCNC matrices N d3 and N u3 and makes it impossible to find viable points.

Discussion and conclusions
Neutral meson oscillations put to severe test all multi-Higgs models that feature Higgsinduced tree-level flavor changing neutral couplings [9]. CP4 3HDM, the three-Higgsdoublet model based on the exotic CP symmetry of order 4, does contain FCNCs with peculiar CP4-driven patterns. Nevertheless, the first phenomenological scan of the model presented in [54] identified parameter space points which pass many experimental constraints including kaon and B/B s -meson oscillations. Later, these parameter space points were ruled out by the light charged Higgs searches [55], but the intrigue remained. Can the FCNC effects be strongly suppressed within CP4 3HDM? In which CP4 Yukawa sectors can it happen and what controls the magnitude of FCNCs? Are there viable examples of the model which pass all the neutral meson oscillation constraints for the new Higgs bosons not heavier than, say, 1 TeV?
In this paper, we systematically investigated all these questions. Following [54], we work in the scalar alignment regime of the CP4 3HDM, which guarantees that the 125 GeV Higgs couplings to quarks are as in the SM, and that the FCNCs can arise only from the heavy new scalars.
First, we constructed the inversion procedure for CP4-invariant Yukawa sectors. That is, we showed how to recover the parameters of the model starting from the physical observables, quark masses and the CKM matrix, and using appropriately shaped quark rotation matrices. This inversion procedure offers a huge advantage over [54] for a phenomenological scan of the entire parameter space of the model because, by construction, each point of the scan agrees with the experimental values of the quark properties.
Our second result is the collection of analytic expressions for the matrices N d2 , N d3 and N u2 , N u3 which describe how the new Higgses from the second and third doublets (in the Higgs basis we use) couple to the down and up-quarks. The off-diagonal entries of these matrices are the FCNC contributions which we want to suppress. These matrices are written in terms of physical quark parameters and quark rotation matrix elements, which offers clear insights which parameters control which FCNC elements.
Third, we numerically studied whether the FCNC couplings can simultaneously pass the constraints imposed by oscillations in the four neutral meson systems: kaons, D-mesons, B-mesons, and B s -mesons. Since the contributions depend on the new scalars masses, we computed the effects for the reference mass of 1 TeV.
In our study, we did not rely on any possible cancellation among several Higgs contributions to meson oscillations. We did this on purpose: we wanted to check whether the FCNC couplings of all neutral scalars can be simultaneously suppressed within any CP4 3HDM scenario without any help of an additional cancellation mechanism. Examples of such mechanisms do exist [13,16] but their effect is not robust as it depends on the details of the scalar sector. If a particular benchmark model leads to destructive interference of Higgs constributions, then the FCNC constraints could be satisfied with scalars masses below 1 TeV. Whether this indeed happens in specific models will become clear after the full phenomenological scan of the model is done, which is delegated to a future work.
In total, there are eight possible CP4-invariant Yukawa sectors. It turns out that most of them lead to FCNC contributions to meson oscillations which cannot be simultaneously kept small in the four neutral meson systems. Often, this comes from structural features of the Yukawa matrices, so that no choice of the free parameters could bring all FCNCs under control.
We find that, out of the eight CP4 3HDM Yukawa sectors, only two benchmark scenarios are capable of producing viable parameter space points passing all four meson constraints.
• Benchmark scenario (A, B 2 ), in which the down-quark sector is completely free from FCNC. The only constraints arise from the D-meson oscillations and can be easily satisfied as the magnitude of FCNCs can be parametrically suppressed.
• Benchmark scenario (B 1 , B 1 ), in which both up and down-quark sectors exhibit FC-NCs but their magnitudes are small if the quark rotation matrices are close to the block-diagonal form.
In both cases, we generated large samples of viable points. All other scenarios fail at least one of the meson constraint by a large margin; we found that they could become viable only for new Higgses as heavy as several TeV. We also noticed the very important role of the D-meson oscillation constraints, which are sufficient to rule out entire scenarios of the CP4 3HDM Yukawa sectors irrespective of their free parameters. Notice that the previous scan [54] did not include the D-meson check. We suspect that all the parameter space points found there would be ruled out by D-mesons.
Our results naturally lead to several follow-up studies. Focusing on either of the two benchmark scenarios, one can now perform a full phenomenological scan of the parameter space of the CP4 3HDM, including the scalar sector. Due to the intrinsic interplay between the scalar and Yukawa sectors, it will be a non-trivial exercise to satisfy simultaneously the experimental constraints on the new scalars and on the fermions. If such benchmark points are found, they can be used to make further collider predictions.
One can also include the lepton sector of the CP4 3HDM. A neutrino mass model based on CP4 was described in [52], but it relied on the unbroken CP4 symmetry. Whether spontaneously broken CP4 can lead to an acceptable charged lepton and neutrino sectors remains an open question.
In summary, the CP4 3HDM achieves remarkably much for a multi-Higgs model based on a single symmetry, and, as we show here, it can survive all the neutral meson oscillation constraints. There remains much to be studied in this exotic but viable model, which will be the subject of the forthcoming papers. , we find the parameters of the correcting matrices C dR and C uR . Since we do not correct the left-handed fields, the CKM matrix is preserved. Thus, we arrive at a viable model with the quark rotation matrices V dL = V (gen.) dL Case B 2 mirrors the previous case; the two conditions to be satisfied involve elements of the third column of matrix M 0 d : Clearly, they can be met through the same two-step procedure, which involves, at the second step, a correcting matrix in left-handed field space. The problem is that, in this case, the relation between V dL and V uL will be broken. Thus, it is desirable to satisfy the conditions in Eq. (A.6) with the aid of right-handed quark rotation only. To do so, we choose the correcting matrix C dR of the following form: In each line, we group terms with s θ and c θ , then divide the two lines, eliminating θ. After that, we regroup the resulting equation separating terms with e iβ and e −iβ : e iβ g 12 g * 33 e iα −g 33 e −iα −g 32 g * 23 e iα −g 13 e −iα = e −iβ g * 22 g * 33 e iα −g 33 e −iα −g * 32 g * 23 e iα −g 13 e −iα . (A.9) The solution to this equation on β exists, if and only if the absolute values of the two long brackets here are non-zero and equal. This equality gives a single equation on α, which can be represented in a compact form as a 1 cos(2α + ψ 1 ) + a 2 cos(2α + ψ 2 ) = c , (A.10) where a 1 = 2|y 1 ||z 1 | , ψ 1 = arg(y 1 z * 1 ) , with y 1 =g 12g * 33 −g 32g * 23 , z 1 =g 12g33 −g 32g13 , a 2 = −2|y 2 ||z 2 | , ψ 2 = arg(y 2 z * 2 ) , with y 2 =g * 22g * 33 −g * 32g * 23 , z 2 =g * 22g 33 −g * 32g 13 , c = |y 1 | 2 + |z 1 | 2 − |y 2 | 2 − |z 2 | 2 . Thus, the procedure for generating a viable case B 2 is as follows. We generate random g ij , compute a, b, c and check if the inequality (A.14) is satisfied. If not, we just interrupt evaluation and pick up another random point. Once we find a point which satisfies condition (A.14), we write the two solutions of Eq. (A.12) as Once the solution for α is found, we insert it into Eq. (A.9) and check whether the long brackets are non-zero. If they are zero, we drop the evaluation and pick another set of randomg ij . Once the long brackets are non-zero, we calculate β. The last step is to insert α and β into (A.8) to get θ. Thus, all parameters of the correcting matrix C dR are determined. Once again, we arrive at a viable model with the quark rotation matrices V dL = V (gen.) dL and V dR = C † dR V (gen.) dR . The same can be done in the up-quark sector without disturbing the CKM matrix.

A.3 Case A
The mass matrix of case A is proportional to Γ 1 , which has a very constrainted structure, see Eq. (2.16). It has 9 real degrees of freedom, much fewer than the other cases. In fact, there exists a transformation R which brings this matrix to the real form: The matrix Γ (r.) 1 is a general real 3 × 3 matrix without any further constraints. Thus, the mass matrix of case A can be produced from D d using arbitrary orthogonal matrices in the left and right-handed quark spaces: This implies V dL = RO dL and V dR = RO dR , solving the inversion problem for case A.
If we combine case A in the down sector with case B 2 in the up sector, then we first compute M 0 d , obtain V dL , which is non-generic, then compute V uL = V dL V † CKM , and follow the above two-step procedure to get the correct B 2 matrix. The same procedure is applied to the combination (B 2 , A). The only difficulty seems to arise for case (A, A), as we are forced to use special matrices V dL and V uL . However, we immediately conclude that in this case V † uL V dL = O T uL O dL is purely real and cannot match the experimentally known CKM matrix. This means that case (A, A) is unphysical as it misses the CP -violating phase of the CKM matrix.

A.4 Case B 3
In case B 3 , the mass matrix constructed from Γ i given in Eqs. (2.19) has a unique feature: its first 2 × 2 block is proportional to a unitary matrix: =ḡ · U , whereḡ = |g 11 | 2 + |g 12 | 2 , U ∈ SU (2) . (A.18) If one tries to follow the above two-step strategy and begins with a generic M (gen.) d , it will be impossible to reach this form by employing block-diagonal V dL and V dR . One must use generic SU (3) rotations.
Instead of trying to find the correcting matrices analytically, we resorted to the fully numerical procedure. Namely, we use the diagonal D d and parametrize the generic unitary V Since the simultaneous shift of all ψ i and all ϕ i does not change the result, we can set one of the rephasing angles to zero, for example ψ 3 = 0, arriving at 9 free parameters per matrix. Treating the 18 angles as independent free parameters, we start with a random seed point and numerically search for their values which satisfy the relations characteristic for case B 3 . With the multidimensional parameter space, this is done easily; the relative accuracy achieved is 10 −10 or better. For case (B 3 , B 3 ), we perform this numerical search simultaneously in the up and down sectors, making sure that the left-handed rotation matrices produce the CKM matrix. Once again, all the conditions are easily satisfied, and in reasonable time we can generate thousands of viable parameter space points for case (B 3 , B 3 ) with all quark masses and mixing matching their experimental values.