Fermion mass hierarchies, large lepton mixing and residual modular symmetries

In modular-invariant models of flavour, hierarchical fermion mass matrices may arise solely due to the proximity of the modulus τ to a point of residual symmetry. This mechanism does not require flavon fields, and modular weights are not analogous to Froggatt-Nielsen charges. Instead, we show that hierarchies depend on the decomposition of field representations under the residual symmetry group. We systematically go through the possible fermion field representation choices which may yield hierarchical structures in the vicinity of symmetric points, for the four smallest finite modular groups, isomorphic to S3, A4, S4, and A5, as well as for their double covers. We find a restricted set of pairs of representations for which the discussed mechanism may produce viable fermion (charged-lepton and quark) mass hierarchies. We present two lepton flavour models in which the charged-lepton mass hierarchies are naturally obtained, while lepton mixing is somewhat fine-tuned. After formulating the conditions for obtaining a viable lepton mixing matrix in the symmetric limit, we construct a model in which both the charged-lepton and neutrino sectors are free from fine-tuning.


Introduction
Understanding the origins of flavour in both the quark and lepton sectors, i.e., the origins of the patterns of quark masses and mixing, of the charged-lepton and neutrino masses, of neutrino mixing and of the CP violation in the two sectors is one of the most challenging unresolved fundamental problems in particle physics [1]. 1 Within the reference three neutrino mixing scheme, the lepton flavour problem consists of three basic elements or sub-problems, namely, understanding: i) the origin of the hierarchical pattern of charged-lepton masses: m e m µ m τ , m e /m µ 1/200, m µ /m τ 1/17; ii) why neutrino masses m ν j are much smaller than the masses of charged leptons and quarks, m ν j ≪ m ,q , q = u, c, t, d, s, b and = e, µ, τ , with m ν j ∼ < 0.5 eV, m ≥ 0.511 MeV, m q ∼ > 2 MeV; iii) the origins of the patterns of neutrino mixing of 2 large and 1 small angles, and of the two independent neutrino mass squared differences, ∆m 2 21 |∆m 2 31 | with ∆m 2 21 /|∆m 2 31 | 1/30, where ∆m 2 ij ≡ m 2 i − m 2 j . Each of these three sub-problems is by itself a formidable problem. As a consequence, individual solutions to each of them have been proposed. The hierarchical pattern of charged-lepton masses can most naturally be understood within the Froggatt-Nielsen mechanism based on the U (1) FN flavour symmetry [2] and its extensions. The enormous disparity between the neutrino masses and the masses of the charged leptons and quarks can be understood within the seesaw or radiative models of neutrino mass generation or else employing the Weinberg effective operator idea [3] (for a concise review see, e.g., [4]). All these approaches lead naturally to massive Majorana neutrinos. Arguably the most elegant and natural explanation of the observed pattern of neutrino (or lepton) mixing of two large and one small mixing angles is obtained within the non-Abelian discrete symmetry approach to the problem (see, e.g., [5][6][7][8][9]).
In the case of the quark sector, the flavour problem similarly has two basic subproblems, namely, understanding: i) the origins of the hierarchies of the masses of the charge 2/3 and of the charge (−1/3) quarks; ii) the origins of the relatively small values of the three quark mixing angles.
The most natural qualitative solution of these two problems is arguably provided by the Froggatt-Nielsen approach [2], although the approach based on non-Abelian discrete symmetries has been applied to the quark flavour problem as well. Solutions to the two flavour problems within the theories with extra dimensions have also been proposed. 2 The specific solutions to the individual lepton flavour sub-problems listed above become problematic when applied to the sub-problems they were not intended to solve. The seesaw and the radiative neutrino mass models do not lead to understanding of the origin of the pattern of neutrino mixing without additional input, consisting typically of imposing specific additional symmetries (of GUT or flavour type) on the relevant constructions. Within the Froggatt-Nielsen approach one most naturally obtains small values of the three neutrino mixing angles and while the charged-lepton and quark mass hierarchies can be qualitatively understood within this approach, the specific predictions suffer from relatively large uncertainties. The symmetry breaking in the lepton and quark flavour models based on non-Abelian discrete symmetries is impressively cumbersome: it requires the introduction of a plethora of 'flavon' scalar fields having elaborate potentials, which in turn require the introduction of large shaping symmetries to ensure the requisite breaking of the symmetry leading to correct mass and mixing patterns.
There have been also attempts to make progress, e.g., on the lepton flavour problem by combining the proposed 'solutions' of the three related sub-problems. In these combined approaches it is difficult, if not impossible, to avoid the drawbacks of each of the 'ingredient' sub-problem 'solutions'. In some cases this can be achieved at the cost of severe fine-tuning. Thus, a universal, elegant, natural and viable theory of flavour that is free from undesired drawback features is still lacking. Constructing such a theory would be a major breakthrough in particle physics.
The unsatisfactory status of the flavour problem and the remarkable progress made in the studies of neutrino oscillations (see, e.g., [10]), which began 22 years ago with the discovery of oscillations of the atmospheric ν µ andν µ by the Super-Kamiokande experiment [11] and lead, in particular, to the determination of the pattern of neutrino mixing, stimulated renewed attempts to seek alternative viable approaches to the lepton as well as to the quark flavour problems. A step in this direction was made in [12] where the idea of using modular invariance as a flavour symmetry was put forward. This new original approach based on modular invariance opened up a new promising direction in the studies of the flavour problem and correspondingly in flavour model building.
The main feature of the approach proposed in [12] is that the elements of the Yukawa coupling and fermion mass matrices in the Lagrangian of the theory are modular forms of a certain level N which are functions of a single complex scalar field τ -the modulus -and have specific transformation properties under the action of the modular group. In addition, both the couplings and the matter fields (supermultiplets) are assumed to transform in representations of an inhomogeneous (homogeneous) finite modular group Γ ( ) N . For N ≤ 5, the finite modular groups Γ N are isomorphic to the permutation groups S 3 , A 4 , S 4 and A 5 (see, e.g., [13]), while the groups Γ N are isomorphic to the double covers of the indicated permutation groups, S 3 ≡ S 3 , A 4 ≡ T , S 4 and A 5 . These discrete groups are widely used in flavour model building. The theory is assumed to possess the modular symmetry described by the finite modular group Γ ( ) N , which plays the role of a flavour symmetry. In the simplest class of such models, the VEV of the modulus τ is the only source of flavour symmetry breaking, such that no flavons are needed. Another appealing feature of the proposed framework is that the VEV of τ can also be the only source of breaking of the CP symmetry [14]. When the flavour symmetry is broken, the elements of the Yukawa coupling and fermion mass matrices get fixed, and a certain flavour structure arises. As a consequence of the modular symmetry, in the lepton sector, for example, the charged-lepton and neutrino masses, neutrino mixing and the leptonic CPV phases are simultaneously determined in terms of a limited number of coupling constant parameters. This together with the fact that they are also functions of a single complex VEV -that of the modulus τ -leads to experimentally testable correlations between, e.g., the neutrino mass and mixing observables. Models of flavour based on modular invariance have then an increased predictive power.
The modular symmetry approach to the flavour problem has been widely implemented so far primarily in theories with global (rigid) supersymmetry. Within the SUSY framework, modular invariance is assumed to be a feature of the Kähler potential and the superpotential of the theory. 3 Bottom-up modular invariance approaches to the lepton flavour problem have been exploited first using the groups Γ 3 A 4 [12,16], Γ 2 S 3 [17], Γ 4 S 4 [18]. After the first studies, the interest in the approach grew significantly and models based on the groups Γ 4 S 4 [19][20][21][22][23][24][25][26] [53] have been constructed and extensively studied. Similarly, attempts have been made to construct viable models of quark flavour [54] and of quark-lepton unification [55][56][57][58][59][60][61][62]. The formalism of the interplay of modular and gCP symmetries has been developed and first applications made in [14]. It was explored further in [63][64][65], as was the possibility of coexistence of multiple moduli [66][67][68][69], considered first phenomenologically in [19,30]. Such bottom-up analyses are expected to eventually connect with top-down results  based on ultraviolet-complete theories. While the aforementioned finite quotients Γ N of the modular group have been widely used in the literature to construct modular-invariant models of flavour from the bottomup perspective, top-down constructions typically lead to their double covers Γ N (see, e.g., [73,75,76,92]). The formalism of such double covers has been developed and viable flavour models constructed in Refs. [93], [94,95] and [96,97] for the cases of Γ 3 T , Γ 4 S 4 and Γ 5 A 5 , respectively.
In almost all phenomenologically viable flavour models based on modular invariance constructed so far the hierarchy of the charged-lepton and quark masses is obtained by fine-tuning some of the constant parameters present in the models. 4 Perhaps, the only notable exceptions are Refs. [98,99], in which modular weights are used as Froggatt-Nielsen charges, and additional scalar fields of non-zero modular weights play the role of flavons.
In the present article we develop the formalism that allows to construct models in which the fermion (e.g. charged-lepton and quark) mass hierarchies follow solely from the properties of the modular forms present in the fermion mass matrices, thus avoiding the fine-tuning without the need to introduce extra fields. We consider theories described by a modular group Γ N with N ≤ 5 (which encompasses the unprimed Γ N ). It was noticed in [19] and further exploited in [27,30,64] that for the three fixed points of the VEV of τ in the modular group fundamental domain, τ sym = i, τ sym = ω ≡ exp(i 2π/3) = − 1/2 + i √ 3/2 (the 'left cusp'), and τ sym = i∞, the theories based on the Γ N invariance have respectively Z S 2 , Z ST 3 , and Z T N residual symmetries. In the case of the double cover groups Γ N , the Z S 2 residual symmetry is replaced by the Z S 4 and there is an additional Z R 2 symmetry that is unbroken for any value of τ (see [94] for further details). The indicated residual symmetries play a crucial role in our analysis.
The fermion mass matrices are strongly constrained in the points of residual symmetries. This suggests that fine-tuning could be avoided in the vicinity of these points if the charged-lepton and quark mass hierarchies follow from the properties of the modular forms present in the corresponding fermion mass matrices rather than being determined by the values of the accompanying constants also present in the matrices. Relatively small deviations of the modulus VEV from the symmetric point might also be needed to ensure the breaking of the CP symmetry [14].
We note that in [100] flavour models in the vicinity of the residual symmetry fixed points, τ sym = i, ω, i∞, have been investigated within the modular invariant A 4 framework (N = 3). The authors find viable lepton (quark) flavour models in the vicinity of each of three residual symmetry values of τ sym (of τ sym = i), in which the mixing arises seemingly without fine-tuning. At the same time, the charged-lepton and quark mass hierarchies are obtained by fine-tuning the values of the constants present in the respective mass matrices.
The aim of this study is to investigate the possibility of obtaining fermion mass hierarchies -and, in models of lepton flavour, large mixing -without fine-tuning. The article is structured as follows. After introducing the necessary tools in section 2, we describe how one can naturally generate hierarchical mass patterns in the vicinity of symmetric points in section 3.1. In section 3.2, the role of decompositions under the residual symmetry groups is highlighted. We perform a systematic scan of attainable hierarchical patterns for N ≤ 5, the results of which are reported in section 3.3. The analysis of two promising lepton flavour models in section 3.4 motivates the discussion, in section 4.1, of necessary conditions to avoid fine-tuned leptonic mixing. We are then driven to a subset of viable models, the most promising of which is explored in section 4.2. We summarise our results and conclude in section 5.

