Multiple modular symmetries as the origin of flavour

We develop a general formalism for multiple moduli and their associated modular symmetries. We apply this formalism to an example based on three moduli with finite modular symmetries SA 4 , S B 4 and S C 4 , associated with two right-handed neutrinos and the charged lepton sector, respectively. The symmetry is broken by two bi-triplet scalars to the diagonal S4 subgroup. The low energy effective theory involves the three independent moduli fields τA, τB and τC , which preserve the residual modular subgroups ZA 3 , Z B 2 and Z C 3 , in their respective sectors, leading to trimaximal TM1 lepton mixing, consistent with current data, without flavons. E-mail: ivo.de@udo.edu E-mail: king@soton.ac.uk E-mail: ye-ling.zhou@soton.ac.uk ar X iv :1 90 6. 02 20 8v 1 [ he pph ] 5 J un 2 01 9


Introduction
The discovery of neutrino mass and mixing implies that the Standard Model (SM) must be extended somehow. An elegant possibility remains the original type Ia seesaw mechanism [1][2][3][4][5][6][7] involving right-handed neutrinos, which, when integrated out, yield the Weinberg operators H u H u L i L j , where H u = H is the Higgs doublet of the SM and L i is a lepton doublet of the ith family 4 . The minimal type Ia seesaw mechanism supplements the particle content of the SM by just two right-handed neutrinos (2RHN) [9,10], and this approach will be followed in the present paper. However, to explain the observed approximate tri-bimaximal lepton mixing, one must go beyond the seesaw mechanism and consider a non-Abelian discrete family symmetry [11,12]. For example, S 4 has been used to account for trimaximal TM 1 lepton mixing [13,14], enforced by a residual Z SU 2 symmetry in the neutrino sector, and a residual Z T 3 in the charged lepton sector 5 . However such realistic models typically involve many flavons.
Recently it has been suggested that finite modular symmetry might be the origin of flavour mixing with neutrino masses as modular forms [45], leading to constraints on the Yukawa couplings. This has led to a revival of the idea that modular symmetries are symmetries of the extra dimensional spacetime with Yukawa couplings determined by their modular weights [46]. The finite modular groups Γ 2 S 3 [47,48], Γ 3 A 4 [45,46,[48][49][50][51], Γ 4 S 4 [52,53] and Γ 5 A 5 [54,55] have been considered, in which special Yukawa structures are consequences of the modular forms. Compared with traditional neutrino models of flavour symmetry, only a minimal set of flavon fields (or no flavons at all) need to be introduced in the new framework 7 , making such an approach very attractive.
Within the framework of finite modular symmetry outlined above, only a single modulus field τ is usually considered, corresponding to a single finite modular symmetry Γ N . It has been pointed out that particular modular forms, corresponding to special values of τ , preserve a residual subgroup of the finite modular symmetry Γ N . For example, such residual symmetries are considered in [51] as subgroups of the modular A 4 symmetry. With the help of two moduli with different residual symmetry Z 3 in the charged lepton sector and Z 2 in the neutrino sector, it was shown how trimaximal TM 2 lepton mixing may 4 Recently an alternative type Ib seesaw mechanism has been proposed yielding the new Weinberg operators H uHd L i L j , whereH d is a charge conjugated second Higgs doublet with opposite hypercharge [8]. 5 We adopt the standard presentation of the S 4 generators S, T, U where S 2 = T 3 = U 2 = (ST ) 3 = (SU ) 2 = (T U ) 2 = (ST U ) 4 = I [11]. 6 The geometric connection between the origin of the family symmetry due to modular symmetry and the orbifolding method with two extra dimensions has recently been discussed, e.g., in [41,42]. 7 Extension to the quark flavour mixing is given in [48,50,56].
be realised [51]. Also brief discussion on residual symmetry after modular S 4 symmetry breaking is given in [53]. However, the formalism for having two or more moduli fields (as necessary for such a scheme) has not so far been developed, providing one of the main motivations for the present paper.
In the present paper, we shall extend the formalism of finite modular symmetry to the case of multiple moduli fields τ J (J = 1, . . . M ) associated with the finite modular symmetry Γ 1 N 1 × Γ 2 N 2 × · · · × Γ M N M . As an example, we shall then present the first consistent example of a flavour model of leptons with multiple modular S 4 symmetries interpreted as a family symmetry. The considered model involves three finite modular symmetries S A 4 , S B 4 and S C 4 , associated with two right-handed neutrinos and the charged lepton sector, respectively, broken by two bi-triplet scalars to their diagonal subgroup 8 . The low energy effective theory consists of a single S 4 modular symmetry with three independent modular fields τ A , τ B and τ C , which preserve the residual modular subgroups Z A 3 , Z B 2 and Z C 3 , in their respective sectors 9 , leading to trimaximal TM 1 lepton mixing, consistent with current data, without requiring any flavons.
The remainder of the paper is organised as follows. In section 2 we show how the formalism of finite modular symmetry with a single modulus field can be extended to include multiple moduli and an extended finite modular group. In section 3 we have focussed on the case of the single finite modular S 4 symmetry, and have analysed its stabilisers and resulting remnant symmetries. In section 4 we have proposed a model based on three moduli fields associated with a high energy finite modular group S 3 4 , which is broken to a single diagonal S 4 with three independent moduli fields at low energies, whose stabilisers leads to different remnant symmetry in the different sectors, which may be used to enforce trimaximal TM 1 mixing, leading to good numerical fits to the data, once right-handed neutrino mixing is taken into account. Section 5 concludes the paper.
2 From single to multiple modular symmetries

A single modular symmetry
The modular group Γ acting on the complex modulus τ (Im(τ ) > 0) as linear fractional transformations: 8 The idea of using multiple non-Abelian discrete symmetries with di-multiplet fields to explain flavour mixing was first introduced in the littlest mu-tau seesaw model [57], whereas [58] also features di-multiplet fields of non-Abelian continuous family symmetries. 9 Having a separate residual symmetry associated with each of the two right-handed neutrinos and the charged lepton sector was also assumed in the tridirect CP approach [59,60], although here we do not assume any (generalised) CP symmetry. An extension of modular symmetry to include general CP symmetries was given in [61].
where a, b, c, d are integers and satisfy ad − bc = 1. It is convenient to represent each element of Γ by a two by two matrix 10 . Then, Γ is expressed as This group is isomorphic to the projective special linear group P SL(2, Z) = SL(2, Z)/Z 2 .
The modular group has two generators, S τ and T τ , which satisfy S 2 τ = (S τ T τ ) 3 = 1. They act on the modulus τ and take the following forms respectively. Representing them by two by two matrices, we obtain Γ is a discrete but infinite group. By requiring a, d = 1 (mod N ) and b, c = 0 (mod N ), N = 2, 3, 4, · · · , i.e., where k a , k b , k c and k d are integers, we obtain a subset of Γ which is also an infinite group and is labelled as The quotient group Γ/Γ(N ), labelled as Γ N , is a finite group, also called the finite modular group. The finite modular group Γ N can be also obtained by imposing an additional condition for T τ , T N τ = 1, which can be achieved to identify τ = τ + N in the upper complex plane 11 . For N taking some small number, Γ N is isomorphic to a permutation group, in particular, Γ 2 S 3 , Γ 3 A 4 , Γ 4 S 4 and Γ 5 A 5 [44].
In a theory satisfying the Γ N modular symmetry, any chiral superfield φ i , as a function of τ (but does not need to be modular forms), non-linearly transforms as [39], where −2k i with k i an integer is the modular weight of φ i , I i is the representation of φ i and ρ I i (γ) denotes a unitary representation matrix of γ with γ an element of Γ N .
Considering an N = 1 supersymmetric model in the finite modular symmetry, the action in general takes the form [39,40] 10 Note that it need not be a unitary matrix. 11 Note that once τ = τ + N is imposed, where h is a positive constant. The Kähler potential K can be changed at most by a Kähler transformation under Γ N , and the superpotential W is required to be invariant, i.e., An example of the Kähler potential satisfying the Kähler transformation takes the following form, After τ gets a vacuum expectation value (VEV), the Kähler potential leaves kinetic terms for the scalar components of the supermultiplets φ i and the modulus field as 12 The superpotential W (φ i ; τ ) is in general a function of the modulus τ and superfieds φ i . Under the modular transformation, the superpotential should be invariant under the modular transformation [39]. Expanding the superpotential W (φ i ; τ ) in powers of φ i , we obtain Here, Y I Y represents a collection of coefficients of the relevant couplings. It transforms as a multiplet modular form of weight 2k Y and representation I Y , where k Y = k i 1 + · · · + k in is required to be a non-negative integral. Its representation and weight are required for the invariance of the operator under the Γ N modular transformation.

Multiple modular symmetries
All lepton flavour models based on finite modular symmetries in the literature so far have been limited to the case of a single modulus field. No theoretical approach or model has so far managed to include more than one modulus fields in a self-consistent approach, although the latter case has been briefly mentioned in some references, e.g. [53]. In this subsection, we will discuss how to include multiple moduli fields consistently. 12 The scalar component of φ i may gain a non-zero VEV, and this VEV also contributes to the kinetic term of τ . We ignore such a contribution by assuming v φi √ h.
For any finite modular transformations γ 1 , ..., γ M in Γ 1 N 1 × Γ 2 N 2 × · · · × Γ M N M , the chiral superfield φ i , as a function of τ 1 , ..., τ M , now transforms as where k i,J and I i,J are the weight and representation of φ i in Γ J N J , respectively, and represents the outer product of the representation matrices for ρ I i,1 , ρ I i,2 , ..., ρ I i,M .
For an N = 1 supersymmetric model in a series of modular symmetries, the action is extended to the form where h is a positive constant. The superpotential W is required to be invariant under all modular transformations and that the Kähler potential K can be changed at most by Kähler transformations.
Including multiple modulus fields, the Kähler potential can be written as, where all h J are positive constants. Since each modular symmetry is independent from each other, one modulus field getting a VEV leaves the rest of the Kähler potential still satisfying the other modular symmetries. For example, after τ 1 gets a VEV, the Kähler potential is left with Once all modulus fields get VEVs, the Kähler potential gives rise to kinetic terms for the scalar components of the supermultiplets φ i and the modulus fields as In this example, the scalar component of each modulus field performs as a scalar field of vanishing weight in the remaining modular symmetries.
The superpotential W (φ i ; τ 1 , ..., τ M ) is in general a function of the modulus fields τ 1 to τ M and superfields φ i . Under the modular transformation, the superpotential should be invariant under the modular transformation [39]. Expanding the superpotential W in powers of φ i , we obtain 3 Modular S 4 symmetry and its remnant symmetries In this section, we temporarily return to the case of a single modular symmetry, focussing on the case of a single modular S 4 symmetry and its remnant symmetries, before generalising the results to the case of multiple S 4 symmetries in the next section.

Modular S 4 symmetry
S 4 is a permutation group of four objects. In the framework of modular symmetry, the S 4 modular group is obtained in the series of Γ N by fixing N = 4. In other word, its generators satisfy S 2 τ = (S τ T τ ) 3 = T 4 τ = I. In previous works, it is common to use three generators S, T and U , which satisfy S 2 = T 3 = U 2 = (ST ) 3 = (SU ) 2 = (T U ) 2 = (ST U ) 4 = I [11], to generate S 4 . These traditional generators are related to the modular generators S τ and T τ as which provides a useful dictionary to relate the two types of generators. In the upper complex plane with the requirement τ = τ + 4, S, T and U can be represented by two by two matrices such as Due to the identification in Eq. (5), these representation matrices are not unique. Using Eq. (23), we write out another three elements of S 4 , namely which are order-three elements of S 4 which will appear in our later discussion, Modular forms of even weights in a modular S 4 symmetry can be explicitly constructed in terms of the Dedekind eta function η(τ ) ≡ q 1/24 ∞ n=1 (1 − q n ), with q = e 2πiτ [52]. At lowest weight 2k = 2, there are five independent modular forms. By defining with a 1 + · · · + a 6 = 0, these five independent modular forms can be constructed to be where ω = e 2πi/3 . These five independent modular forms at lowest weight 2k = 2 form a doublet 2 and a triplet 3 of S 4 , Modular forms with higher even weights (2k = 4, 6, · · · ) can be constructed from these five modular forms. In general, the dimension of the linear space formed by the modular forms of weight 2k and level 4 is 4k + 1 [45]. Namely, the nine independent modular forms of weight 2k = 4, which form one 1, one 2, one 3 and one 3 . Among them, the two triplet modular forms are given by At weight 2k = 6, there are 13 independent forms. They form one 1, one 1 , one 2, one 3 and two 3 s of S 4 . Here we only interested in the two 3 s of S 4 . They are given by These modular forms will be used for our model building in the next section. For modular forms with weights up to 10, a full list can be found in [52].
Extension from a single S 4 modular symmetry to a series of modular S 4 symmetries is straightforwardly achieved by following the procedure in section 2.2 with all levels fixed at N J = 4. In each S J 4 , we denote their generators S, T and U by S J , T J and U J , where the subscript is only used to distinguish groups. Modular forms with weights k Y,1 , ..., k Y,M are multiplets of multiple moduli, namely of of τ 1 , ..., τ M .

Stabilisers and residual symmetries of modular S 4
Although a brief discussion on residual symmetry after modular S 4 symmetry breaking has been given in [53], we note that the essential correlation between the modular field and its residual symmetries has not been discussed. In this section, we will give a thorough analysis of this case, uncovering some new results along the way.
We begin by introducing and reviewing the notion of stabilisers of the symmetry which will play a crucial role in residual symmetries. Given an element γ in the modular group S 4 Γ 4 , a stabiliser of γ corresponds to a fixed point τ γ in the upper complex plane which satisfies γτ γ = τ γ . Once the modular field τ gains a VEV at such a stabiliser, τ = τ γ , an Abelian residual modular symmetry generated by γ is preserved. It is obvious that acting γ on a modular form at its stabiliser leaves the modular form invariant, i.e., Following the standard transformation property in Eq. (13), we obtain This equation lead us to the following important properties for the stabiliser and the modular form: • The stabiliser τ γ satisfies |cτ γ +d| = 1 since (cτ γ +d) −2k is an eigenvalue of a unitary matrix.
A special case is that when (cτ γ + d) −2k = 1 is satisfied, ρ I (γ)Y I (τ γ ) = Y I (τ γ ), and we recover the residual flavour symmetry generated by γ. In general, the eigenvalue does not need to be fixed at 1 in the framework of modular symmetry.
In the follow-up of this subsection, we will consider the following stabilisers, Although τ T and τ T S have been discussed in [52] (identified with τ L and τ R therein, respectively), τ S , τ U and τ ST as stabilisers in the S 4 modular symmetry are discussed here for the first time. Here we apply this notation to take the advantage of modular residual symmetries generated by S, T , U , T S, ST and ST S, respectively. Following Eq. (23), it is straightforward to check that these stabilisers are invariant under the corresponding modular transformations respectively, i.e., It is worthy noting that these stabilisers are some typical examples but not the full list of stabilisers of S 4 .
At the stabiliser, the multiplets formed by the modular form may specify interesting directions. We will discuss how the triplet modular forms Y (for k = 1, 2, 3) gain these directions based on the symmetry argument in Eq. (31).
We begin our discussion from modular forms at the stabiliser γ S . We know that Y (2k) for any weight 2k. Given the well-known representation matrix for S in 3 or 3 , Three eigenvalues are given by 1, −1 and −1. The eigenvector corresponding to the eigenvalue 1 is always fixed at (1, 1, 1) T up to an overall factor. Therefore, we conclude that Y (2k) Here, y 3 ( ) ,S is a overall factor determined by the weight and representation. By taking τ S = i∞ into the exact modular form Y i (τ ) in Eq. (26), we obtain q = 0 and 3,S = 0 and y (4) 3 ,S = −3π 2 /16. For 2k = 6, y 3,S = 0, y At the stabiliser τ S , since the eigenvalue is always fixed at 1 regardless of the weight, the residual Z S 2 modular symmetry is identical to the residual Z S 2 flavour symmetry.
We perform a similar discussion for modular forms at the stabiliser τ T . Eq. (31) is simplified into Thus, the selected eigenvector corresponds to the eigenvalue ω 2k , which is weight-dependent.
In the T -diagonal basis we use in the paper, representation matrix for T is given by The triplet form, as an eigenvalue of T , takes a very simple form where the overall factors are also determined by the weight and representation. These results can be checked numerically by taking τ T into exact formulas of modular forms. It is straightforward to obtain Y 1 (τ T ) = Y 3 (τ T ) = Y 5 (τ T ) = 0, and we are left with only two non-zero modular forms, Y 2 (τ T ) = 2.11219i and Y 4 (τ T ) = −2.43895i. Taking them to Eqs. (27), (28) and (29), we arrive at the same above result with y 3,T = −y In the framework of flavour symmetry, the residual symmetry generated by U is usually called µ-τ symmetry. We discuss the modular form at the stabiliser of U . Y (2k) ρ 3 (U ) has one eigenvalue −1 and the other two degenerate eigenvalues +1. The eigenvector with respective eigenvalue −1 is fixed at (0, 1, −1) T without considering an overall factor. The eigenvector with respective eigenvalue +1 is in principle a linear combination of two independent vectors (2, −1, −1) T and (1, 1, 1) T . For odd and even k, we can express Y The coefficients are determined by the weight. Numerically, 3,U = −i2Y 1 (τ U )a for 2k = 4, and y 3,U = 3 √ 2Y 2 1 (τ U )a for 2k = 6. In 3 representations, ρ 3 (U ) has one eigenvalue +1 and the other two degenerate eigenvalues −1. For an even k the direction of Y for an even k .
In addition, we would like to consider stabilisers for the elements T S, ST and ST S. These elements are order-3 elements and stabiliser for each element preserves a Z 3 symmetry. The representation matrices of T S, ST and ST S take the forms 3 ( ) (τ T S ) corresponds to the eigenvalue (1 − τ T S ) −2k = ω k . Thus, we directly arrive at Taking the explicit formulas of modular forms into account, we obtain the overall factors to be y 3 2 ,T S = 0. We turn to the modular forms at stabilisers τ ST . Y  (45) where y  (46) where y We summarise directions of triplet (3 and 3 ) modular forms for lower weights (2k = 2, 4, 6) at stabilisers (τ = τ S , τ U , τ T , τ T S , τ ST , τ ST S ) in Table 1.
All the above discussion in this subsection is based on a single modular S 4 with a single modulus field. Extending to the case of multiple modular symmetries may allow the theory to have several different residual modular symmetries. Namely, the different moduli fields may take different values at different stabilisers. In the next section, we will apply this property to model building.
3 . Here, we have ignored the overall factor if it is non-zero. 0 represents a vanishing modular form, namely, the one which has a zero overall factor. Combining the results of the previous two sections, we see that the extension from one single modular field to multiple moduli fields, as discussed in section 2, opens a window into a new type of modular model building, in which several moduli fields can appear, with each one having a different modular form with a different residual symmetry, of the kind discussed in section 3.

A modular S 3 4 model
As a concrete example, we will show how the results of the previous sections can lead to a consistent model of trimaximal TM 1 mixing, analogous to the traditional approach [13,14]. At high energies, the model in Table 2 is based on three modular symmetries, S A 4 , S B 4 and S C 4 , with moduli fields labelled by τ A , τ B and τ C , respectively. After the moduli fields gain different VEVs, different textures of mass matrices are realised in charged lepton and neutrino sectors.
The transformation properties of the leptons are given in Table 2. We arrange that each lepton has no more than one non-vanishing modular weight in either S A 4 , S B 4 or S C 4 . We note that: 1) The lepton doublets L form a triplet of S C 4 with zero weight; 2) the right-handed leptons e c , µ c and τ c are singlets of S C 4 but have different weights 2k C = −6, −4, −2, respectively; 3) We introduce only two right-handed neutrinos N c A and N c B , which are all singlets but have weights 2k A = −6 and 2k B = −4 in S A 4 and S B 4 , respectively. It is in principle possible to arrange one field with non-vanishing weights in more than one modular symmetry, so our choice here is just for simplicity.
In addition, we introduce two scalars Φ AC and Φ BC . These scalars are assumed to be bi-triplets in the flavour space, arranged in S A 4 × S B 4 × S C 4 as Φ AC ∼ (3, 1, 3) and Φ BC ∼ (1, 3, 3) with zero weights. As bi-triplets, they transform as for any elements γ A , γ B and γ C of S A 4 , S B 4 and S C 4 , respectively. These scalars are introduced to connect three S 4 's together as shown in the superpotential below, where the leptonic superpotential includes the terms responsible for generating lepton masses.
To be invariant under the modular transformation, Y e,µ,τ are 3-plet modular forms in the modular space S C 4 with weight 2k C = 2, 4, 6, respectively, Y A and Y B are 3-plet modular forms in the modular space S A 4 , S B 4 with weights 2k A = 6 and 2k B = 4, respectively. A term, e.g., LΦ AC Y A (τ A )N c A is explicitly written as where m A , m B,1 , m B,2 and m AB are complex free parameters with a mass dimension.

