A method to explore flavor symmetries of the 3HDM and their implications on lepton masses and mixing

We present a method to identify symmetry groups of the Yukawa sector of the three-Higgs-doublet model and to determine the implication that the symmetry has on the lepton masses and mixing. The method can accommodate different hypotheses about the group representation assignments, and thus support the exploration of candidate symmetry groups. For one particular representation selection scheme we apply the computer-implemented method to scan all discrete groups of order less than 1035. It can be proven that none of these groups defines a flavor symmetry that implies masses and neutrino mixing angles consistent with the experimental lepton data, although several cases are found that are partially or approximately consistent.


Introduction
Many important questions still remain unanswered in the Standard Model (SM), which on the other hand agrees with experiment with great accuracy [1,2]. In this situation, various extensions of the SM have been proposed to address some unresolved questions, for example, the origin of the three generations of fermions, with such different mass hierarchies and radically different mixing patterns of leptons compared to quarks. In the lepton sector, which is the focus of this work, various attempts have been made to devise a theory that could predict the lepton masses and their mixing (e.g. see [3][4][5][6][7]). The experimental data on these quantities have gained better precision in recent years (e.g. see [8][9][10]), and it is getting increasingly hard to explain them satisfactorily. An often encountered approach to this problem is to look for some flavor symmetry of the lepton interaction Lagrangian that can, in a way consistent with the experimental data, predict the values of the lepton masses and their mixing angles. It is well known that a continuous nontrivial flavor symmetry in the lepton sector does not exist, which follows from the fact that the lepton masses are distinct [11]. Discrete symmetries remain possible in the Yukawa sector. However, due to Schur's first lemma, the acting of two inequivalent representations in flavor space (one on the lepton doublets, one on the lepton singlets) implies that any mass matrix is either proportional to the unit matrix or vanishes. One strategy to avoid this problem is to break the flavor symmetry explicitly in such a way that the charged lepton mass matrix and the neutrino mass matrix are separately invariant under two different subgroups of a larger discrete symmetry group G [12]. The subgroups are usually kept small, although in [13], groups G of order up to 1000 have been investigated. It was shown that the possible lepton mixing patterns then depend on how the two subgroups are embedded within G. Another approach is based on a Lagrangian including mass terms constructed with lepton, Higgs doublet and scalar flavon fields. Group-invariance of the terms is looked for by systematic probing all plausible representation assignments [14], where dynamic parameters (VEVs) affect the predicted mixing angles [15] [16].
The aforementioned explicit flavor symmetry breakdown in the Yukawa sector can be avoided in the presence of Higgs fields (causing spontaneous symmetry breaking) that transform under G. Two-Higgs doublet models (2HDM) [17] are obvious candidates to investigate that principle. In [18] and [19] non-Abelian groups that have two-and three-dimensional irreducible representations have been investigated. The lepton masses and the angles of the mixing matrix resulting from the assumed discrete symmetry were found, but the agreement with the experimental data was not satisfactory. Although the implied lepton masses of the model were nondegenerate, the obtained mixing matrices appeared to be always monomial. However, the results are different when one more Higgs doublet is added, as for example in the three-Higgs-doublet model (3HDM) [20,21]. This model contains neither flavons nor other additional fields and it allows the treatment of neutrinos as Dirac particles or as Majorana particles. We will assume that the Higgs doublets form a flavor vector transforming under a group G, as do the charged-lepton and neutrino flavor vectors. In this model ( [18], [22]), the mass-squared matrices are not affected by Schur's first lemma: the masses can be nondegenerate and the neutrino mixing nontrivial, even when the assigned 3D irreps are inequivalent and G is non-Abelian.
The 3HDM symmetry depends on how group representations are assigned to the flavor vectors. In order to either prove the existence of a proper symmetry or to exclude any symmetry, different hypotheses about the representation assignments must be tested. For this purpose we have built a computer tool that analyses discrete groups of arbitrary order. The theoretical foundation, the design of the tool and the initial outcomes are presented in this article.
In the following section we define the extension of the SM to 3HDM and describe the construction of the Yukawa Lagrangian. In Chapter 3 we define one particular scheme for the representation assignments and summarize the scan results for groups with order of less than 1035. Chapter 4 describes the method for finding those transformations (defined by finite subgroups of U(3)) that impose non-trivial relationships between lepton masses and mixing angles. In Chapter 5 we list alternative assignments of the representations, and directions for improvement of the exploration method.