Modular symmetries as flavour symmetries
We start by briefly reviewing the modular invariance approach to flavour. In this supersymmetric (SUSY) framework, one introduces a chiral superfield, the modulus τ , transforming non-trivially under the modular group Γ ≡ SL(2, Z). The group Γ is generated by the matrices Elements γ of the modular group act on τ via fractional linear transformations, while matter superfields transform as 'weighted' multiplets [12,92,101], where k ∈ Z is the so-called modular weight 5 and ρ is a unitary representation of Γ. In using modular symmetry as a flavour symmetry, an integer level N ≥ 2 is fixed and one assumes that ρ(γ) = 1 for elements γ of the principal congruence subgroup Hence, ρ is effectively a representation of the (homogeneous) finite modular group Γ N ≡ Γ Γ(N ) SL(2, Z N ). For N ≤ 5, this group admits the presentation The (lowest component of the) modulus τ acquires a VEV which is restricted to the upper half-plane and plays the role of a spurion, parameterising the breaking of modular invariance. Additional flavon fields are not required, and we do not consider them here. Since τ does not transform under the R generator, a Z R 2 symmetry is preserved in such scenarios. If also matter fields transform trivially under R, one may identify the matrices γ and −γ, thereby restricting oneself to the inhomogeneous modular group Γ ≡ P SL(2, Z) ≡ SL(2, Z) / Z R 2 . In such a case, ρ is effectively a representation of a smaller (inhomogeneous) finite modular group Γ N ≡ Γ Γ(N ) ∪ Z R 2 . For N ≤ 5, this group admits the presentation In general, however, R-odd fields may be present in the theory and Γ and Γ N are then the relevant symmetry groups. As shown in Table 1, the finite modular groups Γ N N 2 3 4 5 dim M k (Γ(N )) k/2 + 1 k + 1 2k + 1 5k + 1 Table 1: Finite modular groups and dimensionality of the corresponding spaces of modular forms, for N ≤ 5. Note that for N = 2 only even-weighted modular forms exist. and Γ N are isomorphic to permutation groups and to their double covers for small N . Group-theoretical results for the Γ N groups are collected in appendix B of [14], while for the double cover groups Γ N they can be found in Refs. [93,94,96].
Finally, to understand how modular symmetry may constrain the Yukawa couplings and mass structures of a model in a predictive way, we turn to the Lagrangian -which for an N = 1 global supersymmetric theory is given by Here K is the Kähler potential, while the superpotential W can be expanded in powers of matter superfields ψ I , where one has summed over all possible field combinations and independent singlets of the finite modular group. By requiring the invariance of the superpotential under modular transformations, 6 one finds that the field couplings Y I 1 ...In (τ ) have to be modular forms of level N . These are severely constrained holomorphic functions of τ , which under modular transformations obey Modular forms carry weights k = k I 1 +. . .+k In and furnish unitary representations ρ Y of the finite modular group such that ρ Y ⊗ ρ I 1 ⊗ . . . ⊗ ρ In ⊃ 1. Non-trivial modular forms of a given level exist only for k ∈ N, span finite-dimensional linear spaces M k (Γ(N )), and can be arranged into multiplets of Γ ( ) N . The fact that these spaces have low dimensionalities for small values of k and N (as shown in Table 1) is at the root of the predictive power of the described setup, since only a restricted number of τ -dependent Yukawa textures are allowed in the superpotential.
Note that modular forms are functions of τ and are thus invariant under R. In order to compensate the (−1) k factor in eq. (2.9), odd-weighted forms must furnish representations with ρ Y (R) = −1 (we use hats to denote such representations). For even-weighted modular forms, one has instead ρ Y (R) = 1.