Symmetry breaking of S 3 4 to the diagonal S 4 subgroup
The modular symmetries are broken after the bi-triplet scalars τ A , τ B and τ C gain VEVs. Unlike the flavons introduced in most flavour models in the literature, the VEVs of these scalars are not responsible for special Yukawa textures for leptons, but rather their purpose is to break three modular S 4 's to a single modular S 4 symmetry, identified as the diagonal subgroup and denoted as S D 4 , as depicted in Fig. 1.
The VEVs of Φ AC and Φ BC take the following forms Here again, P 23 represents the (2,3) row/column-switching transformation matrix, and α = 1, 2, 3 corresponds the entries of the triplet of S C 4 , while i = 1, 2, 3 (m = 1, 2, 3) corresponds to those of S A 4 (S B 4 ). These VEV structures are not arbitrarily assumed, but can be simply achieved following the standard driving field method. They are essentially related to the group structure of S 4 and its explicit form is basis-dependent 13 . For details of how to derive them without loss of generality, we refer the reader to Appendix B.
Although S A 4 , S B 4 and S C 4 are broken by these VEVs, the diagonal subgroup S D 4 survives below the symmetry breaking scale, corresponding to the associated transformation γ A = γ B = γ C . In more detail, the S D 4 survives since, given any γ A of S A 4 , there always exists an element γ C of S C 4 which is identical to γ A , and the VEV of Φ AC is invariant under this "contravariant" transformation. Furthermore, there also exists an element γ B of S B 4 which is identical to γ C , and the VEV of Φ BC is also invariant under the transformation. Thus, the modular S D 4 symmetry corresponds to a universal transformation.