Flavor symmetry of the 3HDM
The SM comprises one electroweak SU(2) Higgs doublet, that generates the masses of the gauge bosons and the fermions. The fermions that we study in this work are the three charged leptons electron, muon and tau, denoted l e , l µ and l τ , and three neutrinos, that we name ν e , ν µ and ν τ . The subscripts e, µ and τ are called lepton flavors.
Suppose for a moment that the three charged leptons have equal mass and also the three neutrinos have equal mass. Then the charged leptons would be indistinguishable and the same goes for the neutrinos. Any Lagrangian constructed from the lepton fields would be invariant under interchanging of charged fermions or interchanging neutrinos or, more generally, under any unitary transformation of the triplets (or flavor vectors) (l e , l µ , l τ ) T and (ν e , ν µ , ν τ ) T . We would have a trivial continuous U(3) flavor symmetry. If any two charged leptons (or any two neutrinos) have different masses, the Lagrangian will in general lose its U(3) flavor symmetry. We will see that for the Lagrangian of the SM we can apply the reverse argument: if the Yukawa term is invariant under certain continuous U(3) transformations of the flavor vectors then the masses of the charged leptons, and of the neutrinos, are degenerate. The key point here is that unlike the SM the 3HDM does allow nondegenerate lepton masses under specific discrete symmetries and, in addition, predict relations among mass and mixing quantities that can be compared against the experimental data.
In the following we consider two scenarios. In the first we treat all six fermions as Dirac particles. In the second we treat the three neutrinos as Majorana particles instead. A difference between the Dirac neutrinos and Majorana neutrinos is the way they acquire mass through the Higgs mechanism. As an extension of the SM the 3HDM has three electroweak SU(2) Higgs doublets (instead of only one), denoted Φ i , i=1,2,3. Doublet i contains an electrically neutral Higgs field φ 0 i and a negative Higgs field φ − i . The left-handed components of the leptons form electroweak SU(2) doublets L α = (ν αL , l αL ) T for each of three flavors α. Lepton masses arise from Yukawa couplings between the Higgs fields and the (left-and right-handed components of the) lepton fields.
To define flavor symmetry we first need to know the Lagrangian of the system responsible for the lepton masses and mixing. In this work we account for the Yukawa sector of the Lagrangian only, not for the vacuum alignments imposed by symmetry of the potential energy term of the Higgs bosons. The Yukawa Lagrangian L l of the charged leptons has 27 terms (plus their Hermitian conjugates), for the occurring combinations (i, α, β), compactly denoted Here summation over the lepton flavor indices α, β = e, µ, τ and over i = 1, 2, 3 is understood.
Yukawa matrix defining the couplings due to Higgs doublet Φ i . L αL = (ν αL , l αL ). If we treat neutrinos as Dirac particles, i.e. there exist right-handed neutrino fields, then the Yukawa Lagrangian of the neutrinos is (similar to L l ): If we assume that neutrinos are Majorana particles instead then separate right-handed neutrino fields are not needed. For each flavor α the charge-conjugate lepton doublet L c αL = CL T αL is a right-handed doublet. Then an effective dimension-five operator can be constructed to define the Yukawa Lagrangian of Majorana neutrinos as [23] where g and M are constants and h M ij is a 3 × 3 matrix of coupling constants. The number of terms of this Lagrangian is 81 (not counting the H.c. terms), which is the number of combinations (i, j, α, β).
After spontaneous electroweak symmetry breaking, the vacuum expectation values v i for the Φ i generate the mass matrices M l , M ν and M M as The mass-Lagrangian terms for the charged leptons, and for the Dirac neutrinos and the Majorana neutrinos will read where l L , l R , ν L , ν R are flavor vectors for the left/right-handed charged leptons and neutrinos, respectively. For example l L = (l eL , l µL , l τL ) T . We also define v = (v 1 , v 2 , v 3 ) T as the flavor vector of the Higgs doublets.
We are now ready to explore symmetries of the three Lagrangian terms; for that purpose we can leave out the Hermitian conjugate terms, so it is sufficient to consider the quantities l L , l R , ν L , ν R and v will each be assigned a three-dimensional representation of some finite group G. In the present study we restrict to irreducible representations; other choices of representations will be discussed later. All representation matrices must be unitary in order to conserve the total lepton number and to ensure that where G F is the Fermi coupling constant. We should remember that l l and ν l reside in the same SU(2) doublet and hence should be assigned the same representation.
As said in the introduction, we want to find such a discrete symmetry of the full Lagrangian of leptons that will give mass matrices of charged leptons and neutrinos, where several free parameters can be selected so, that after diagonalization, we obtain lepton masses, mixing angles and the CP symmetry breaking phases consistent with experimental observations. For this purpose, we should study the symmetry of two terms, 1) the Yukawa Lagrangian which, after spontaneous symmetry breaking, produces lepton mass matrices, and 2) the Higgs potential. In the present study, we will only deal with the symmetry of a well-defined Yukawa Lagrangian, leaving for further analysis the selection of a Higgs potential the form of which is more freely defined.
The first step is to identify distinct groups G, isomorphic to a U(3) subgroup, that have one or more unitary three-dimensional irreducible representations, and assign these representations to the flavor vectors so that some or all of the three transformed terms of Eq. (6) remain invariant under simultaneous operation of matrices of the irreducible representations A L (g), A lR (g), A νR (g) and A Φ (g) for all g in G and any vector l L , l R and v. Note that different flavor vectors can be assigned different representations of G.
As an example let us consider the Yukawa Lagrangian for the charged leptons (7). It is Since all matrix operators are linear, we are dealing with a set of |G| equations with 27 unknown complex numbers (the entries of the three Yukawa matrices h l i ). |G| is the order of G. So if Eq. (10) can be solved, it defines the h l i . Here we realize that if Eq. (10) is solved for a particular group element g then the system of equations obtained after transformation of all flavor vectors with the operators A L (g), A lR (g), A νR (g) and A Φ (g), will also be solvable for g, giving exactly the same solution(s) h l i . So if Eq. (10) is solved for g, it will be solved for g n for all n. In that case we can remove the mutually equivalent equations. Further, if Eq. (10) is solved for group elements g and h it will be solvable for g n h m too. More generally, we need only to take into account the generators of G, rather than all its elements. In a similar way we can define the set of equations to solve for h ν i and h M ij : where for Eq. (12) the number of unknowns is 81 (the entries of nine 3 × 3 matrices). A discrete symmetry of lepton interactions given by some group G requires that Eqs. (10) and (11) are satisfied for Dirac neutrinos and Eqs. (10) and (12) for Majorana neutrinos for all generators of G.