Residual symmetries
The breakdown of modular symmetry is parameterised by the VEV of the modulus and there is no value of τ which preserves the full symmetry. Nevertheless, at certain socalled symmetric points τ = τ sym the modular group is only partially broken, with the unbroken generators giving rise to residual symmetries. Recall that the R generator is unbroken for any value of τ , so that a Z R 2 symmetry is always preserved. The fundamental domain D of the modular group is shown in Figure 1, along with its symmetric points. There are only three inequivalent symmetric points, namely [19]: Finally, it is worth noting that these symmetric values preserve the CP symmetry of a CP-and modular-invariant theory (e.g. a modular theory where the couplings satisfy a reality condition) [14,94]. A Z CP 2 symmetry is preserved for Re τ = 0 or for τ lying on the border of D, but is broken at generic values of τ .
3 Mass hierarchies without fine-tuning

Mass matrices close to symmetric points
In theories where modular invariance is broken only by the modulus, the flavour structure of mass matrices in the limit of unbroken supersymmetry is determined by the value of τ and by the couplings in the superpotential. At a symmetric point τ = τ sym , flavour textures can be severely constrained by the residual symmetry group, which may enforce the presence of multiple zero entries in the mass matrices. As τ moves away from its symmetric value, these entries will generically become non-zero. The magnitudes of such (residual-)symmetry-breaking entries will be controlled by the size of the departure from τ sym and by the field transformation properties under the residual symmetry group (which may depend on the modular weights). This is shown in what follows.
Consider a modular-invariant bilinear where the superfields ψ and ψ c transform under the modular group as 7 so that each M (τ ) ij is a modular form of level N and weight K ≡ k + k c . Modular invariance requires M (τ ) to transform as Taking τ to be close to the symmetric point, and setting γ to the residual symmetry generator, one can use this transformation rule to constrain the form of the mass matrix M (τ ). We consider each of the three symmetric points in turn.

τ sym = i∞
The representation basis for the group generators S and T typically found in the literature is the T -diagonal basis, in which ρ (c) (T ) = diag(ρ (c) i ). This basis is particularly useful for the analysis of models where τ is 'close' to τ sym = i∞, i.e. models with large Im τ . By setting γ = T in eq. (3.3), one finds It is convenient to treat the M ij as a function of q ≡ exp (2πiτ /N ), so that parameterises the deviation of τ from the symmetric point. Note that the entries M ij (q) depend analytically on q and that q Expanding both sides in powers of q, one finds where M (n) ij denotes the n-th derivative of M ij with respect to q. This means that in the vicinity of the symmetric point. It crucially follows that the entry M ij is expected to be O( l ) whenever Im τ is large. The power l only depends on how the representations of ψ and ψ c decompose under the residual symmetry group Z T N . This point will be made explicit in section 3.2.

τ sym = i
For the analysis of models where τ is in the vicinity of τ sym = i, it is convenient to switch to the basis where the S generator is represented by a diagonal matrix. In this which not only simplify the algebra, but also correspond to representations of the residual symmetry group, see eq. (3.23). By setting γ = S in eq. (3.3), one finds We now treat the M ij as functions of so that, in this context, ≡ |s| parameterises the deviation of τ from the symmetric point. Note that the entries M ij (s) depend analytically on s and that s S − → −s. Thus, in terms of s, eq. (3.10) reads ij denotes the n-th derivative ofM ij with respect to s. It should be clear from eq. (3.13) that for τ i the mass matrix entry M ij ∼M ij is only allowed to be O(1) whenρ c iρ j = 1. If insteadρ c iρ j = −1, the entry M ij ∼M ij is expected to be O( ), with = |s|. Note that, unlike in the previous section, the relevant factorsρ (c) i depend on the weights k (c) via eq. (3.9).

τ sym = ω
Finally, for the analysis of models where τ is in the vicinity of τ sym = ω, we consider the basis where the product ST is represented by a diagonal matrix.
which are representations under the residual symmetry group, see eq. (3.24). By setting It is now convenient to treat the M ij as functions of . Expanding both sides in powers of u, one obtains ij denotes the n-th derivative ofM ij with respect to u. It follows that for τ ω the mass matrix entry M ij ∼M ij is only allowed to be O(1) whenρ c iρ j = 1. More generally, ifρ c iρ j = ω l with l = 0, 1, 2, then the entry M ij ∼M ij is expected to be O( l ) in the vicinity of τ = ω, with = |u|. Like in the previous section, the factorsρ (c) i depend on the weights k (c) , see eq. (3.14).

Decomposition under residual symmetries
We have just shown that, as τ departs from a symmetric value τ sym , entries of the mass matrices are subject to corrections of O( l ), where parameterises the deviation of τ from τ sym . The powers l are extracted from products of factors which, in this section, are shown to correspond to representations of the residual symmetry group. One can systematically identify these residual symmetry representations for the different possible choices of Γ N representations of matter fields. This knowledge will later be exploited to construct hierarchical mass matrices via controlled corrections to entries which are zero in the symmetric limit.
We start by noting that matter fields ψ furnish 'weighted' representations (r, k) of the finite modular group Γ N . Whenever a residual symmetry is preserved by the value of τ , matter fields decompose into unitary representations of the residual symmetry group. Modulo a possible Z R 2 factor, 9 these groups are the cyclic groups Z T N , Z S 4 , and Z ST 3 (cf. section 2.2). A cyclic group Z n ≡ a | a n = 1 has n inequivalent 1-dimensional irreducible representations (irreps) 1 k , where k = 0, . . . , n − 1 is sometimes referred to as a charge. The group generator a is represented by one of the n-th roots of unity, For odd n, the only real irrep of Z n is the trivial one, 1 0 (the reality of an irrep is indicated by removing the boldface). For even n, there is one more real irrep, 1 n/2 . All other irreps are complex, and split into pairs of conjugated irreps: To illustrate the aforementioned decomposition of representations at symmetric points, we take as an example a (3, k) triplet ψ of S 4 . It transforms under the unbroken One can check that the eigenvalues of ρ 3 (ST ) are 1, ω and ω 2 , and so in a suitable (ST -diagonal) basis the transformation rule explicitly reads which means that ψ decomposes as ψ 1 k ⊕ 1 k+1 ⊕ 1 k+2 under the residual Z ST 3 . One can find the residual symmetry representations for any other 'weighted' multiplet of a finite modular group in a similar fashion. For a given level N , the decompositions of fields under a certain residual symmetry group only depend on the pair (r, k). In general: where for the last equality we have assumed to be in a T -diagonal basis (no sum over i). The phase factors ρ i correspond to the Z T N irreps into which ψ decomposes. It follows that each ρ i is a power of ζ = exp(2πi/N ), depending on r but not on k.
• At τ = i, ψ ∼ (r, k) transforms under the unbroken γ = S as where for the last equality we have assumed to be in an S-diagonal basis (no sum over i). The phase factorsρ i = i k ρ i correspond to the Z S 4 irreps into which ψ decomposes. It follows that eachρ i is a power of i which depends both on r and on k (mod 4).
where for the last equality we have assumed to be in an ST -diagonal basis (no sum over i), as in the explicit example of eq. (3.21). The phase factorsρ i = ω k ρ i correspond to the Z ST 3 irreps into which ψ decomposes. It follows that eachρ i is a power of ω which depends both on r and on k (mod 3).
After identifying the ( ) ρ i and ( ) ρ c i factors for the fields ψ and ψ c entering a bilinear (equivalently, their irrep decompositions), one can apply the results of the previous section to determine the structure of a mass matrix in the vicinity of a symmetric point in terms of powers of , in the appropriate basis. It follows from the above that, in the analysis with large Im τ , the product (ρ c i ρ j ) * matches some power ζ l with 0 ≤ l < N , while in the analysis corresponding to τ ω one necessarily hasρ c iρ j = ω l with l = 0, 1, 2. These were tacitly taken as the most general possibilities in sections 3.1.1 and 3.1.3. The same reasoning implies that, in the τ i context,ρ c iρ j is some integer power i l , with l = 0, 1, 2, 3. It turns out that only two out of the four possibilities are viable, namely l = 0, 2 so thatρ c iρ j = ±1, as considered in section 3.1.2. This is due to the fact that M (τ ) ij is R-even and thus the fields ψ c i and ψ j need to carry the same R-parity (see also appendix A).
We list in Tables 6 -9 of appendix A the decompositions of the weighted representations of Γ N (N ≤ 5) under the three residual symmetry groups, i.e. the residual decompositions of the irreps of Γ 2