The effective low energy theory with modular S 4 symmetry
The effective low energy superpotential, below the S 3 4 breaking scale, involves only a single surviving modular S 4 symmetry, and may be written as, where terms such as e.g., LY A (τ A )N c A may be explicitly written as which is straightforwardly obtained from Eq. (49). This superpotential involves only the single residual S D 4 , and three modular fields τ A , τ B and τ C at the same time.
The above superpotential may be taken as a starting point for models based on a single modular S 4 symmetry, where the three moduli fields introduced in an ad hoc way and taken to be independent fields. However, we have shown that such a model can consistently arise from a high energy model involving three modular groups S 3 4 . The key point of such a model is that, in the low energy effective theory, the three moduli transform under the same S D 4 , i.e., for any γ D ∈ S D 4 , τ A , τ B and τ C transform in the following way, for J = A, B, C. We also write out transformation properties of leptons and those for modular forms where α = e, µ, τ and k e,µ,τ = 3, 2, 1.
To summarise, we have derived a low energy effective flavon-less leptonic flavour model with one modular S 4 symmetry and three independent moduli fields. The importance of this for model building is that, as we shall see shortly, by making use of the different moduli fields, we can access different sets of triplet modular forms, corresponding to having different residual symmetries in different sectors of the theory. This is similar to the traditional approach to model building based on S 4 , but of course is achieved now without having to introduce flavons with certain vacuum alignments.