Solving the invariance equations
We take Eq. (10) as an example again. Working out the matrix multiplications it can be expressed as The Kronecker product gives a 27×27 matrix, and h l is the 27-dimensional vector built from the Yukawa matrices h l 1 , h l 2 and h l 3 , in that order, row-wise. If h l is an invariant eigenvector satisfying Eq.(13) for all generators of G then L l is G-invariant. The invariance equations for the terms L ν and L M are  (14) and (15) are equivalent to the (Clebsch-Gordan) tensor product decomposition equations [24] still assuming that all matrix operators are unitary. Eqs. (16) (or (17)) are fulfilled if and only if the 9-dimensional tensor product operator allows a decomposition containing at least one three-dimensional matrix operator. Eq. (18) requires that the 27-dimensional tensor product operator contains at least one three-dimensional matrix operator. To investigate which groups and which representation assignments can accomplish such invariance one could expect that it is necessary to explicitly solve the invariance equations for each instance. However, in the next section it will be shown that a significant part of those calculations can be skipped. The computational complexity of the calculations needs to be considered further because groups with |G| >∼ 1000 can have 8 or more generators and hundreds of inequivalent three-dimensional representations, leading to up to 10 million equation sets Eqs. (13)-(15) for a single group. We therefore put some restrictions to the analysis. Only irreducible representations are considered and a candidate group must have at least one faithful three-dimensional representation. Our precise selection criterion is that at least one of the irreducible representations assigned to a mass term is faithful. As mentioned, different selection criteria can be tested too.
From a scan of groups with |G|<1035, while applying the selection criteria, we find 749 groups that provide one or more solutions h l to Eq. (13), and in each of these cases also a solution h ν to Eq. (14). In the entire set of (over 6 million) solutions h l (or h ν ), 2130 are linearly independent. 216 groups provide one or more solutions to Eq. (15), resulting in 70 linearly independent vectors h M .
It follows from Eq. (4) that for given solutions h l , h ν and h M the mass matrices (and hence the mixing matrices) are functions of the VEVs. Since we do not know the absolute scale of the Higgs couplings nor the flavor assignments we have v2/v1 and v3/v1 as 4 real parameters, and we can at most determine mass ratios m i /m j of the leptons. Also any PMNS matrix can be determined up to permutation of its rows and columns and (for the Majorana case) up to a phase for two rows.
The low number of solutions that, despite the 4 parameters, imply nontrivial neutrino mixing is somewhat surprising. The majority of the obtained mixing matrices is monomial and exclude mixing. For none of the tested groups we find a symmetry that is entirely compatible simultaneously with the experimental masses and mixing data, neither when the neutrinos are assumed Dirac or Majorana particles. However, there are symmetries (e.g. for ∆(96) and S 4 ) that imply partially or approximately correct quantities.