From entries to masses
We are in a position to use the results found so far and construct hierarchical mass matrices in the vicinity of a symmetric point. We have seen that in an appropriate basis M (τ ( )) ij ∼ O( l ). For each (i, j) pair, the power l can be obtained from the residual symmetry group decompositions of Tables 6 -9.
Note that a modular-symmetric mass matrix M (τ ( )) depends analytically on the small real parameter , defined in section 3.1 for each symmetric point. Physical masses are the singular values of M (τ ) and are also analytic functions of . 10 After SUSYbreaking, the leading superpotential contribution to each fermion mass is thus expected to be proportional to a power of which depends on the hierarchical structure of the entries of M . To find out which, one can make use of the following set of relations, valid for any n × n complex matrix M [106]: where p = 1, . . . , n is fixed, m i are the singular values of M , and the sum on the right-hand side goes over all possible p × p submatrices M p×p of M . In the particular case of n = 3, we denote the masses by (m 1 , m 2 , m 3 ) such that their leading terms are respectively of order where ∼ refers to power counting in and not necessarily to a reliable approximation. Note that so far the considered mass spectrum is generic. This is to be contrasted with the special case of a hierarchical 3 × 3 mass matrix, for which d 1 > d 2 > d 3 ≥ 0 and thus m 1 m 2 m 3 . In this case, eqs. (3.26) turn into useful approximations, and lead to reliable expressions for m 3 , m 2 and m 1 = | det M |/(m 2 m 3 ).

Example and results
As an example of application of our results, consider a model at level N = 5 with τ having a large imaginary part and with matter fields in weighted representations ψ ∼ (3, k) and ψ c ∼ (3 , k c ). From Table 9 one sees that ψ 1 0 ⊕ 1 1 ⊕ 1 4 and ψ c 1 0 ⊕ 1 2 ⊕ 1 3 under the residual group at the symmetric point τ sym = i∞. One can then identify Resorting to (3.26), one finds that spectrum is hierarchical, with (m 3 , m 2 , m 1 ) ∼ (1, , 4 ). Note that to have a non-zero mass matrix one needs the sum K = k + k c to be even (in this case), since matter fields furnish unhatted representations of the finite modular group and should carry the same R-parity (see appendix A). Furthermore, in order to obtain the full structure of eq. (3.28) and the expected hierarchical spectrum, K must be large enough that sufficient modular forms contribute to M (τ ). For instance, for K = 2 the superpotential may turn out to include a unique contribution: denotes the modular form multiplet of level N , weight K and irrep r, with µ possibly labelling linearly independent multiplets of the same type. The rightmost matrix subscript indicates the multiplet to which the Y i components belong. We have considered the T -diagonal basis for A 5 . One can explicitly check that, at leading order in = |q|, the components of Y 4 , with q = exp (2πiτ /5) and a common normalisation N . The power structure matches that of eq. (3.28) and naïvely this corresponds to the desired (1, , 4 ) spectrum. Upon closer inspection, however, one realises that the determinant of M vanishes identically for any value of τ , meaning that at least one fermion is massless. In the vicinity of τ sym = i∞, we have (m 3 , m 2 , m 1 ) ∼ (1, , 0). This issue is resolved already at weight K = 4, for which the modular multiplets Y Let us pause and describe our philosophy going forward. We are interested in identifying 3 × 3 hierarchical mass matrices where the hierarchical pattern is a result of the proximity of the modulus to a point of residual symmetry and no massless fermions are present in the spectrum. We assume to effectively be dealing with bilinears of the type (3.1) and consider all possible 3-dimensional representations for the fields ψ and ψ c . While the representations r and r c are in general reducible, we focus on the case where the same weight is shared between the irreps into which they decompose. 11 Thus, in our search, we take r (c) to be either irreducible or a direct sum of irreps sharing the same ρ(R). While it is possible for r (c) to be a direct sum of hatted and unhatted representations, the requirement of a common weight k (c) would result in the co-existence of R-odd and R-even fields within ψ (c) . The fact that M (τ ) is R-even would then imply the isolation of these sectors by the Z R 2 symmetry and the vanishing of some mixing angles.
Finally, it is straightforward to recognise that if all mass matrix entries are either O(1) or O( ), then leading contributions to the masses themselves are not expected to be smaller than O( ), unless one resorts to cancellations. Therefore, for τ i one cannot produce the desired hierarchical patterns solely as a consequence of the smallness of .
The result of our analysis is given in Tables 10 -13 of appendix B. These tables summarise, for each of the levels N ≤ 5, the patterns which may arise in the vicinity of the two potentially viable symmetric points, τ sym = ω and i∞, for all (r, r c ) pairs of 3-dimensional representations and all weights k (c) . One finds that it is only possible to obtain hierarchical spectra for a small list of representation pairs, the most promising of which are collected here, in Table 2. Table 2: Hierarchical mass patterns which can be realised in the vicinity of symmetric points. These patterns are unaffected by the exchange r ↔ r c and may only be viable for certain weights (see appendix B). Subscripts run over irreps of a certain dimension, and 1 a = 1 a for N = 3, while 1 a = 1 a for N = 4. Primes in parenthesis are uncorrelated.
We have excluded from this summary table reducible representations made up of three copies of the same singlet, since in those cases at least three independent modular multiplets of the same type must contribute to the mass matrix to avoid a massless fermion, and the number of superpotential parameters is unappealingly high. Still, such cases may result in other interesting hierarchical patterns such as (1, 2 , 3 ) and ( , 2 , 3 ) and can be found in the tables of appendix B.
Inspired by these results, we now turn to the construction of realistic models of lepton flavour where mass hierarchies are a consequence of the described mechanism.  Table 3: Best-fit values and 1σ ranges for neutrino oscillation parameters obtained from the global analysis [107], and for charged-lepton mass ratios, given at the scale 2 × 10 16 GeV with the tan β averaging described in [12], obtained from Ref. [108]. The parameters entering the definition of r are δm 2 ≡ m 2 2 − m 2 1 and ∆m 2 ≡ m 2 3 − (m 2 1 + m 2 2 )/2.