Flavour structure in the charged lepton sector
In the charged lepton sector, only S C 4 plays a role. We assume the VEV of τ C fixed at τ C = τ T = ω. Following Eq. (39), we obtain for weights 2k C = 6, 4, 2, respectively. These modular forms will lead to diagonal Yukawa couplings for the charged leptons, where all lepton mixing arises from the neutrino sector. Although the diagonal Yukawa couplings are independent, we do not gain any understanding of the charged lepton mass hierarchy in this model.

Flavour structure in the neutrino sector
In the neutrino sector, by selecting τ A = τ T S = 1 2 + i √ 3 2 and τ B = τ U = 1 2 + i 2 we have residual modular symmetries Z T S 3 and Z U 2 , respectively. Following the discussion in section 3.2, we obtain the modular form for the Yukawa coupling whereĉ R ≡ cos θ R e α 1 andŝ R ≡ sin θ R e iα 2 . The Dirac mass matrix M D in the basis where charged lepton and right-handed neutrino mass matrices are diagonal is obtained through V acting on the right of M D , which mixes the columns: Applying seesaw formula, we obtain where µ 1 and µ 2 are real inspect of an overall phase. There are five physical parameters µ 1 , µ 2 , θ R , α 1 and α 2 .
The PMNS matrix is obtained by diagonalising the neutrino mass matrix, U T M ν U = diag{0, m 2 , m 3 }. Since both Y A and Y B are orthogonal to (2, −1, −1) T , we directly arrive at the TM 1 form of lepton mixing matrix [62][63][64][65], TM 1 lepton mixing implies three equivalent relations: leading to a prediction θ 12 ≈ 34 • , in excellent agreement with current global fits, assuming θ 13 ≈ 8.5 • . By contrast, the corresponding TM 2 relations imply θ 12 ≈ 36 • [64], which is on the edge of the three sigma region, and hence disfavoured by current data. TM 1 mixing also leads to an exact sum rule relation relation for cos δ in terms of the other lepton mixing angles [64],