Method to find relations among the lepton mass and mixing parameters
The selection and processing of groups is fully automated using the computer-algebra system GAP [25]. To determine which representation assignments would solve a particular invariance equation it is sufficient to observe the group's character table, which is readily provided by GAP. By doing so a tremendous amount of redundant computation is avoided; we therefore go into some detail of this point. The character of g ∈ G in representation A, denoted χ A (g), is defined as the trace of matrix A(g). The mapping χ A is called the character of A. Let A and B be representations of G, then defines the inner product of characters χ A and χ B . It can be proven that Let also H be an irreducible representation of G. Then χ A⊗B , χ H is the multiplicity of H occurring in the tensor product decomposition of A ⊗ B. In Eq. (16), can take a value from 0 to 3. This is the number of linearly independent solutions h l to Eq. The process of identifying 3HDM symmetry groups is schematically given in Algorithm 1. The first if-statement filters out the groups with |G| not divisible by 3. It can be proven that groups with 3D irreps (and those are the ones we consider) must have order divisible by 3. In the next statement the character table of G is retrieved, using GAP. Let us take group S 4 (the group of all permutations of 4 objects) as an example. S 4 has order |G| = 24. Its character table is shown in Table 1. The top row in the table lists how many group elements are member of a particular conjugacy class of the group. S 4 has 5 conjugacy classes and the list adds up to 24. The first conjugacy class has one member, which is the identity element of the group. The next 5 rows list the vectors of characters for all irreducible representations of S 4 . Each row defines the 5 entries of a character χ i , where each entry equals the trace of a representation matrix. For example χ 3 = (2, 0, −1, 2, 0). Its 4th component specifies that the three group elements in the 4th conjugacy class each correspond to a representation matrix with trace equal to 2. The first numerical column (belonging to the identity element) indicates the dimensions of the matrix representations, since the trace of any identity matrix equals its dimension.
In the character table of S 4 the algorithm finds two three-dimensional irreducible characters namely χ 4 and χ 5 , listed in the last two lines of 1. So the set E (introduced in line 5 of the algorithm) is going to contain these two characters. The algorithm also signifies that both are faithful, since all their character entries (except the first) are smaller than 3, and thus none of the corresponding representation matrices will be an identity matrix. This guarantees that any two different group elements will be represented by two different matrices.
where i = 1, . . . , n G and n G is the number of generators, then BaseFixedSpace takes the matrices K l i as input and returns their n 1 -dimensional common eigenspace. In our application n 1 =1, so the eigenspace defines h l .
In line 13 M l is calculated as a function of v, using Eq. (4). From line 15 on it is checked whether assignment (A, C, H) can solve decomposition Eq. (17) and if so, the Dirac neutrino mass matrix M ν is calculated. Finally, if both M l and M ν are available, U PMNS is calculated (using Mathematica [27]) for the duplet ((A, B, H), (A, C, H)) using the corresponding representations ((3 A , 3 B , 3 H ), (3 A , 3 C , 3 H )). For group S 4 it turns out that an invariant U PMNS exists, as will be described later.
For a given (A, B, C, H), Algorithm 1 outputs either nothing or M l or M ν or both mass matrices, and in the latter case U PMNS as well. From this output the invariance of the individual terms L l or L ν can be studied, as well as the invariance of (L l +L ν ). Since U PMNS is extracted from observations of SU(2) L -symmetric charged current interactions the characters A and H both appear twice in duplet ((A, B, H), (A, C, H)); l αL and ν αL both reside in the same SU(2) L doublet and thus l L and ν L should transform equally. Likewise the transformations of v and v differ by complex conjugation, since the two corresponding Higgs doublets differ by complex conjugation.
In case the neutrinos have Majorana nature, the computations are similar to Algorithm 1. The character assignments are then of the form ((A, B, H), (A, H)).
From M l , M ν (or M M ) and U PMNS the 3HDM predicts lepton mass ratios and mixing angles as functions of v 2 /v 1 and v 3 /v 1 . As an example we take group S 4 again, which happens to be the smallest group allowing nontrivial flavor mixing in case neutrinos are Dirac particles. The two representation assignments that provide such neutrino mixing symmetry are ((3 1 , 3 1 , 3 1 ), (3 1 , 3 1 , 3 1 )) and ((3 1 , 3 1 , 3 2 ), (3 1 , 3 1 , 3 2 )), where 3 1 and 3 2 are the two inequivalent three-dimensional representations of S 4 . The first representation assignment leads to anti-symmetric mass matrices and a monomial mixing matrix. The second implies the mass matrices: where c ν is an arbitrary constant. Based on these two mass matrices, we obtain fits to the experimental neutrino mixing data with sin 2 θ 12 and sin 2 θ 23 ending up approximately 3σ larger than the observed values; sin 2 θ 13 and δ CP deviate less than 1σ from the observed values. Fig. 1 shows mixing angle distributions calculated at uniform sampling within a limited box of the four-dimensional parameter space. Points within 1σ from the experimental sin 2 θ 13 and δ CP are red colored. The best fit is obtained for the red points closest to the intersection of the shaded areas in the left-most subplot. Higher-order groups, such as C i × S 4 and C 4 × A 4 give similar distributions.