Promising models
To ascertain the viability of modular-invariant models of lepton flavour one must confront their predictions with experimental data on ratios of charged-lepton masses, neutrino mass-squared differences and leptonic mixing angles, see Table 3 (the constraints on the Dirac CPV phase δ are ignored in our fit). The reader is referred to [19] for details on our numerical procedure. The minimal-form Kähler potential is here considered, with Λ 0 of mass dimension one. We further take Higgs doublets H u and H d to be singlets under the modular group. Charged-lepton masses are obtained from their Yukawa interactions, where L and E c denote lepton doublet and charged-lepton singlet superfields with weights k L and k E , respectively. Neutrino masses are generated either by the Weinberg operator, 33) or within a type-I seesaw UV completion, where at least 2 neutrino gauge-singlet superfields N c of weight k N are present in the model. To compensate the modular weights of field monomials, the modular forms entering the Weinberg and Majorana terms need to have weights k W = 2k L and k M = 2k N , while those in Yukawa terms need instead k Y = k L + k E and k Y = k L + k N . The relevant superpotentials can be cast in the form . .

(3.37)
Finally, aiming at minimal and predictive examples, we impose a generalised CP symmetry on the models, enforcing the reality of coupling constants [14].

A 5 models with L ∼ 3, E c ∼ 3
We start by considering the most 'structured' series of hierarchical models, i.e. the case with both fields L, E c furnishing complete irreducible representations of the finite modular group. According to Table 2, the only such possibility arises at level N = 5 in the vicinity of τ = i∞ when L and E c are different triplets of the finite modular group A 5 . For definiteness, we choose L ∼ (3, k L ), E c ∼ (3 , k E ). The predicted charged-lepton mass pattern is (m τ , m µ , m e ) ∼ (1, , 4 ). We have performed a systematic scan restricting ourselves to models requiring modular forms of weight not higher than k = 6, and involving at most 8 effective parameters (including Re τ and Im τ ). Models producing a massless electron are rejected. For neutrino masses generated via a type I seesaw, we have considered gauge-singlet superfields N c furnishing a complete irrep of dimension 2 or 3. Out of the 36 models thus identified, we have selected the only one which i) is viable in the regime of interest, ii) produces a charged-lepton spectrum which is not fine-tuned, 12 and iii) is consistent with the experimental bound on the Dirac CPV phase δ. For this model, k L = 3, k E = 1 and N c ∼ (2 , 2). The corresponding superpotential reads: The modular forms entering W are obtained from the lowest weight (k = 1) sextet [97], where θ and ε are functions of the modulus, θ(τ ) = 1 + 3q 5 /5 + 2q 10 /25 − 28q 15 /125 + . . . and ε(τ ) = q 1 − 2q 5 /5 + 12q 10 /25 + 37q 15 /125 + . . . , with q = exp (2πiτ /5). In our regime of interest, |q| is small and the sextet modular multiplet approximately reads Y (1, 2, 5q, 0, 0, 0), which motivates an alternative definition of expansion parameter, | | with ≡ 5ε/θ 5q, used only in the context of this section. The charged-lepton mass matrix is then approximated by at leading order in | |. These expressions alone isolate viable ( -independent) regions in the plane ofα 2 andα 3 . Taking the 1σ ranges for charged-lepton mass ratios from Table 3, we plot these regions in Figure 2. The superimposed contours refer to the Barbieri-Giudice measure of fine-tuning [109] in the charged-lepton sector, max(BG), corresponding to the largest of four quantities |∂ ln(mass ratio)/∂ lnα 2,3 |. Expansions similar to (3.40) for the neutrino Yukawa and mass matrices are not useful since | | is not the only small parameter in the neutrino sector. In particular, some entries of M ν are proportional to 1 + (5/12)g 2 /g 1 − √ 3 g 3 /g 1 which is forced by the fit to be O(10 −2 ) in the viable region. Fine-tuning in the neutrino sector is expected (see section 4.1) and is related to the fact that the residual symmetry constrains not only the charged-lepton masses, but also the lepton mixing pattern. By sending Im τ → ∞ while keeping the couplings fixed to their best-fit values, one can check that the (13), (23), (31) and (32) entries of the PMNS matrix become zero. In particular, sin 2 θ 23 → 0, and therefore the symmetric limit does not provide an approximate description of leptonic mixing. The result of the fit of this A 5 model is summarised in the first column of Table 5, which includes best-fit values and the corresponding 3σ ranges. The viable region in the τ plane is shown on the left side of Figure 5 and corresponds to a neutrino spectrum with inverted ordering (IO).