Numerical Fit
As described in previous subsections, we obtain through the use of modular symmetries a flavon-less effective theory which fulfils TM 1 lepton mixing. In this section, we make use of the above analytical sum rules for TM 1 lepton mixing as well as the diagonalisation of the 2 × 2 symmetric matrices which result from the rotation of the neutrino mass matrix by the TB mixing matrix, following the analytic methods presented in [66]. We are thus able to express each observable (the 3 mixing angles, the squared mass ratio and the CPviolating phase δ) in terms of the model parameters ({x}) = ({α 1 , α 2 , θ R , µ 1 , µ 2 }), i.e. the phases α 1 and α 2 , the angle parametrizing the rotation originating from RH neutrino sector, and the parameters governing the contribution from Y A and Y B , µ 1 and µ 2 . These formulas are somewhat complicated and not particularly illustrative, but enable us to easily run a numerical minimisation procedure on a χ 2 function: where P i are the model predictions, BF i the current best-fit values, and the errors σ i correspond here to the average of the 1σ ranges for each observable. We use the bestfit values and 1σ ranges from NuFit 4.0 [67,68]  parameters and the observables tested are the 3 PMNS mixing angles, the phase δ, and the absolute masses obtained from the square roots of the squared mass differences (taking into account that we have only 2 RH neutrinos and normal mass ordering, m 1 = 0).
The obtained best-fit point (BF) corresponds to a χ 2 = 0.74, with the model parameters shown in Table 3, together with the respective predictions for the observables, including mixing parameters, neutrino masses, and the effective neutrino mass parameter in neutrino-less double beta decay m ee = |µ 1ĉ 2 R + µ 2ŝ * 2 R |. These observables (predicted by the analytical formulas for the specific point in parameter space) completely match with the values obtained by performing an entirely numerical diagonalisation for the same point in parameter space. For the best-fit point the observables are all within the 1σ range except δ = 290 • .
For comparison we present also two other benchmark points. In Benchmark 1 (B1), the observables are all within the 1σ range except δ, which is slightly smaller (285 • ) than in the best fit point. Conversely, θ 23 deviates slightly from its best-fit point. The total χ 2 = 1.6 is slightly worse. In Benchmark 2 (B2), δ = 254 • is within the 1σ range. Conversely, θ 13 is slightly deviated from its best fit point and θ 23 = 41.5 • is strongly deviated from its best-fit point and is indeed outside the 1σ range. The total χ 2 = 55 is much worse, although we note that this value is somewhat spurious, given that the expression in Eq. (67) is based on Gaussian distributions, which is not the case for θ 23 , which for B2 contributes 0.99 of the total χ 2 . We are taking the best-fit point for θ 23 = 49.6 • from NuFit 4.0 [67,68].
It is worth emphasizing that these predictions originate from the special directions Y A and Y B obtained from the fixed points in the respective modular symmetries. The best-fit point observables all lie within the 1σ range except δ, which nevertheless lies within its 3σ range and takes a value close to maximal (290 • ).