Conclusions
We have presented a method to determine discrete symmetries for the Yukawa sector in the 3HDM model, which gives the relationship between the masses of leptons and the elements of the PMNS mixing matrix. The implied relations among the lepton masses and mixing can be compared to the experimental data. Using the method we have searched for discrete flavour symmetry groups with order less than 1035. With the applied representation selection criteria, none of the studied groups provides a symmetry which predicts lepton masses and mixing simultaneously compatible with the experimental data. Partial and approximate agreement with the data has been explored and appeared to exist for a small number of groups.
The method can accommodate different selection criteria by altering the presented algorithm. For example the restriction to U(3) subgroups can be relaxed, and unfaithful, reducible and/or lower-dimensional representations can be admitted. Also, the requirement that the tensor product decomposition is unique could be omitted, allowing multi-dimensional solutions. Presumably it will thus be feasible for a given group, to either detect a symmetry compatible with the experimental data or to rule out any such symmetry even for the most relaxed hypothesis.
There are several possible improvements and extensions to the described method. The present results are based on numerical sampling in four-dimensional VEV-space, where the choice of search interval and the sampling density is limited for practical reasons. The analysis would be highly enhanced when analytic expressions for the eigenvalues of mass matrices are used to find bounds on solutions of physical quantities.

Data Availability Statement:
The data presented in this study are available upon reasonable request from the corresponding authors.

Conflicts of Interest:
The authors declare no conflict of interest.