S 4 models with L ∼2 ⊕1, E c ∼3
In the second most 'structured' case, one of the fields L, E c is an irreducible triplet, while the other decomposes into a doublet and a singlet of the finite modular group. This possibility is realised at level N = 4 in the vicinity of τ = i∞, see Table 2. For defi-niteness, we take L = L 12 ⊕ L 3 with L 12 ∼ (2, k L ), L 3 ∼ (1, k L ), and E c ∼ (3 , k E ). The charged-lepton mass pattern in this regime is predicted to be (m τ , m µ , m e ) ∼ (1, , 3 ).
We have performed a systematic scan restricting ourselves to models involving at most 8 effective parameters (including Re τ and Im τ ; no limit on modular form weights). Once again, models predicting a massless electron are rejected, while the N c (when present) furnish a complete irrep of dimension 2 or 3. Out of the 60 models thus identified, we have selected the only one which i) is viable in the regime of interest and ii) produces a charged-lepton spectrum which is not fine-tuned. This model turns out to be consistent with the experimental bound on the Dirac CPV phase. It corresponds to k L = k E = 2 and N c ∼ (3, 1) and the superpotential reads: The modular forms entering W are obtained from the lowest weight (k = 1) triplet [94], where θ(τ ) = 1 + 2q 4 + 2q 16 + . . . and ε(τ ) = 2q + 2q 9 + 2q 25 + . . ., with q = exp (iπτ /2). Using in the context of this section the expansion parameter obtained from the ratio of these functions, | | with ≡ ε/θ 2q, the charged-lepton mass matrix is approximately given by withα 2 ≡ α 2 /α 1 andα 3 ≡ α 3 /α 1 . It matches the expected power structure in |q|, as one can check. 13 One can also find approximate expressions for the charged-lepton mass ratios, which read m e m µ 18 The expressions (3.45) isolate viable ( -independent) regions in the plane of coupling constants, sayα −1 2 = α 1 /α 2 andα 3 /α 2 = α 3 /α 2 . We plot these regions in Figure 3, 13 Aside from consulting the third column of Table 8, one must keep in mind the ordering of the ρ (c) i , which depends on the representation basis for the group generators (we are using that of [94]). including contours quantifying the degree of fine-tuning involved in the relation between charged-lepton mass ratios and superpotential parameters (as described in the previous section). Note that the model best-fit point in particular corresponds to a small value of max(BG) 0.74.
Regarding the neutrino sector, one can check that close to the symmetric limit the neutrino masses follow the pattern (a| | 6 , b 0 − b 2 | | 2 , b 0 + b 2 | | 2 ), which naturally leads to IO. However, we were unable to find viable regions with IO and, instead, the model predicts a neutrino spectrum with normal ordering (NO). Within the viable region, the approximate pattern (a| is not accurate since | | is not the only small parameter in the neutrino sector. In particular, some of the entries of M ν are proportional to (1+ √ 6g 2 /g 3 ) which is forced to be O(10 −3 ) in our fit. As in the previous section, the fine-tuning in the neutrino sector is explained by the necessity to introduce large corrections to the symmetric-limit PMNS matrix, which has a zero either in the (31) or (33) entry in the case of NO or IO, respectively.
The result of the fit of this S 4 model is summarised in the second column of Table 5. The viable region in the τ plane is located near the point τ = 2.7i and is shown on the left side of Figure 5. 4 Large mixing angles without fine-tuning 4

.1 Viable PMNS matrix in the symmetric limit
We have seen in the previous sections that a slightly broken residual modular symmetry allows to accommodate hierarchical charged-lepton masses without fine-tuning of the corresponding couplings. However, the resulting models are still subject to fine-tuning in the neutrino sector, since residual symmetries typically constrain not only the chargedlepton masses, but also the form of the PMNS matrix by forcing some of its entries to be zeros. This raises the question of whether it is possible to have a PMNS matrix which is close to the observed one even in the symmetric limit, i.e. such that either none of its entries vanish, or only the (13) entry vanishes as → 0.
This possibility has been investigated in Ref. [110] for arbitrary flavour symmetry groups. In particular, this analysis directly applies to the case of the flavour symmetry being a residual modular symmetry. One of the main conclusions of Ref. [110] is that only a limited number of flavour symmetry representation choices for L and E c give rise to a PMNS matrix which is viable in the symmetric limit (as defined above). Most notably, there are only two such cases consistent with hierarchical charged-lepton masses: where 1 is some real singlet of the flavour symmetry, and r is some (possibly reducible) representation such that r ⊃ 1; where 1 is some complex singlet of the flavour symmetry, 1 * is its conjugate, and r is some (possibly reducible) representation such that r ⊃ 1, 1 * .
The above original result makes use of the assumption that one charged-lepton mass and at least one neutrino mass does not vanish in the symmetric limit. However, one can also deduce from the analysis performed in [110] that the PMNS matrix is generically unconstrained in the symmetric limit when the opposite is true. Therefore, we extend the list of viable cases with the following two: 3. all charged-lepton masses vanish in the symmetric limit, i.e. the corresponding hierarchical pattern involves only positive powers of , e.g. ( , 2 , 3 ); 4. all light neutrino masses vanish in the symmetric limit, i.e. L decomposes into three (possibly identical) complex singlets none of which are conjugated to each other.
It follows that a modular-symmetric model of lepton flavour with hierarchical chargedlepton masses may be free of fine-tuning if it satisfies any of the properties 1-4. Applying this filter to the promising hierarchical cases of Table 2, one is left with the representation pairs listed here, in Table 4. In this summary table, we have once again disregarded reducible representations made up of three copies of the same singlet.
We now proceed by constructing such a model in the following section (clearly, the models described in section 3.4 do not satisfy any of the properties 1-4). As a final Table 4: Hierarchical charged-lepton mass patterns which may be realised in the vicinity of symmetric points without fine-tuned mixing (PMNS close to the observed one in the symmetric limit). The property which is satisfied (from 1-4, see text) is given in the last column and may depend on the weights k and k c . The case N = 3 with τ ω is the only one in the table for which r E c ↔ r L may be required, and for which not all k (c) choices are viable. For other notation, see the caption of Table 2.
remark, we note that the argument of Ref. [110] is only valid in the case when the flavour symmetry analysis can be applied directly to the light neutrino mass matrix. In our setup, this corresponds to the situation when light neutrino masses arise either directly from a modular-invariant Weinberg operator, or via a type-I seesaw UV completion such that none of the gauge-singlet neutrinos N c becomes massless in the symmetric limit (so that they can be integrated out). This is the case in the two models considered so far and in the model described in the following section.

S 4 models with τ ω
We finally turn to the most 'structured' cases within the surviving lepton flavour models of Table 4. These arise at level N = 4 in the vicinity of τ = ω and correspond to E c and L being a triplet and the direct sum of three singlets of the finite modular group S 4 , respectively. The expected charged-lepton mass pattern is (m τ , m µ , m e ) ∼ (1, , 2 ).
We have performed a systematic scan restricting ourselves to promising models involving the minimal number of effective parameters (9, including Re τ and Im τ ). Once again, models predicting a massless electron are rejected, while the N c furnish a complete irrep of dimension 2 or 3 (N c are present since Weinberg models require more parameters). Out of 48 models, we have identified the only one which i) is viable in the regime of interest, ii) is not fine-tuned in this regime, and iii) is consistent with the 2σ range for the Dirac CPV phase, predicting δ π while other models lead to δ 0. For this model, L = L 1 ⊕ L 2 ⊕ L 3 with L 1 , L 2 ∼ (1, 2), L 3 ∼ (1 , 2), and E c ∼ (3,4) and 1). The corresponding superpotential reads: Since L 1 and L 2 are indistinguishable, one of the constants α i , with i = 1, . . . , 4, is effectively not an independent parameter and can be set to zero by a suitable rotation without loss of generality. We choose to set α 2 = 0.
At leading order in the small parameter | |, with ≡ 1 − 1+ in the context of this section, the charged-lepton mass matrix reads while the charged-lepton mass ratios are given by withα i ≡ α i /α 1 , i = 3, 4, 5. With respect to charged-lepton mass ratios, the model best-fit point is found to correspond to max(BG) 0.85. Up to an overall normalisation K, the light neutrino mass matrix is instead given by: at leading order in | |, whereg i ≡ g i /g 1 , i = 2, 3. Note that the smallness of | | does not constrain the M ν contribution to the mixing matrix, which depends only on the couplings g i , and large mixing angles are allowed. From the form of M ν it is clear that, in the limit of unbroken SUSY, there is a massless neutrino, even though N c is a triplet. This follows from the modular-symmetric superpotential, which implies the proportionality of the first two columns of Y, reducing its rank and therefore the rank of M ν . The neutrino masses thus read and imply the -independent prediction which, by taking into account the 1σ range for r in Table 3, isolates a viable region in the plane of coupling constants. We show this region in the plane of g 2 /g 3 and g 1 /g 3 in Figure 4. Contours refer to a Barbieri-Giudice measure of fine-tuning for the ratio of neutrino mass-squared differences, given by max{|∂ ln r/∂ lng 2 |, |∂ ln r/∂ lng 3 |}. At the model best-fit point, it has an acceptable value of 2.9. Additionally, the 3σ ranges for g 2,3 are not especially narrow. The result of the fit of this S 4 model is summarised in the last column of Table 5. The viable region in the τ plane corresponds to a neutrino spectrum with NO and is located very close to τ sym = ω, as can be seen from the magnified plot on the right side of Figure 5.
In summary, in the vicinity of the symmetric point, i.e. for small | |, this model can naturally lead to the observed charged-lepton mass hierarchies, see eq. (4.3). The neutrino mass-squared difference ratio r is, in this region, insensitive to and depends only on the two ratiosg 2,3 of neutrino couplings, see eq. (4.6). Furthermore, it is not especially sensitive to these couplings. Finally, since light neutrino masses vanish in the symmetric point, the symmetric limit allows for a generic mixing matrix (case 4 of section 4.1). Therefore, the fit is not expected to be tuned in a way that compensates some 'wrong PMNS' symmetric prediction. In fact, we have numerically verified that sending τ → ω ( → 0) has almost no effect on the values of mixing angles. This can be understood by considering, in turn, each of the contributions to the mixing matrix. The rotation to the mass basis in the neutrino sector, on the one hand, is independent of in the region of interest, see eq. (4.4), and thus has a well-defined limit as → 0 (it is unchanged) even though light neutrinos become massless. This rotation depends only on the ratiosg 2,3 of neutrino couplings. On the other hand, one can check that the charged-lepton rotation arising from the diagonalisation of M e M † e , with M † e given in eq. (4.2), also has a well-defined limit as → 0 even though two of the three charged leptons become massless. This limiting form closely matches the rotation obtained for finite, non-zero , and depends only on the ratiosα 3,4,5 of charged-lepton couplings.

