Dirac neutrinos in the 2HDM with restrictive Abelian symmetries

Recently, there has been a growing interest in extensions of the Standard Model in which naturally small Dirac neutrino masses arise due to existence of a symmetry which protects neutrino's Diracness. Motivated by this, we consider an extension of the Standard Model with a second Higgs doublet (2HDM) and three right-handed neutrinos where lepton number is conserved and, thus, neutrinos are Dirac particles. In this framework, we identify the most restrictive texture-zero combinations for the Dirac-neutrino and charged-lepton mass matrices that lead to masses and mixings compatible with current experimental data. We then investigate, in a systematic way, which of these combinations can be realized by Abelian continuous U(1) or discrete $\mathbb{Z}_N$ symmetries. We conclude that, from the 28 initially possible sets of maximally-restricted lepton mass matrices, only 5 have a symmetry realization in the 2HDM. For these cases, one-to-one relations among the Yukawa couplings and the neutrino mass and mixing parameters are established, and the fermion interactions with the neutral and charged scalars of the 2HDM are also determined. Consequences for lepton universality in $\tau$ decays and rare lepton-flavor-violating processes are also discussed.


I. INTRODUCTION
The discovery of neutrino masses and mixing through the observation of neutrino oscillations provided so far the only laboratory evidence for physics beyond the Standard Model (SM). In the last decades, several experiments using neutrinos from various sources (the Sun, the atmosphere, reactors, and accelerators, among others) have been measuring most of the parameters involved in their flavor oscillations with very good precision. In spite of this remarkable achievement, there remain several questions to be answered about neutrinos. Perhaps the most fundamental one is related with their nature, namely, whether neutrinos are Dirac or Majorana particles, which is of utmost importance when constructing SM extensions with massive neutrinos. Unfortunately, this question cannot be addressed by neutrino oscillation experiments and the experimental data presently available is compatible with the existence of either Dirac or Majorana massive neutrinos. Establishing the nature of neutrinos has proved very challenging and, among the several proposals to solve this riddle [1][2][3][4][5][6][7], the most promising one seems to rely on neutrinoless double beta decays [8]. Meantime, in the absence of a solid evidence in favor (or against) the existence of Dirac and/or Majorana neutrinos, both scenarios should be equally considered.
Dirac neutrino masses require adding singlet righthanded (RH) neutrino fields ν R to the SM field content. At the same time, Majorana mass terms m R ν c R ν R must be forbidden to ensure lepton number conservation. The main objection that may be raised against this simple Dirac-neutrino scenario is that extremely small Yukawa couplings ∼ O(10 −11 ) are required to generate sub-eV neutrino masses via the usual Higgs mechanism. This is in contrast with the Majorana neutrino case, in which neutrino masses are naturally suppressed by the existence of a large scale Λ commonly related with the mass of some new (heavy) particles. From the effective theory viewpoint, integrating out these heavy mediators gives rise to the dimension-five Weinberg operator ν L ν L φ 0 φ 0 /Λ [9], from which naturally small Majorana masses are generated after electroweak symmetry breaking (EWSB). The ultraviolet completions of the SM that realize this operator at tree level lead to the well-known type I [10][11][12][13][14][15], II [16][17][18][19][20] and III [21][22][23] seesaw mechanisms.
If neutrinos are indeed Dirac particles, the smallness of their masses should be natural [24] in the context of new symmetries beyond those of the SM. 1 To avoid dealing with very tiny couplings, several proposals have been put forward to explain Dirac neutrino mass suppression from an effective theory viewpoint [27][28][29]. Obviously, any of such solutions must contemplate a symmetry which forbids the dimension-four term L ν R Φ. In this way, small Yukawa couplings may originate from higher-dimensional operators which can be realized as in the Majorana neutrino case through, for instance, a Dirac seesaw mechanism [30,31]. At lowest order in the effective theory, one can define the dimension-five operator where L and ν R are a left-handed (LH) lepton doublet and RH neutrino singlet, respectively, and Φ is the SM Higgs doublet. Notice that, in this case, the generic couplings Y do not need to be very small. The scalar singlet S, after acquiring a vacuum expectation value (VEV) S = v S , leads to small effective couplings Y ν L ν R Φ with Y ν ≡ Yv S /Λ 1 as long as Λ v S .
In the above scenario, SM-like Dirac neutrino masses through a dimension-four term can be forbidden with a simple Z 2 symmetry under which ν R and S are odd, while the SM fields are even. 2 To further prevent the appearance of the ν R ν R Majorana mass term, at least a Z 3 symmetry should be considered [33,34]. For instance, the charge assignment ν R ∼ ω, with ω 3 = 1 and all remaining fields transforming trivially, would guarantee the absence of Y ν couplings and m R masses at the renormalizable level. A more economical scenario relies on a U(1) B−L symmetry with unconventional charges chosen in such a way that the gauge-invariant operators L ν R Φ and ν R ν R are forbidden [35].
The identification of new directions towards naturally small Dirac neutrino mass generation, has triggered a growing interest on this subject. Recently, several models have been proposed in the context of tree-level/radiative Dirac neutrino mass generation, and their classification with respect to the dimensionality of the corresponding generating operators has been put forward [24,27,32,[36][37][38]. Links to the dark matter problem have also been established within frameworks where the symmetry protecting neutrino Diracness play the additional role of stabilizing the dark matter particle [37][38][39][40][41][42][43][44][45]. Moreover, mechanisms for the generation of the matter-antimatter asymmetry can also be envisaged within scenarios with Dirac neutrinos [46,47] and, ultimately, both these problems can be related [48][49][50].
Independently of the mechanism for neutrino mass generation, one is always confronted with the problem of explaining the observed neutrino mass and mixing pattern. While the general frameworks described above provide an explanation for the smallness of (Majorana or Dirac) neutrino masses, in general they do not address the flavour problem per se. Thus, one is compelled to consider more sophisticated realizations of certain neutrino mass models in which flavour symmetries are considered. One of the approaches is to explore the existence of vanishing elements (texture zeros) in the Yukawa and mass matrices which reflect the violation of a symmetry by a certain interaction [51][52][53][54][55][56][57][58][59][60][61][62]. The simplest of these symmetries are those based on continuous U(1) and discrete Z N transformations.
As it will be shown later, in the SM extended with RH neutrinos, U (1) or Z N -motivated texture zeros turn out to be incompatible with data since, in general, they lead to massless neutrinos or charged-leptons and/or vanishing lepton mixing angles (which are already excluded by data). This is a direct consequence of the fact that all fermions in the SM couple to the same Higgs field. However, this is not the case in the 2HDM [63].
Abelian symmetries have been used in the context of the 2HDM with controlled flavor-changing neutral cur-rents (FCNC). Namely, in the so-called Branco-Grimus-Lavoura (BGL) scenarios [64][65][66], flavor-changing couplings (FCC) depend only on the Cabibbo-Kobayashi-Maskawa (CKM) quark mixing matrix or on the lepton mixing matrix. In these cases, the symmetry realizes a minimal flavour violation scenario without addressing the question of reconstructing the parameters of the Lagrangian (e.g. Yukawa couplings) in terms of the observable fermion mass and mixing parameters. In the most economical generalized BGL realization for quarks, in which all Yukawa couplings are real and CP violation arises from the relative phase between the VEVs of the neutral components of the two Higgs fields [67], there are nineteen parameters to be confronted with ten physical quark mass and mixing observables (six masses, three mixing angles and one Dirac CP-violating phase). However, if CP is broken explicitly in the Yukawa interactions, that number increases.
In this work, we adopt a different approach by considering restrictive Abelian symmetries in the context of the 2HDM, extended with RH neutrino fields such that neutrinos are Dirac particles. We call these symmetries restrictive since the number of relevant flavor (Yukawa coupling) parameters in the lepton sector is the same as the number of observables, i.e. nine (ten) in the case of two (three) massive Dirac neutrinos.
The paper is organized as follows. After a brief introduction on the 2HDM and general considerations regarding Dirac neutrinos presented in Section II, some general aspects of Abelian symmetries in the 2HDM with Dirac neutrinos are described in Section III. We then identify in Section IV the maximally-restrictive textures for the lepton mass matrices compatible with present neutrino oscillation data. Section V is devoted to investigate whether those textures are realizable by Abelian flavor symmetries in the 2HDM. Here, we apply the procedure developed in Ref. [68], together with the Smith normal form (SNF) and canonical methods for the identification of the U(1) flavor symmetry [69,70]. We conclude that, among the 28 combinations of maximally-restrictive texture-zero mass matrices compatible with data, only 5 can be realized by Abelian symmetries in the 2HDM. In Section VI, we discuss the reconstruction of the Yukawa parameters in terms of lepton masses and mixing, paying special attention to leptonic CP violation. The constraints coming from lepton-flavour universality in leptonic τ decays, δ are also discussed. Finally, our concluding remarks are given in Section VII.

II. DIRAC NEUTRINOS IN THE 2HDM
As in the SM, Dirac neutrino masses can be generated in the 2HDM by adding RH neutrino singlet fields ν R , which couple to the SM lepton doublets L and the two Higgses Φ a defined, as usual, by being φ + a and φ 0 a the charged and neutral components of the scalar doublets. In this framework, the lepton Yukawa interactions can be written as where e R are the charged-lepton RH singlets, and Φ a = iσ 2 Φ * a . The general 3 × 3 complex matrices Y a and Y ν a encode the charged-lepton and Dirac neutrino Yukawa interactions, respectively. In line with the discussion presented in Section I, very small Yukawa couplings Y ν a may originate from dimension-five operators: such that Y ν a ≡ Y a v S /Λ are sufficiently suppressed to generate sub-eV Dirac neutrino masses upon EWSB, i.e. when φ 0 a acquire VEVs: with v 246 GeV. Thus, from now on we will consider that there is such a mechanism responsible for the smallness of Y ν a .
The charged-lepton and Dirac-neutrino mass matrices can be diagonalized by a set of appropriate unitary matrices U ,ν L,R so that where m e,µ,τ and m 1,2,3 denote the charged-lepton and neutrino masses, respectively, all being real and positive.
To extract the LH rotation matrices, one diagonalizes the Hermitian matrices in the following way: The unitary transformations U ,ν L define the lepton mixing matrix U appearing in lepton charged-current interactions as In the case of massive Dirac neutrinos, U can be parametrized by three mixing angles θ ij and a single CPviolating Dirac phase δ, such that [71] with c ij ≡ cos θ ij and s ij ≡ sin θ ij . Thus, the lepton sector is characterized by ten parameters: three chargedlepton and three neutrino masses, three mixing angles and one phase. Global analyses of all available neutrino oscillation data constrain the parameters of the matrix U as shown in Table I [ [72][73][74], for both normal (m 1 < m 2 < m 3 ) and inverted (m 3 < m 1 < m 2 ) ordering of neutrino masses. In the case of m 1 = 0 (m 3 = 0), m 2 and m 3 (m 1 and m 2 ) are fully determined by the mass-squared differences ∆m 2 21 = m 2 2 − m 2 1 and ∆m 2 31 = m 2 3 − m 2 1 as NO : m 2 = ∆m 2 21 , m 3 = ∆m 2 31 , so that the number of measurable quantities is reduced to nine. The cosmological constraint on the sum of the neutrino masses, i m i < 0.12 eV at 95% CL [75], should also be taken into consideration. If one neutrino is massless, then i m i 0.059 eV for NO and i m i 0.099 eV for IO, thus being the cosmological limit automatically obeyed.

III. ABELIAN SYMMETRIES IN THE 2HDM WITH DIRAC NEUTRINOS
As shown in Appendix A, the most constrained mass matrices M and M ν (compatible with experiment) achievable by imposing Abelian symmetries in the SM with RH neutrinos contain 23 real parameters to be compared with 9 measurable quantities (lepton masses and mixings). This is not the case in models with extra scalar fields like the 2HDM, where Abelian symmetries may lead to interesting constraints in the flavor sector [51,76,77]. Thus, in this section we will go through some general aspects on this subject as a setup for the forthcoming analysis presented in this work. Denoting Φ ≡ (Φ 1 Φ 2 ) T , and requiring that the full Lagrangian is invariant under the field transformations where S Φ ∈ U(2) and {S , S e , S ν } ∈ U(3) are unitary matrices, yields the following constraints on the Yukawa couplings where a sum over b = 1, 2 is implicitly assumed. By performing basis transformations identical to those in Eq. (13), with the appropriate choices of unitary matrices V ∈ U(2) and {V , V e , V ν } ∈ U(3), one can bring the matrices S into the form [77]: where θ i , α i , β i and γ i are continuous phases. Under the general transformation (15), the invariance condition (14) reads where i, j = 1, 2, 3 are flavor indices and x = , ν. The phase matrices Θ x a , which encode the transformation properties of each Yukawa interaction, are defined as These phases can be written in terms of charges (α , β , γ , θ ) and a parameter ϕ ∈ [0, 2π[ such that The particular case of ϕ = 2π/N , with N = 2, 3, ..., corresponds to a discrete Z N symmetry. With these redefinitions, Eq. (14) can be interpreted in terms of charge relations. Namely, invariance of (Y x a ) ij under the U(1) symmetry implies (Θ x a ) ij = 0 (mod 2π). It is also straightforward to conclude that, as a consequence of having θ 1 − θ 2 = 0 (mod 2π), 3 a non-zero entry in Y x 1 will automatically imply a zero entry in Y x 2 , and vice versa. Moreover, one can see from Eq. (14) that overall rephasings of the type with e i(α−β−θ) = e i(α−γ+θ) = 1, do not alter the invariance condition on the Lagrangian. Thus, without loss of generality, we can set one of the Higgs and fermion transformation phases to zero. Here, we choose θ 1 = α 1 = 0. The aforementioned U(1) charges determine the presence (or absence) of zero entries in the Yukawa and mass matrices Y x a and M x , defined in Eqs. (3) and (6), respectively. In particular, Vanishing elements in a mass matrix or Yukawa interaction matrix are usually dubbed as "texture zeros". In this work, whenever a general matrix structure contains texture zeros, we will refer to it as a "texture". Before concluding this section, it is worth comparing the lepton and quark sectors in the framework considered in this work, i.e. the 2HDM with RH neutrinos. The Yukawa interactions for quarks are the same as in Eq. (3) after replacing L by Q L (LH quark doublets), e R by d R (RH down-quark singlets) and ν R by u R (RH up-quark singlets). The charged-lepton and Dirac neutrino Yukawa matrices must be also replaced by the down and up-quark ones, Y d a and Y u a , respectively. Textures for quarks can, in principle, be implemented in the same way as for leptons by imposing Abelian symmetries as the ones discussed above. However, as explained in Appendix A, this cannot be done in the SM. As in the lepton 3 Otherwise one recovers, in terms of a phase transformation like in Eq. (17), the SM context, since the same texture would be enforced on the Yukawa matrices Y x 1 and Y x 2 and, hence, in the resulting mass matrix Mx. sector, quark textures can be implemented by Abelian symmetries in the 2HDM, being the main difference in the fact that all six quarks must be massive. This is in contrast with leptons, for which the existence of a massless neutrino is allowed by current experimental data.

IV. MAXIMALLY-RESTRICTIVE TEXTURES FOR LEPTONS
In Ref. [78], all possible textures for M and M ν were identified and grouped into equivalence classes, considering both Majorana and Dirac massive neutrinos.
For charged leptons, two textures M and M are equivalent if they can be transformed onto each other by performing permutations of the L and e R fields, i.e. if where P ,e can be any two matrices of the 3-dimensional representation of the S 3 permutation group, namely, Said otherwise, two M textures are equivalent when they are equal up to permutations of rows and columns.
In the case of Dirac neutrinos, two textures M ν and M ν are considered equivalent if where P ν is also a permutation matrix. Thus, two neutrino mass matrix textures are equivalent if they can be transformed onto each other by column permutations. We shall combine the above M and M ν classes into all possible (M , M ν ) pairs, keeping only one representative texture of each M and M ν equivalence class. Pairs leading to the same leptonic mixing matrix U are equivalent and, thus, redundant. This is the case when the mass matrices can be related by for any two texture pairs (M , M ν ) and (M , M ν ). Notice that, in order to leave U invariant, P must be the same in both transformations. Therefore, two texture pairs are equivalent if they can be obtained from each other through column and row permutations, being the row permutation identical for both mass matrices in the pair. This is why in Eq. (23) only column permutations are considered, avoiding the possibility of excluding relevant cases. The outlined procedure aims at eliminating redundant cases that lead to the same physics, i.e., that reproduce the same mass and mixing parameters. Since in this work we are interested in the most predictive M and M ν , it is crucial to introduce the concept of maximally-restrictive textures [78,79]: A texture pair (M , M ν ) is said to be maximally restrictive if the predicted values for the lepton masses, mixing angles and CP phase are compatible with the experimental data, and the addition of one more texture zero in either M or M ν makes the pair incompatible with data.
Essentially, these are the pairs with least parameters, which are viable when confronted with observations. In order to identify them, we shall perform an analysis similar to the one of Refs. [78,79], considering the updated neutrino oscillation parameters and including the current ranges for the Dirac phase δ (see Table I). We will require compatibility at 3σ confidence level (CL).

A. Compatibility with data
In order to identify the maximally-restrictive texture pairs (M , M ν ) we performed a standard χ 2 -analysis with the function where x denotes the matrix elements of M and M ν , P i (x) is the model prediction for the observable O i , O i is the corresponding best-fit value, and σ i denotes its 1σ error. In our search for viable pairs (M , M ν ), we require the charged-lepton masses to be at their central values [71], so that the χ 2 -function is minimized only with respect to the six neutrino observables O i (the two neutrino masssquared differences ∆m 2 21,31 , the three mixing angles θ ij and the Dirac phase δ) following the numerical method presented in Refs. [56,79]. If the deviation of each neutrino observable from its experimental value is at most 3σ at the χ 2 minimum for a given (M , M ν ) pair, the corresponding lepton textures are said to be compatible with data. In such cases, we test compatibility at the 1σ as well.
Our results show that the maximally-restrictive pairs (M , M ν ) compatible with data are those presented in Table II, where the labeling follows the notation of Ref. [78]. 4 A representative texture of each equivalence   Table I). The pair (6 1 , 4 ν 17 ) was found to be consistent with experimental data only at 3σ and for IO.  Table II. class is presented in Tables III and IV for M and M ν , respectively. All pairs in Table II were found to be consistent with neutrino oscillation data at 1σ, for both NO and IO, except for the pair (6 1 , 4 ν 17 ) which is consistent with data only at 3σ CL and for a NO mass spectrum. Also notice that, with the exception of texture 4 ν 17 , any representative of M ν given in Table IV features a massless neutrino, since it contains a full column of zeros.
We emphasize that these maximally-restrictive texture pairs cannot be implemented in the SM by imposing Abelian symmetries (see Appendix A). Hence, in the next section we will address the question of whether (or which of) the texture pairs in Table II can be implemented in the context of the 2HDM with Abelian flavor symmetries.  Table II.

V. ABELIAN SYMMETRY REALIZATION OF COMPATIBLE TEXTURES
Having identified the maximally-restrictive texture pairs that are compatible with data, next we aim at ascertaining the pairs in Table II that can be realized in the 2HDM by imposing continuous or discrete Abelian symmetries. At the same time, we wish to identify the corresponding transformation properties of the various fields according to Eq. (13). Keeping this in mind, two methods shall be employed, namely the canonical and SNF methods, which we briefly review below.
In the canonical method, the phases in Eq. (18) are considered as variables and, for each texture pair of Table II, a system of equations corresponding to the conditions in Eq. (20) is defined. Subsequently, one checks if the system has a non-trivial solution for the phases, considering all possible decompositions of M and M ν into Yukawa matrices Y 1,2 and Y ν 1,2 . If no solution is found, the texture (or texture pair) is not realizable by means of an Abelian symmetry (see, for example, Refs. [68,77] for more details).
The SNF method is used to identify rephasing symmetries of the Lagrangian for specific decompositions of (M , M ν ) texture pairs into Yukawa matrices. With n f complex fields, and in the absence of phase-sensitive terms, the Lagrangian is invariant under a [U(1)] n f symmetry. In the presence of terms transforming nontrivially under rephasing, the symmetry group is altered and the SNF method can be used to identify it. With this purpose, one starts by defining a k × n f integer matrix D where each row d l corresponds to one of the k interactions allowed in the Lagrangian, with where the number of d f components is equal to the number of fields of type f , i.e. two for d Φ and three for the remaining fields. Specifically, d fi = n (d fi = −n) when the (conjugate of) field f i appears n times in the interaction. For instance, d l = (0, 1, −1, 0, 0, 0, 0, 1, 0, 0, 0) corresponds to the coupling e L τ R Φ 2 .
The matrix D can then be brought to its unique SNF [80] such that D = R D SNF C with where d i is a positive integer, divisor of d i+1 , and r is the rank of D. The matrices R and C encode the addition, sign flip and permutation operations on the rows and columns, respectively. The minimal rephasing symmetry of the Lagrangian can be identified as where d i = 1 and d i = 0 correspond to a continuous U(1) symmetry and the trivial group, respectively. For the case of discrete groups, the symmetry charges can be determined from the columns of C −1 , although no information regarding the charges associated to continuous symmetries can be obtained. We should emphasize that the SNF method identifies symmetries under which all interactions included in D are invariant, without guaranteeing that all remaining terms are absent from the Lagrangian. A more detailed discussion on this point can be found in Ref. [80]. We first apply the canonical method to find which of the (M , M ν ) pairs given in Table II cannot be implemented with Abelian flavor symmetries in the 2HDM. If both M and M ν textures can be implemented but not simultaneously, then the corresponding pair is not realizable by the symmetry. From this analysis, which is detailed in Appendix B, we conclude that 23 of the 28 maximally-restrictive pairs (M , M ν ) appearing in Table II cannot be realized through Abelian symmetries in the present context. Moreover, all redundant pairs (not presented in the table) which are equivalent to these 23 pairs cannot be implemented either. 5 In conclusion, only the 5 pairs (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ) have no inconsistencies and can be implemented by an Abelian symmetry. In order to identify it, we now study all possible decompositions of these pairs into Yukawa matrices.

A. Decomposition into Yukawa textures
In Ref. [68], a thorough analysis on the realization of Yukawa textures through Abelian symmetries in multi-Higgs-doublet models was performed. In particular, it was established that any Yukawa texture Y x realizable by a transformation of the type (15) can be expressed as where A k is one of the 3 × 3 texture matrices defined in Eq. (31) of Ref. [68], and P and P are permutation matrices. The A k textures are classified according to the degeneracy of the LH and RH symmetry phases α i and β j in Eq. (15). Those A k that are realized through transformations with the same phase degeneracy are grouped together in classes (i,j), where i and j indicates the number of non-degenerate α i and β j , respectively. Within each class (i,j), all P A k P which can be simultaneously realized without common non-zero elements are grouped into "chains", being C (i,j) n the nth chain of the (i,j) class. In the 2HDM framework, a given decomposition of a mass matrix M x into Yukawa textures is realizable through Abelian symmetries if there is at least one chain containing the corresponding Y x 1 and Y x 2 . Since none of the textures belonging to the pairs (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ) has identical columns, and neither 4 3 nor 5 1 has identical rows, the implementation of any of these textures require non-degenerate transformations, i.e three distinct α i , three distinct β i and three distinct γ i phases. Thus, the search for decompositions is limited to the class (3,3), generated by with the corresponding chains C given in Ref. [68]. By inspecting these chains, we determine the viable decompositions for each matrix of the (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ) pairs (see Table V). The associated chains and their corresponding building matrices are shown in Tables VI and VII, for the charged-lepton and neutrino textures, respectively.
Based on Table V, we generate all possible pairs of (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ) decompositions. Notice that one must consider both Yukawa matrix orderings, since swapping Y ν 1 and Y ν 2 , while maintaining the charged-lepton ordering of Y 1 and Y 2 , will affect the corresponding implementation. In what follows, we choose to fix the charged-lepton Yukawa ordering and consider two possible orderings for the neutrino Yukawa textures.
At this point, it remains to determine which of the decomposition pairs are viable, i.e those in which all four Yukawa matrices can be simultaneously realized by the same symmetry. Obviously, at least one pair of decompositions must be realizable for at least one Y ν 1,2 ordering since, otherwise, the corresponding mass matrix pair would have been already excluded by the canonical method. Given that we have reduced the number of possible combinations to only 16, including the two possible orderings (12 and 4 for the 4 3 and 5 1 textures, respectively), the SNF method becomes now more practical to identify the minimal Abelian symmetry realization.

B. Minimal Abelian symmetries
We apply the SNF method to all the (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ) decomposition pairs in Table V, considering both orderings of Y ν 1,2 . We find that D SNF can only take one of the following two forms: where, according to the discussion around Eq. (27), the two cases correspond to the Abelian symmetry groups G This result can be interpreted as follows. Irrespective of the type of Yukawa textures, the Lagrangian is invariant under the global continuous symmetries U(1) Y and U(1) L , where Y is the SM hypercharge and L is the lepton number. These global symmetries simply appear because the 2HDM is invariant under the global version of the U(1) Y gauge symmetry, while U(1) L ensures lepton number conservation, which is required when Dirac neutrino Yukawa interactions are considered. Notice that neither of these symmetries imposes texture zeros in the mass matrices. As such, the flavor symmetry identified by the SNF method is, in fact, For all Yukawa matrices in Table V, couplings with the field ν R1 are forbidden, leading to the existence of a massless Weyl neutrino. The absence of phase-sensitive terms involving ν R1 implies a continuous U(1) symmetry associated to that field, which is preserved from the original rephasing symmetry [U(1)] n F of the Lagrangian without phase-sensitive terms. Therefore, the actual SNF symmetry which realizes the non-zero elements of the decomposition pairs is Texture decomposition Building matrices Associated chains where U(1) ν Ri follows from the invariance under ν Ri → e iα ν Ri , with all remaining fields transforming trivially. Notice that we take ν Ri with i = 1, 2, 3 since row permutations of M ν generate an equivalent pair (M , M ν ), as seen from Eq. (24). After applying the SNF we have concluded that the non-zero elements of all Yukawa matrices in each decomposition can be implemented (i.e. are compatible) by a U(1) or a Z 2 symmetry. Moreover, a U(1) ν Ri is present due to the absence of ν Ri interactions. It now remains to check whether the texture zeros can be implemented by the same symmetries.
Using the canonical method, we determine the decompositions and orderings of the texture pairs (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ) that are exactly realized by the symmetries (34) (see Table VIII). We find that all 5 mass matrix pairs have one decomposition which can be implemented imposing a G F = U(1) × U(1) ν Ri symmetry, being U(1) alone able to reproduce the entire corresponding texture structure. This means that U(1) ν Ri is, in fact, an accidental symmetry arising from the particular form of the M ν texture. The viable ordering in each case corresponds to that indicated in Table V. The three decompositions which cannot be reproduced (and are omitted in Table VIII) exhibit, for one ordering each, a rephasing symmetry G F = U(1) × U(1) ν Ri that is unable to enforce one of their respective zero textures. On the other hand, the G F = Z 2 × U(1) ν Ri cases correspond to the alternative Yukawa orderings of all 8 decompositions, with the Z 2 rephasing symmetry simply enforcing the corresponding ordering, while it fails in imposing any texture zero on the mass matrices.
From now on, we will only refer to the U(1) symmetry, corresponding to g = 0 in Eq. (34). The U(1) implementation of the realizable (M , M ν ) pairs is summarized in Table IX, were the continuous phases and discrete charges are presented. Out of the 28 maximallyrestrictive (M , M ν ) texture pairs compatible with data given in Table II, only 5 can be realized through Abelian symmetries in a 2HDM, each admitting a single decomposition into Yukawa matrices, realizable by a U(1) flavor symmetry. As it turns out, in all cases, the minimal set of discrete charges corresponds to a Z 5 symmetry.

VI. PHENOMENOLOGY
In the previous section we have selected among all maximally restricted texture-zero pairs (M , M ν ) those which are simultaneously compatible with neutrino data and realizable in the context of the 2HDM with a U(1) flavor symmetry. In the following, we will focus on their phenomenology in what concerns leptonic CP violation, rare lepton decays and lepton universality. More specifically, we will establish the relation between the only complex phase which remains in the Yukawa sector of the theory with that appearing in the lepton mixing matrix U of Eq. (11), i.e. the Dirac phase δ. We also apply the constraints coming from universality tests in τ decays and from rare (two-and three-body) lepton decay searches.

A. Leptonic CP violation
Although, in general, all elements of the Yukawa matrices Y ,ν a in Eq. (3) are complex, some phases have no physical significance and, thus, can be removed by rephasing the fermion fields as ψ j → e iϕj ψ j . It is straightforward to show that, for the maximallyrestricted pairs of matrices considered in this work, all elements of M and M ν (or Y ,ν 1,2 ) can be made real and positive, except one. This is due to the fact that a single unremovable phase α, which will be necessarily correlated with the Dirac CP-violating phase δ in Eq. (11), remains after exhausting the rephasing freedom. Our convention for the position of α is given in Table VIII, where all parameters a i , b i , x i and y i are real and positive.
Given that for all (M , M ν ) pairs there are nine real parameters in total, four of them remain undefined after ensuring a mass spectrum that reproduces the observed charged-lepton masses m e,µ,τ and neutrino mass-squared differences ∆m 2 21,31 (since one neutrino is massless, the two Dirac neutrino masses can be written in terms of ∆m 2 21,31 ). In both 4 3,I and 6 ν i,I two free parameters remain, while in 5 1,II and 5 ν 8,I we are left with one and three free parameters, respectively. For the texture 5 1,II (4 3,I ) we choose to write a 1 and b 1,2 (a 1,2 and b 1 ) in terms of a 2 (a 3 and b 2 ) and m e,µ,τ . As for 6 ν i,I (5 ν 8,I ), we express x 2 and y (x 1 and y 1 ) in terms of x 1 (x 2 , y 2 ) and ∆m 2 21,31 , as shown in Table XI of Appendix C. For the texture pairs in Table VIII, we perform a scan of the free parameters in their validity ranges and, for each input set, H and H ν are defined as in Eq. (8). After diagonalizing these two matrices, U gets determined by Eq. (10), and θ ij and  Table VIII. NH1,2 and IH1,2 solutions are classified according to which neutrino is considered to be the lightest (see Appendix C for more details). In all points, the mixing angles θij lie within the 3σ ranges given in Table I. δ appearing in Eq. (11) can be extracted. Demanding agreement with the 3σ ranges given in Table I, we plot δ  Table I.  Table V. We have took advantage of fieldrephasing freedom to eliminate all unphysical phases, placing the only irremovable phase α in Mν .
as a function of α. Our results are shown in Figs. 1 and 2 for the (4 3 , 6 ν i ) and (5 1 , 5 ν 8 ) pairs, respectively. In all cases, both the NH and IH mass spectra are considered.
By looking at Fig. 1 we see that the results are similar for different pairs of (M , M ν ) textures. This could suggest that those pairs are equivalent in the sense that they can be transformed into each other by permutation transformations (see discussion at the beginning of Section IV). However, we advocated that only nonequivalent pairs of textures were kept. The similarity between some of the results have to do with the fact that up to a very small parameter, which does not have much impact in the results for the mixing angles and phases, those pairs are indeed equivalent. Let us illustrate this with one example, namely (4 3 , 6 ν 3 ) and (4 3 , 6 ν 7 ). By looking at Table VIII we notice that 6 ν 3 = P 13 6 ν 7 P 23 , 4 3 = a3→0 P 13 4 3 P 23 .
The last relation simply indicates that 4 3 = P 13 4 3 P 23 in the limit a 3 → 0 and, thus, the two pairs would be equivalent in this case. 6 Obviously, a 3 = 0 would lead to a massless charged lepton, which is not acceptable. However, it turns out that a 3 is very small and, indeed, 4 3 P 13 4 3 P 23 , leading to a similar U for both pairs. A similar reasoning can be applied to the remaining pairs for which the results look the same. We emphasize that the existence in the mass matrices of such small parameters is a consequence of the hierarchical nature of the charged-lepton masses. We now turn our attention to the results for the (5 1 , 5 ν 8 ) pair shown in Fig. 2. As will become clear in the next section, this case is the most interesting one and deserves a more detailed attention. The first to notice is that in the matrix 5 1 one of the charged-lepton states is decoupled. In particular, the Hermitian matrices H = M M † and H = M † M and their diagonalizing unitary transformations U L and U R , are given by where c L,R ≡ cos θ L,R and s L,R ≡ sin θ L,R with As explained in Appendix B, m 2, 3 correspond to the masses of the two non-decoupled charged-lepton states. Thus, depending on which charged lepton 1 is identified as decoupled, three different cases must be considered: This explains the notation used in the plots of Fig. 1. Taking into account that the unitary matrix U L must be such that Eq. (9) is verified with the correct chargedlepton mass ordering, we have When M ν is of the type 5 ν 8 , the unitary matrix which diagonalizes H ν as in Eq. (9) is presented in Eqs. (C1) and (C2) of Appendix C for NH and IH, respectively. After some algebra, we conclude that the charged-lepton rotation set by θ L given in Eq. (36) is crucial since compatibility with the measured neutrino mixing angles θ ij would not be possible from U ν L alone. 7 In particular, if one considers U = U ν L , then the following conditions must be satisfied when M is of the 5 e 1 type: where r ≡ ∆m 2 21 /∆m 2 31 0.03, according to the data given in Table I. The above relations imply tan 2 θ 23 = 4s 2 13 (1 ± rs 2 12 ) 2 r 2 sin 2 (2θ 12 ) 4s 2 13 r 2 sin 2 (2θ 12 ) 110 , (46) where the −(+) sign in ± is for the NH (IH) case, and in the final estimate the values of Table I have been considered. As for the 5 µ,τ 1 , relations (44) and (45) are replaced by leading to θ 23 θ, π/2 − θ with which implies θ 23 0, π/2. From these results we conclude that the contribution to the mixing coming from the charged-lepton sector is crucial to get compatibility with data. Moreover, it can be shown that the Jarlskog invariant J CP , which signals Dirac-type CP violation, obeys confirming the fact that, for CP violation to occur in the lepton sector, θ L = nπ/2 and α = nπ (n is integer) must hold. In the following, we will obtain relations among the parameters in M and M ν (a 2 , x 2 , y 2 and α) and the three mixing angles θ ij and the CP phase δ. Again, we focus on the (5 e 1 , 5 ν 8 ) pair. From Eqs. (10), (37), (41), (C1) and (C2), the lepton mixing matrix U is computed and the mixing angles and the phase δ are extracted. Notice that, for the case (5 e 1 , 5 ν 8 ), one has U 1j = (U ν L ) 2j . Therefore, given the parametrization (11), x 2 and y 2 in M ν depend only on It now remains to express θ L (or b 2 ) appearing in Eq. (37) and the phase α in terms of the measurable neutrino parameters. Including the charged-lepton corrections to the mixing we have where t L ≡ tan θ L , c α ≡ cos α and s α ≡ sin α. In the limit θ L , α → 0, we recover the result (46) where the − (+) sign corresponds to the NH (IH) case.
The above relation provides a very good approximation for the behavior of the charged-lepton mixing angle θ L in terms of θ ij , r and α. From Eqs. (38), (54) and the defining conditions for a 2 1 and b 2 1,2 given in Table XI, one can determine the parameters a 1,2 and b 1,2 in terms of m e , m µ , m τ and θ 23 , namely, In order to relate δ with α, we notice that arg(U 23 ) arctan 2s 13 s α rt L sin(2θ 12 arg(U 33 ) − arctan 2s 13 t L s α r sin(2θ 12 ) + 2s 13 for both NH and IH. These relations imply from which, after performing some rephasing transformations to bring U to the form given in Eq. (11), we obtain confirming the results plotted in Fig. 2. Following the same procedure for the 5 µ 1 (5 τ 1 ) we obtain δ α (δ π + α), which also agrees with the numerical output shown in Fig. 2. In conclusion, all parameters in the mass matrices M and M ν can be determined in terms of the chargedlepton and neutrino masses and mixing angles through Eqs. (51)- (58) and (62).

B. Lepton universality and rare LFV decays
In the 2HDM, Yukawa interactions may induce flavorchanging neutral currents (FCNC) at the tree and loop levels. Therefore, the viable maximally restrictive textures previously obtained (cf. Table VIII) should be confronted with the current experimental bounds on such processes. In particular, the constraints on universality in purely leptonic decays, lepton-flavor-violating decays − α → − β + γ − δ , and α → β γ should be considered. 8 With this purpose, we first briefly review the interactions among leptons and the neutral and charged scalars in the 2HDM. Since we are considering scenarios with a U(1) symmetry under which one of the Higgs doublets is charged, there is no CP violation in the scalar potential and, thus, no mixing between CP-even and CP-odd scalars. It is convenient to rotate (Φ 1 , Φ 2 ) to the Higgs basis (H 1 , H 2 ) through such that one can write with H 2 = 0. The neutral (charged) Goldstone boson is denoted G 0 (G + ), while I is the U(1) Goldstone boson, which is massless in the exact U(1) symmetric limit. In order to avoid this massless particle, a soft U(1) symmetry breaking term of the type m 2 12 Φ † 1 Φ 2 can be included in the scalar potential originating a mass m 2 I ∝ m 2 12 for the decoupled CP-odd scalar.
Throughout this work we will also assume that there is no mixing between R and H 0 and identify H 0 with the SM Higgs boson discovered by the ATLAS and CMS collaborations at the LHC. Thus, we consider the limit where the physical mass state h is identified with H 0 , which implies m H 0 ≡ m h 125 GeV. Under these premises, the CP-even scalars are h and R, while I is CP-odd. As already said, in our framework there is no CP-violation in the scalar sector and, thus, I is decoupled from h and R, which allows us to take m R , m I and the mass of the charged scalar m H + as independent parameters.
The relevant scalar-fermion interactions can be read off from the Lagrangian (3) which, after appropriate transformations, takes the form with In Eq. (65) all the fermion fields are mass eigenstates (see Section II for details and definitions regarding the unitary transformations and Yukawa matrix conventions). Lepton universality tests aim at probing the SM prediction that all leptons couple with the same strength to the charged weak current. A relevant quantity to test universality in purely leptonic τ decays is where x αβ ≡ m 2 α /m 2 β . In the presence of scalar and vector interactions, the α → β νν branching ratio (BR) is [81] Br( α → β νν) = 1 + 1 4 g S RR,αβ where are the phase space functions and Finally, the flavor-dependent coefficients g iαjβ are modelspecific, being in our framework given by with N e defined as in (66). Current experimental constraints yield [82] |g µ /g e | − 1 = 0.0019 ± 0.0014 (73) and g S RR,µe < 0.035, g S RR,τ e < 0.70, g S RR,τ µ < 0.72, (74) at 95% CL [71].
Currently, the experimental upper limits on the branching ratios of the 3-body LFV decays are [71] Br(τ − → e − e + e − ) < 2.7 × 10 −8 , at 90% CL. Finally, neglecting contributions proportional to the neutrino masses and sub-leading terms in m 2 /m 2 R,I , the decay width of the radiative lepton-flavor violating process α → β γ is given, up to one-loop level, by [81] Br( α → β γ) where α e = e 2 /(4π) and the amplitudes A L,R are in the present framework given by where a sum over i = e, µ, τ is implicitly assumed and chirally suppressed terms proportional to m β /m α were neglected. Current experimental upper bounds at 90% CL are [71] Br(µ → eγ) < 4.2 × 10 −13 , We now aim at studying the compatibility of the texture pairs given in Table VIII with the constraints coming from lepton universality and rare LFV decays discussed above. Simultaneously to the analysis performed in the previous section, we randomly vary tan β in the range 0.01 to 100 (these values ensure that the Yukawa couplings are always 1), the charged-Higgs scalar masses m H ± 80 GeV [83], and the neutral scalar masses m R,I 100 GeV [71]. We limit our search to cases where the m H ± 1 TeV and m R,I 10 TeV. For each input parameter set compatible with neutrino data, we compute |g µ /g e | − 1, g S RR,ij , and the BRs of all LFV 3-body and radiative charged-lepton decays. In all cases, we keep only those points obeying Eqs. (74), (78) and (81). As for |g µ /g e |, we demand |g µ /g e | − 1 ≥ 10 −4 , keeping in mind the result (73).
In Fig. 3 we show the results for the (5 e 1 , 5 ν 8 ) and (5 µ,τ 1 , 5 ν 8 ) texture pairs in the two upper and lower panels, respectively. In the left (middle) column we plot |g µ /g e |−1 as a function of m H ± (tan β), while in the right column the same points are shown in the (m R , m I )-plane. We conclude that the (5 µ,τ 1 , 5 ν 8 ) cases are disfavored by the |g µ /g e | − 1 constraint (73) (indicated by the horizontal gray bands in the plots). Instead, for (5 e 1 , 5 ν 8 ) the deviation from universality is in agreement with Eq. (73) for 80 GeV m H ± 200 GeV and tan β 0.03 or tan β 30, for both NH and IH. Notice that for large (small) tan β the Yukawa couplings in Y 1 (Y 2 ) are enhanced, leading to an enhancement of |g µ /g e | − 1.
Similar results are presented in Fig. 4 for the (4 3 , 6 ν k ) texture pairs given in Table VIII. We do not present the results in terms of tan β since the behavior is similar to that of the (5 1 , 5 ν 8 ) cases i.e., in general, there is a small and large tan β region. The main difference between the results in Figs. 3 and 4 is evident from the comparison of the (m R , m I ) plots. While for the (5 1 , 5 ν 8 ) texture pair all constraints are verified for non-correlated m R,I masses (see the left column plots in Fig. 3), for the texture sets TABLE X. Allowed α → β γ and α → β γ for 5e,µ,τ , indicated in each case by particle flavor indices (α, β) and (α, βγδ).
(4 3 , 6 ν k ) the mass tuning m R /m I 1 is needed to pass all the constraints, as shown in the second and forth column plots of Fig. 4. This is easy to understand. First we notice that the matrix N e defined in Eq. (66) has the following forms for the 5 1 and 4 3 textures: Ultimately, for a certain value of tan β, the non-zero entries marked with a × could be expressed in terms of the charged-lepton and neutrino masses and lepton mixing angles, as illustrated for the case of the 5 1 texture discussed in the previous section. Taking into account Eqs. (75), (76), (79) and (80), we can immediately conclude that most of the − α → − β + γ − δ and α → β γ are forbidden at the one loop level for the 5 1 textures. This is due to the coupling structure imposed by the U(1) flavor symmetry which, in the case of charged leptons, only allows mixing between two flavors. For instance, the decay µ → eγ only occurs when M ∼ 5 τ 1 , since for 5 e 1 (5 µ 1 ) the electron (muon) is decoupled and, thus, µ − e transitions are not allowed. On the other hand, τ radiative decays are forbidden in that case since the τ is decoupled. Applying the same reasoning to the 3-body decays − α → − β + γ − δ , we conclude that, at most, only two of these processes are allowed for each of the 5 1 case (see Table X). Thus, the stringent constraints coming from the µ decays, are naturally satisfied in the 5 e,µ 1 case, while for the 5 τ 1 texture the rates are suppressed by the small couplings in the µ − e sector. Summarizing, for 5 1 the constraints from LFV decays are respected without requiring any special relation among the scalar masses m R,I , as can be seen from the plots in Fig. 3.
The natural suppression of LFV decays does not however occur when M ∼ 4 3 . As can be seen from (83), in these cases the couplings N e do not exhibit any decoupling behavior and, thus, the decay rates are not naturally suppressed. In the particular case of µ → eγ, the terms enhanced by m τ /m µ are potentially large and the experimental bound on that decay is respected only In the left (middle) columns we plot |gµ/ge| − 1 as a function of m H ± (tan β). The horizontal grey bands correspond to the constraint (73). In the right column, the same points as in the corresponding |gµ/ge| − 1 plots are shown in the (mR, mI )-plane. All points obey the constraints (78) and (81), and the mixing angles θij lie within the 3σ ranges given in Table I for NH and IH. when there is a cancellation between the two terms proportional to m τ /m µ in A R , i.e. when m R m I . This is reflected in the (m R , m I ) plots of Fig. 4, where that correlation is clear. Notice that in the 5 1 case those terms were absent since (N e ) µτ (N e ) τ e = 0, which is not the case for the 4 3 textures. For illustration, we show in Fig. 5 the dependence of Br(µ → eγ) on the mass ratio m I /m R for the texture pair (4 3 , 6 ν 7 ), which confirms the fact that quasi-degenerate m R,I masses are required to respect the MEG µ → eγ bound. Notice also that, in all cases, sizable deviations from lepton universality require m H ± 200 − 300 GeV.  4. Results for the (4 3 , 6 ν k ) for k = 1, 3, 7, 9 (from top to bottom rows). In the first (third) column we plot |gµ/ge| − 1 as a function of m H ± for NH (IH). The horizontal grey bands correspond to the experimental constraint (73). The same points as in the corresponding |gµ/ge| − 1 plots are shown in the (mR, mI )-plane for NH (IH) in the second (forth) column. In all points the constraints (78) and (81) are satisfied, and the mixing angles θij lie within the 3σ ranges given in Table I for NH or IH.

VII. CONCLUSIONS
At present, no tangible experimental evidence exists in favor of Dirac or Majorana massive neutrinos. In pure theoretical grounds, both scenarios can be implemented in a natural way from the point of view of an effective SM in which neutrino masses are suppressed by a large scale. In the case of Dirac neutrinos, once RH neutrino singlets are added, their direct couplings with the LH lepton and the SM Higgs doublets must be forbidden by a symmetry. In this way, the window is open for seesawlike suppression of Dirac neutrino masses. This principle can be extended to the 2HDM to suppress the effective Dirac couplings to both Higgs doublets. The presence of the second Higgs doublet opens the possibility of implementing texture-zero structures in the charged-lepton and neutrino mass matrices compatible with data. In this work, we have addressed this problem by considering the maximally-restricted scenarios that can be implemented by imposing U(1) or Z N Abelian symmetries in the 2HDM Lagrangian. Our approach differs from usual BGL symmetries in the sense that the number of relevant flavor parameters in the Lagrangian is the same as the number of lepton masses, mixing angles and CP phases. We stress that these are the most restrictive symmetries in the sense that the number of parameters cannot be further reduced.
Our analysis shows that from the 28 initial pairs (M , M ν ) compatible with data, only 5 can be implemented in the 2HDM by imposing a U(1) flavor symmetry. Clearly, this does not preclude the possibility of implementing the remaining pairs in models with more than FIG. 5. Dependence of Br(µ → eγ) on the mass ratio mI /mR. Although the results are shown for the (4 3 , 6 ν 7 ) texture pair, for the remaining (4 3 , 6 ν k ) cases shown in Table VIII the results are similar. In all points the mixing angles θij lie within the 3σ ranges given in Table I and |gµ/ge| − 1 ≥ 10 −4 . The red points are excluded by the µ → eγ MEG bound given in (77).
two Higgs doublets. For the realizable pairs we established the relation between the only complex phase parameter α in the Lagrangian and the Dirac CP-violating phase δ to which neutrino oscillation experiments are sensitive. In particular, for the most interesting case of M ∼ 5 e 1 , we have δ −α. Imposing the constraints coming from lepton universality in τ decays and β → α γ and − α → − β + γ − δ searches, we have shown that 5 e 1 is the only case which does not require a tuning among the masses of the R and I scalars, due to a flavor suppression in the µ − e channel stemming from the U(1) flavor symmetry. The lepton universality constraint coming from the ratio |g µ /g e | selects m H + to be at most 300 GeV, with small and large tan β values. Nevertheless, it is worth emphasizing that these limits are loosen and wider intervals for m H + and tan β are allowed if the universality constraint is relaxed.
Among all the cases considered, only the pair (5 e 1 , 5 ν 8 ) leads to naturally small rates for LFV decays (in particular µ → eγ) without the need of cancellations between the amplitudes of R and I mediated contributions, only possible if m I m R . Such a condition among the scalar masses, together with a light H + , does not arise naturally in the 2HDM. In fact, m I m R can be easily accomplished in the 2HDM decoupling limit but, in general, H + is degenerate with R and I. Therefore, all texture combinations of Table VIII would be phenomenologically viable if one drops the requirement of having a non-negligible deviation from universality in τ decays.
In conclusion, we have shown that maximallyrestrictive lepton mass matrices realizable in the 2HDM with U(1) Abelian symmetries are phenomenologically viable from the point of view of lepton masses and mix-ing and the constraints imposed on flavour-changing processes. The results obtained in this work pose the natural question of what happens if we apply the same principle of maximally-restrictive textures to the quark sector. In this case, besides the symmetry implementation and compatibility with the observed quark masses and mixing, more severe constraints have to be checked, such as those coming from universality tests in τ and meson semileptonic decays, and from B → X s γ and meson µ + µ − decays. The extension of the analysis presented here to the quark sector is under preparation [85], as well as the generalization to the case of seesaw-generated Majorana neutrino masses.
To further motivate our analysis of texture-zero Yukawa and mass matrices in the 2HDM, let us analyze the consequences of imposing a U(1) symmetry in the SM such that the matrix M x contains texture zeros. First we note that, in the SM, Eq. (14) takes the form where we have set the SM Higgs U(1) phase to zero. Unlike in the 2HDM, the lepton mass matrices will be proportional to a single Yukawa matrix, leading to the following transformation properties of the Hermitian matrices H and H ν defined in Eq. (8): Using Eq. (10), and defining the unitary matrices A and B as one obtains that the lepton mixing matrix U must satisfy the relation On the other hand, the invariance condition (A2) reads which, in the case of non-degenerate charged-lepton and neutrino masses (as required by experiment), imply where A i and B i are general phases. Together with Eq. (A4), this leads to the relations Current experimental data (see Table I) indicate that all elements of U are nonzero. Therefore, compatibility of U with data forces the solution of the above equation to be and, consequently, Equation (A3) can be interpreted as the diagonalization of S . However, for fully degenerate eigenvalues, as demanded by Eq. (A9), S has a unique diagonal form, up to a rephasing of the columns. This leads to the lepton mixing matrix U = diag{e iχ1 , e iχ2 , e iχ3 }, which is obviously incompatible with the experimental observations.
On the other hand, one could require that S commutes with both U L and U ν L . Then A = B = S , implying that the charges A i = B i are in direct correspondence with the charges of S . Equation (A8) would then yield which is the only transformation that allows for a viable mixing matrix U in the context of the SM [84]. One can always set A 1 to zero so that S = 1 1, S e = diag(e iβ1 , e iβ2 , e iβ3 ), and S ν = diag(e iγ1 , e iγ2 , e iγ3 ). The phase transformation matrices would now be showing that, in the SM, imposing a single texture zero in M or M ν and requiring U to be compatible with experiment implies an entire column of zeros, resulting in a massless particle. This is clearly not acceptable for charged leptons. However, for neutrinos, such a possibility is not excluded, since a massless neutrino is compatible with current data. In conclusion, in the context of the SM, one cannot impose texture zeros in the chargedlepton mass matrix, while maintaining compatibility with all experimental observations. 9 The most economical scenario compatible with data would be the one with textures 9 This is in agreement with the conclusions of Ref. [84].
where P is a 3×3 permutation matrix (see Eq. (22)), and × denotes a general non-zero entry. In this case there are 23 independent parameters in the mass matrices, to be compared with 9 measurable quantities (5 masses, 3 mixing angles and 1 phase). The above analysis can be trivially extended to the SM quark sector being the CKM matrix the analogue of U. Since all quarks are known to be massive, neither the up nor the down quark mass matrix can have a zero eigenvalue. This means that, in the SM, no texture zeros can be imposed by Abelian flavor symmetries in the quark sector while ensuring compatibility with all experimental observations. or (B2) for distinct indices i, j leads to incompatible decompositions.
Symmetry transformations that allow for a full texture row or column in a mass matrix must respect Eq. (B1) or (B2), respectively, for the same pair of indices i, j. Since a set of transformations with α 1 = α 2 cannot generate a texture where the first and second rows are nonidentical, our conclusion regarding the lack of an implementation for 4 1 , 4 2 , 4 ν 17 , still holds for several symmetry transformations. The case where θ 1 = θ 2 in (15) should also be considered since, out of the several transformations, only one is required to act on the scalar fields. In this case, a full texture row or column leads to β 1 = β 2 = β 3 and α 1 = α 2 = α 3 , respectively. Thus, the imposition of a single texture zero leads to the requirement of an entire row or column of zeros, which is not compatible with the matrices 4 1 , 4 2 and 4 ν 17 . The above analysis allows us to generalize, to any number of consecutive symmetry transformations, the conclusion that these three textures cannot be realized through any continuous or discrete Abelian symmetries in the 2HDM. We remark that this result only holds in the context of the 2HDM, since the addition of a third Higgs doublet could lift these constraints. From Table II, we then conclude that 11 maximally-restrictive pairs of leptonic mass matrix textures, associated to 4 1 , 4 2 and 4 ν 17 , can be excluded from our model implementation perspective.
Next, we use the canonical method to eliminate the matrix pairs that are composed of textures that are individually realizable through Abelian symmetries, but for which the symmetry transformation cannot be implemented simultaneously in both matrices. First we note that, for each of the texture pairs (3 2 , 7 ν 1,3 ), (6 1 , 4 ν 1 ) and (5 1 , 5 ν 1,6 ), one of the textures is characterized by a column full of non-zero entries and it has two identical texture rows (generated by two identical left-handed continuous phases α i ). Since α i are common to the charged-lepton and Dirac-type neutrino sectors, Eq. (17) demands that both M and M ν have two identical texture rows i and j when α i = α j for some pair i, j. However, one finds that, for all these pairs, the second texture has no identical rows. As such, we conclude that these 5 pairs cannot be implemented through a continuous or discrete Abelian symmetry in a 2HDM.
Additional symmetry transformations with θ 1 = θ 2 do not change this conclusion since realizing the texture with a full column of non-zero entries would require α 1 = α 2 = α 3 , exacerbating the issue found with the implementation of these mass matrix pairs. All symmetry transformations attempting at generating these pairs must respect Eq. (B2) for the same pair i, j. Therefore, there is always at least one texture zero in these pairs which cannot be imposed, while allowing all of the non-zero entries. Thus, we can generalize to any number of symmetry transformations (15) the conclusion that (3 2 , 7 ν 1,3 ), (6 1 , 4 ν 1 ) and (5 1 , 5 ν 1,6 ) cannot be realized through any continuous or discrete Abelian symmetries in the 2HDM. Consequently, another 5 maximallyrestrictive pairs of leptonic mass matrix textures can be eliminated from Table II for the purpose of our model implementation. From the application of the canonical method one can conclude that, in the 2HDM, a mass matrix texture including a column with two non-zero entries associated to non-identical rows i and j can only be realized by a symmetry transformation which respects an appropriate version of In each of the pairs (4 3 , 6 ν 2,5,8 ) and (5 1 , 5 ν 5 ), the constituent textures can only be realized by transformations which obey non-compatible versions of Eq. (B3). Therefore, we conclude that these pairs cannot be realized through a continuous or discrete Abelian symmetry in the 2HDM.
Let us now consider the case of consecutive symmetry transformations. First we notice that, in the 2HDM, a texture including a column with two non-zero entries, in rows i and j, can only be implemented through a transformation which respects either or Eq. (B3). Two symmetry transformations that obey Eq. (B3) and Eq. (B4), respectively, for the same rows i, j generate incompatible decompositions into Yukawa matrices and cannot be used together to implement such a texture. If the rows i and j are non-identical then its implementation requires at least one symmetry transformation which obeys Eq. (B3) and, hence, all transformations used must do so. Furthermore, they must respect the same sign in ±(θ 1 − θ 2 ), as it determines the ordering of the decomposition into Yukawa matrices, which must always be the same. This is essentially the same conclusion as for the case of a single transformation. A symmetry transformation with θ 1 = θ 2 can only implement a texture with a column with two non-zero entries if it obeys Eq. (B4) and is, therefore, not useful in this context.
We conclude that the texture pairs (4 3 , 6 ν 2,5,8 ) and (5 1 , 5 ν 5 ) cannot be realized through any continuous or discrete Abelian symmetries, in the context of the 2HDM, regardless of the number of consecutive symmetry transformations of type (15) imposed. As a result, we exclude another 4 maximally-restrictive pairs of leptonic mass matrix textures from Table II. 4. Texture pairs (4 3 , 6 ν 4,6 ) and (5 1 , 5 ν 4 ) From the application of the canonical method one can conclude that if a mass matrix texture T 1 is characterized by two non-identical columns c 1 and c 2 with non-zero entries in some row i and zero entries in some row j, then it cannot be realized through Abelian symmetries in the 2HDM, together with a texture T 2 generated by a symmetry transformation obeying Eq. (B3) for the same rows i and j.
The above statement arises because a symmetry transformation as in Eq. (15) is unable to impose a zero in both (T 1 ) jc1 and (T 1 ) jc2 . Despite this, it is possible, in principle, to find transformations able to generate each zero separately, which could be an ideal setup for multiple symmetry transformations generating a texture pair unattainable through a single one. However, for the texture pairs considered here, this is not the case. In the context of the statement, one of the columns c 1 /c 2 of T 1 already has two non-zero entries in the rows i and k = j. Thus, a symmetry transformation imposing a zero in (T 1 ) jc2 /(T 1 ) jc1 generates a column c 1 /c 2 full of non-zero entries in T 1 . Such a transformation must have at least two identical charges α i . The only viable case is α k = α j (mod 2π). For each of the three texture pairs (4 3 , 6 ν 4,6 ) and (5 1 , 5 ν 4 ), we can identify the indices i, j, k and easily arrive at the conclusion that α k = α j (mod 2π) is not compatible with one of the mass matrix textures in the pair. The case θ 1 = θ 2 does not bear any consequence for the same arguments given in the previous section.
In conclusion, none of the pairs (4 3 , 6 ν 4,6 ) and (5 1 , 5 ν 4 ) can be realized through a continuous or discrete Abelian in the 2HDM, regardless of the number of consecutive symmetry transformations. Another 3 of the maximallyrestrictive pairs of leptonic mass matrix textures of Table II must be excluded from our model implementation.
In this appendix, we have determined that 23 out of the 28 maximally restrictive leptonic mass matrix texture pairs have no implementation through Abelian symmetries in the 2HDM. In the process of identifying nonrealizable texture pairs, we have also discussed the possibility of imposing several consecutive symmetry transformations of the form given in Eq. (15). It is worth emphasizing here that the same analysis can also be ex-tended to the realizable pairs, namely, (4 3 , 6 ν 1,3,7,9 ) and (5 1 , 5 ν 8 ). We find that no other decompositions (besides those given in Section V A) arise from the application of such consecutive transformations.
Finally, we stress that all our conclusions were drawn in the context of the 2HDM. The addition of more Higgs doublets would relax most of the constraints, even for the case of a single symmetry transformation.

Appendix C: Relations for M and Mν parameters
We have seen that for the maximally-restrictive textures the total number of parameters in the mass matrices is nine, equalling that of lepton mass and mixing parameters, namely three charged-lepton and two neutrino masses, three mixing angles and one CP-violating phase. For the mass matrices M and M ν , the number of free parameters not fixed by requiring the lepton masses to be compatible with experiment is two (one) for the 4 3 and 5 ν 8 (5 1 and 6 ν k ) textures. Those free parameters are determined by requiring compatibility with the measured values of the mixing angles and CP phase. In Table XI, we provide the defining relations of the parameters that can be written as functions of the masses and free parameters for the textures realizable with Abelian symmetries in the 2HDM (see Table VIII). The choice of the free parameters in the mass-defining relations is not unique but, obviously, the results do not depend on it. We also give the ranges within which the free parameters can vary. For each texture combination known to be compatible with data, the results shown in Section VI are obtained as follows.
The free parameters are varied within the ranges given in Table XI, considering all possible cases of lepton mass ordering and both the NH and IH cases. This gives rise to possible different scenarios labelled as i) and ii). Notice that, since the masses do not depend on the phase α, we consider α ∈ [0, 2π] . For each set of input parameters, the lepton mixing angles θ ij and the Dirac CP-violating phase δ in Eq. (11) are determined by computing U via Eqs. (9) and (10). For the texture 5 1 , U is determined by Eqs. (36)- (43). In the case 4 3 , an analytical form can also be obtained for U , but since the relations are more involved it is not presented here. As for U ν , we have in the 5 ν 8 case U ν = K α U ν , where K α = diag(e −iα , −iα , 1) and  ii) x 2 2 = ∆21 + ∆31 , y 2 = ∆31 − x 2 1 x 2 1 < ∆21 x 2 1 < ∆31 6 ν 9 : TABLE XI. Defining conditions for the charged-lepton and neutrino Yukawa textures given in Table VIII. For the 4 3 (5 1 ) texture we write a 2 1 , a 2 2 and b 2 1 (a 2 1 , b 2 1 and b 2 2 ) in terms of a 2 3 and b 2 2 (a 2 2 ) and the charged-lepton masses. In 5 1 , the state 1 can be identified with e, µ or τ leading, respectively, to the cases 5 e 1 , 5 µ 1 and 5 τ 1 discussed in Section VI. For the 6 ν k (5 ν 8 ) textures we write x 2 2 and y (x 2 1 and y 2 1 ) in terms of x 2 1 (x 2 2 and y 2 2 ) and the neutrino mass-squared differences ∆21 and ∆31. The Yukawa couplings in Y ,ν 1 (Y ,ν 2 ) entering Eq. (3) are determined by dividing ai and xi (bi and yi) by v cos β (v sin β). For simplicity, we use the notation ∆21 ≡ ∆m 2 21 and ∆31 ≡ ∆m 2 31 > 0 (∆m 2 31 corresponds to |∆m 2 31 | in the IH case).
for NH and IH neutrino mass spectra, respectively. In the above equations and in Table XI, we have used the notation ∆ ± = ∆ 31 ± ∆ 21 . Moreover, the −(+) sign in U ν corresponds to the case i (ii) of the validity ranges for x 2 and y 2 . The numerical procedure starts by varying the free parameters in their validity ranges. For each input set (a i , b i , x i , y i , α), the lepton mixing matrix U is obtained and compatibility with the 3σ ranges for θ ij given in Table I is checked. If compatibility is found, the phase δ is plotted as a function of α. This procedure fixes all mass parameters a i , b i , x i , y i and α in the mass matrices M and M ν . The Yukawa couplings in Y ,ν 1 (Y ,ν 2 ) are determined dividing the corresponding mass-matrix elements by v cos β (v sin β). Thus, for a given tan β, the couplings N e are completely determined in terms of lepton masses and mixing parameters, as illustrated in Section VI for the case (5 e 1 , 5 ν 8 ). Finally, the β → α γ and − α → − β + γ − δ BRs, as well as the value of |g e /g µ |−1, are computed.