Conclusion
In this paper we have considered, for the first time, leptonic flavour models based on multiple moduli fields with an extended finite modular symmetry. We reviewed the case of a single modular symmetry Γ with a single modulus field τ and N = 1 supersymmetry, then extended the formalism to include a series of modular groups Γ 1 , Γ 2 , ..., Γ M , where the modulus field for each modular symmetry Γ J is denoted as τ J , where J = 1, ..., M , resulting in the finite modular symmetry Γ 1 We then returned to the case of a single modular symmetry, focussing on the case of modular S 4 symmetry and its remnant symmetries, exploring relations of stabilisers of modular transformations, residual symmetries and modular forms in the framework of finite modular symmetry. In the case of modular S 4 symmetry, several new stabilisers of residual symmetries were identified, where each stabiliser preserves a Z 2 or Z 3 residual symmetry. We discovered a strong correlation between the modular transformation and the modular form at its stabiliser, namely that a modular form at a stabiliser of any modular transformation is an eigenvector of the representation matrix of the modular transformation. Based on this correlation, we were able to determine some new types of modular forms without knowing exact expressions for those modular forms.
As an application of the preceding results, we constructed a flavour model of leptons involving two right-handed neutrinos and three finite modular symmetries S A 4 × S B 4 × S C 4 . Here, S A 4 and S B 4 are modular symmetries for two right-handed neutrinos, respectively, while S C 4 is the modular symmetry in charged lepton sector. They are connected by two bi-triplet scalars. After they gain VEVs, three S 4 's are broken to a single S D 4 , i.e., Independent fixed points in the extra dimensions associated with S A 4 and S B 4 specify (flavon-less) special directions that preserve subgroups of the respective symmetries, whereas a scalar transforming as a triplet of both S A 4 and S C 4 and another scalar transforming one of both as a triplet of both S B 4 and S C 4 acquire vacuum expectation values that break S A 4 × S B 4 × S C 4 to its diagonal subgroup S D 4 . We emphasise that these scalars do not carry any information about flavour.
After the three S 4 's are broken, we arrive at an effective low energy flavour mixing model with a single S 4 modular symmetry but three independent modular fields τ A , τ B and τ C . The independence of these modular fields allows us to assign different VEVs for them which determine the flavour structure. We fix the VEV of τ C at a stabiliser which satisfies a modular Z C 3 symmetry. A diagonal charged lepton mass matrix is obtained. VEVs of τ A and τ B are fixed at other two stabilisers which preserve a different Z A 3 symmetry and a Z B 2 symmetry respectively. The residual modular symmetries justify the special directions that lead to TM 1 mixing. This is similar to the traditional approach to model building based on S 4 , but of course is achieved now without having to introduce flavons with certain vacuum alignments.
Finally, we performed an analysis of the predictions of the model taking into account the existence of RH neutrino mixing (in the model-building basis). When this is taken into account, the 5 observables depend on 4 real model parameters, and we obtain an excellent fit to experiment, with all 3 mixing angles and the squared mass ratio within 1σ of their experimental values and a near-maximal value for δ = 290 degrees. Having two right-handed neutrinos, the model predicts the absolute neutrino mass scale m 1 = 0.
In conclusion, we have developed a general formalism for multiple modular symmetries, analysed the residual symmetries of modular S 4 symmetry, and proposed a realistic model based on modular S 3 4 symmetry, which yields the successful trimaximal TM 1 lepton mixing, without requiring any flavons. 1 The generators of S 4 in different irreducible representations are listed in Table 4, in the basis we used in the main text. This basis is widely used in the literature since the charged lepton mass matrix invariant under T is diagonal in this basis. The following basis-dependent property is satisfied, ρ T 3 ( ) (γ)P 23 = P 23 ρ 3 ( ) (γ −1 ) for any γ of S 4 . The products of two 3 dimensional irreducible representations a and b can be expressed as where The products of two doublets a = (a 1 , a 2 ) T and b = (b 1 , b 2 ) T are divided into  In Appendix B, we apply another basis to calculate the vacuum alignment. Both bases are widely used and have no physical difference. However, we apply this basis because it is simpler to carry out the respective calculations. Representation matrices for generators in this basis are listed in Table 5. Representation matrix for any γ of S 4 satisfiesρ T 3 ( ) (γ) = ρ 3 ( ) (γ −1 ). Basis transformation between the first and second basis are given bỹ for any γ ∈ S 4 , where Irreducible products of two triplet representationsã andb are simply written as where 1 i , 3 i and 3 j are given the same as in Eq. (70).
The full solution for the above equation is not hard to obtained. There are 24 solutions in total. It is convenient to write them as 3 × 3 unitary matrices, The matrices in Eq. (78) are identical to the 3-plet representation matrices of the 24 elements of S 4 in basis in Table 5. Using the basis transformation, we obtain solutions for Φ AC iα = j,β (U † ω ) ij (U † ω ) αβ Φ AC jβ in the basis used in the main text, i.e., that in Table 4. With the help of U † ω = U ω P 23 , we can express these solutions as Φ AC = ρ 3 (γ)P 23 v AC with γ any element of S 4 .
In the main text, we achieved the breaking of modular symmetry S A 4 × S C 4 to the flavour symmetry S D 4 by assuming the VEV Φ AC iα = (P 23 ) iα v AC (corresponding to the first solution in Eq. (78), Φ AC iα ∝ δ iα ). We comment in the following that this VEV is not special and any VEV with the form ρ 3 (γ)P 23 v AC can lead to the breaking of two S 4 's to a single S 4 .
To see this feature more clearly, let us pick up the second solution in Eq. (78) as an example. This solution is represented as a 3 × 3 matrix to be Φ AC = ρ 3 (S)P 23 v AC , i.e., corresponding to the element S. The operator LΦ AC Y A N c A H u , for instance, after Φ AC gains the VEV ρ 3 (S)P 23 v AC , is effectively expressed as For any element γ D of S D 4 , which leads to (L 1 , L 2 , L 3 ) → (L 1 , L 2 , L 3 )ρ T 3 (γ D ), with the help of basis-dependent property ρ T 3 (γ D )P 23 = P 23 ρ 3 (γ −1 D ), one can always require the transformation properties for Y A and N c A in the following, such that the effective operator L Φ AC Y A N A H u is invariant in S D 4 . The transformation property for Y A is equivalent to the following transformation property for [ρ 3 (S)Y A ], Similarly, S B 4 × S C 4 can be broken to S D 4 by the VEV Φ BC mα = (P 23 ) mα v AC . Eventually, we realise the breaking of three S 4 's to a single S 4 .