Summary and Conclusions
We have investigated the possibility to obtain fermion mass hierarchies without finetuning in modular-invariant theories of flavour, which do not include flavons. In these theories, hierarchical fermion mass matrices may arise solely due to the proximity of the VEV of the modulus τ to a point of residual symmetry. In our analysis we have considered theories with flavour symmetry described by a finite inhomogeneous or homogeneous modular group, Γ N or Γ N , with N ≤ 5. For N ≤ 5, the finite modular groups Γ N are isomorphic to the permutation groups Γ 2 S 3 , Γ 3 A 4 , Γ 4 S 4 and Γ 5 A 5 , while the groups Γ N are isomorphic to their double covers S 3 ≡ S 3 , A 4 ≡ T , S 4 and A 5 .
In the simplest class of such models considered by us, the VEV of the modulus τ is the only source of flavour symmetry breaking, such that no flavons are needed. Another appealing feature of the proposed framework is that the VEV of τ can also be the only source of CP symmetry breaking in the theory [14]. There is no value of τ which preserves the full modular symmetry. Nevertheless, at certain so-called symmetric points τ = τ sym the modular group is only partially broken, with the unbroken generators giving rise to residual symmetries. There are only three inequivalent symmetric points in the fundamental domain of the modular group [19]: τ sym = i, τ sym = ω ≡ exp(i 2π/3) = − 1/2 + i √ 3/2 (the 'left cusp'), and τ sym = i∞. In these three points, the theories based on Γ N invariance have respectively Z S 2 , Z ST 3 and Z T N residual symmetries. In the case of the double cover groups Γ N , there is an additional Z R 2 symmetry that is unbroken for any value of τ [94], thus enlarging the residual symmetries Z ST 3 and Z T N by the factor Z R 2 , while the Z S 2 symmetry is enlarged to a Z S 4 one. In each of the three symmetric points the standard Z CP 2 symmetry may also be conserved [14]. The indicated residual symmetries play a crucial role in our analysis. In Γ ( ) N modular invariant theories of flavour the fermion mass matrices are modular forms of a given level N . As we show, the mass matrices can be strongly constrained in the vicinity of points of residual symmetries by the properties of the respective modular forms. For each of the three symmetric points, we have developed the formalism which allows to determine the degree of suppression of the elements of the fermion mass matrices, and correspondingly, of their singular values -the fermion masses -in the vicinity of a given symmetric point. More specifically, our analysis showed that, if parameterises the deviation of τ from a given symmetric point, | | 1, the degree of suppression is given by | | l , where l is an integer and can take values i) l = 0, 1, ..., N − 1 in the case of τ sym = i∞, ii) l = 0, 1, 2 if τ sym = ω, and iii) l = 0, 1 when τ sym = i.
These results show, in particular, that it is impossible to obtain the charged-lepton and quark mass hierarchies in the vicinity of the symmetric point τ sym = i as a consequence only of the smallness of | |. As we have proven, the specific value of the power l depends only on how the representations of the fermion fields in the mass term bilinear, denoted for brevity as ψ i and ψ c j , decompose under the considered residual symmetry group. We have derived the decompositions of the weighted irreducible representations of Γ N (N ≤ 5) under the three residual symmetry groups, i.e., the residual decompositions of the irreducible representations (irreps) of Γ 2 S 3 , Γ 3 A 4 = T , Γ 4 S 4 = SL(2, Z 4 ), and Γ 5 A 5 = SL(2, Z 5 ) (they are listed in Tables 6 -9 in Appendix A). The results include also the case of irreps of Γ N , since they represent a subset of the irreps of Γ N .
Having these results we proceeded to identify 3×3 hierarchical fermion mass matrices where the hierarchical pattern is a result of the proximity of the modulus to a point of residual symmetry and no massless fermions are present in the spectrum. We analysed bilinears of the type ψ c i M (τ ) ij ψ j , M being the mass matrix, and considered all possible 3-dimensional representations for the fields ψ i and ψ c j , i, j = 1, 2, 3. While the representations of these fields r and r c are in general reducible, we focused on the case where the same weight is shared between the irreps into which they decompose. The results of this analysis are given in Tables 10 -13 of Appendix B. These tables summarise, for each of the levels N ≤ 5, the patterns of the three fermion masses which may arise in the vicinity of the two potentially viable symmetric points, τ sym = ω and i∞, for all (r, r c ) pairs of 3-dimensional representations and all weights k (c) . We have found that it is only possible to obtain hierarchical spectra for a small list of representation pairs, the most promising of which are collected in Table 2 and correspond to the patterns (1, , 2 ), (1, , 3 ) and (1, , 4 ).
Using the developed formalism and the aforementioned results we have performed a scan searching for phenomenologically viable models of lepton flavour where the chargedlepton mass hierarchies are a consequence of the described mechanism. The chargedlepton masses were obtained from their Yukawa interactions, while neutrino masses are generated either by the Weinberg operator or within a type-I seesaw UV completion. Aiming at minimal and predictive models, we have imposed a generalised CP symmetry on the models, enforcing the reality of coupling constants [14] and restricted the scan to models involving at most 8 constant parameters (including Re τ and Im τ ). Models producing a massless electron were rejected. Out of the many models thus identified, we have selected only those which i) are phenomenologically viable in the regime of interest, and ii) produce a charged-lepton spectrum which is not fine-tuned.
We have found two viable models and both are in the vicinity of τ = i∞. The first is an A 5 model where the lepton doublets L and charged-lepton singlets E c are different triplets of A 5 . For this model, the charged-lepton, Dirac-neutrino and N c Majorana mass terms involve modular forms of weights 4, 5 and 4, respectively. The neutrino masses are generated by the seesaw mechanism with the gauge-singlet neutrino fields N c furnishing a doublet of A 5 . The predicted charged-lepton mass pattern is (m τ , m µ , m e ) ∼ (1, , 4 ). The best description of the input data on the charged-lepton and neutrino masses and mixing (N σ = 0.43) was found to be obtained for Re τ = −0.47, Im τ = 3.11. The second viable model found by us is based on S 4 modular symmetry. In this model the charged-lepton, Dirac-neutrino and N c Majorana mass terms involve modular forms of weights 4, 3 and 2, respectively. The charged-lepton mass pattern is predicted to be (m τ , m µ , m e ) ∼ (1, , 3 ). The viable region in the τ plane is centred around τ = 2.65 i.
Both the A 5 and S 4 viable models were found to require a certain amount of finetuning when describing the neutrino masses and mixing. The presence of fine-tuning in the neutrino sector is explained by the necessity to introduce large corrections to the symmetric-limit PMNS matrix. Addressing the problem of fine-tuning in the neutrino sector, we have found that a modular-symmetric model of lepton flavour with hierarchical charged-lepton masses is expected to be free of fine-tuning 14 if it satisfies at least one of four conditions (see section 4.1). Two of the conditions were formulated earlier in Ref. [110] for arbitrary flavour symmetry groups.
We have constructed a viable model based on S 4 modular symmetry in the vicinity of τ = ω, which is free of fine-tuning in both the charged-lepton and neutrino sectors. It has altogether nine parameters. The neutrino masses are generated via the seesaw mechanism. The charged-lepton mass pattern is predicted to be (m τ , m µ , m e ) ∼ (1, , 2 ). The model predicts, in particular, δ π and m 1 0. We have found also other viable non-fine-tuned S 4 models which, however, predict δ 0.
For the three lepton flavour models constructed, we collect in Table 5 the best-fit values and 3σ allowed ranges of i) Re τ , Im τ and the superpotential parameters, of ii) the charged-lepton masses and neutrino mass and mixing observables used as input in the statistical analysis of the models, and of iii) the predicted lightest neutrino mass, the Dirac and Majorana CPV phases, the sum of the neutrino masses and effective neutrinoless double beta decay Majorana mass.
The results obtained in the present article show, in particular, that the requirement of absence of fine-tuning in both the charged lepton and neutrino sectors in lepton flavour models based on modular invariance is remarkably restrictive. It is hoped that using this requirement it might be possible to identify not more than a few, if not just one, modular-invariant models providing a simultaneous, viable and appealing solution to both the lepton and quark flavour problems.
Note Added. While this work was in its conclusion, Ref. [111] appeared on the arXiv in which the authors investigated the possibility to generate the charged-lepton mass hierarchy in the vicinity of the symmetric point τ sym = i. The two models presented in Ref. [111] are restricted to level N = 3. In these scenarios, the electron is massless by construction, m e = 0, and m e = 0 is expected to arise either through SUSY breaking or from a dim-6 operator. The ratios m e /m τ and m µ /m τ are then associated to independent parameters and not to different powers of the same expansion parameter, as is the case of our work.

A Residual group decompositions
The multiplets of Γ N are 'weighted', i.e. are described by a pair (r, k). 15 At a symmetric point these multiplets decompose into 1-dimensional representations of the corresponding residual symmetry group. In this appendix we present the decompositions of Γ N multiplets (N ≤ 5) under the three residual groups of interest (Tables 6 -9). As seen in section 2.2, these are Z S 4 , Z ST 3 × Z R 2 and Z T N × Z R 2 . Before proceeding, let us comment on the Z R 2 factors in Z ST 3 × Z R 2 and in Z T N × Z R 2 . While kept as part of the residual symmetry group definition in this appendix, they have been omitted in the main text of section 3.2. To understand why they can be ignored without loss of generality, note that a direct product Z n ×Z 2 ≡ a, b |a n = b 2 = 1, ab = ba has 2n irreps 1 ± k , k = 0, . . . , n − 1, which are simply given as products of the Z n and Z 2 irreps: (A.1) In this notation, 1 + 0 is the trivial irrep. The representation under Z 2 is just a sign and does not affect the reality/complexity of a representation. Hence real irreps are 1 + 0 , 1 − 0 and, for even n, 1 + n/2 , 1 − n/2 (one also has (1 ± k ) * = 1 ± n−k ). Since M (τ ) in the bilinear of eq. (3.1) is a function of τ alone, it is R-even. The fields ψ and ψ c are then constrained to carry the same R-parity, i.e. transform with the same sign under Z R 2 . Fields in unhatted representations r -for which ρ r (R) = 1 -are even (odd) under Z R 2 if k is even (odd), while the opposite happens for hatted representations. Keeping this in mind, one can omit the Z R 2 factor and ignore the superscript signs in the following tables. Finally, notice that a Z R 2 factor is hidden in the residual Z S 4 , as S 2 = R. Fields transforming under Z S 4 as 1 0 or 1 2 are R-even while fields transforming as 1 1 or 1 3 are R-odd. Requiring that ψ and ψ c carry the same R-parity implies that one effectively works with Z S 4 /Z R 2 Z 2 , which is why it is generic to considerρ c iρ j = ±1 in section 3.1.2. Table 6: Decompositions of 'weighted' (r, k) multiplets of Γ 2 S 3 under the residual symmetry groups. Irrep subscripts should be understood modulo n, where n = 4, 3 in the first and second columns, respectively. Upper (lower) signs correspond to even (odd) values of k.   Table 9: Decompositions of 'weighted' (r, k) multiplets of Γ 5 A 5 = SL(2, Z 5 ) under the residual symmetry groups. Irrep subscripts should be understood modulo n, where n = 4, 3 in the first and second columns, respectively. Upper (lower) signs correspond to even (odd) values of k.

B Possible hierarchical patterns
In this appendix we list the hierarchical patterns which may arise in the vicinities of the two symmetric points of interest (see main text). We consider in turn the finite modular groups Γ 2 S 3 , Γ 3 A 4 = T , Γ 4 S 4 = SL(2, Z 4 ), and Γ 5 A 5 = SL(2, Z 5 ) (Tables 10 -13). We have focused on 3-dimensional (possibly reducible) representations (r, r c ) entering the bilinear (3.1). Dependence on the weights k (c) may only arise for τ ∼ ω and through the combination K = k + k c , modulo 3.
Note that for N = 2 the residual symmetry group at τ sym = i∞ is Z T 2 . Mass matrix entries are then expected to be either O(1) or O( ) and, as was the case for τ i, one cannot obtain the sought-after hierarchical patterns from the smallness of alone. As such, only τ ω is considered in Table 10. Table 10: Leading-order mass spectra patterns of bilinears ψ c ψ in the vicinity of the symmetric point ω, for 3d multiplets ψ ∼ (r, k) and ψ c ∼ (r c , k c ) of the finite modular group Γ 2 S 3 . Spectra are insensitive to transposition, i.e. to the exchange ψ ↔ ψ c . Congruence relations for k + k c are modulo 3. r r c τ ω k + k c ≡ 0 k + k c ≡ 1 k + k c ≡ 2 Table 13: Leading-order mass spectra patterns of bilinears ψ c ψ in the vicinity of the symmetric points ω and i∞, for 3d multiplets ψ ∼ (r, k) and ψ c ∼ (r c , k c ) of the finite modular group Γ 5 A 5 = SL(2, Z 5 ). Spectra are insensitive to transposition, i.e. to the exchange ψ ↔ ψ c . Congruence relations for k + k c are modulo 3. r r c τ ω τ i∞ k + k c ≡ 0 k + k c ≡ 1 k + k c ≡ 2 3 3