SUSY SU(5) $\times S_{4}$ GUT Flavor Model for Fermion Masses and Mixings with Adjoint, Large $\theta^{PMNS}_{13}$

We propose an $S_{4}$ flavor model based on supersymmetric (SUSY) SU(5) GUT. The first and third generations of \textbf{10} dimensional representations in SU(5) are all assigned to be $1_{1}$ of $S_{4}$. The second generation of \textbf{10} is to be $1_{2}$ of $S_{4}$. Right-handed neutrinos of singlet \textbf{1} and three generations of $\bar{\textbf{5}}$ are all assigned to be $3_{1}$ of $S_{4}$. The VEVs of two sets of flavon fields are allowed a moderate hierarchy, that is $\langle\Phi^{\nu}\rangle \sim \lambda_{c}\langle\Phi^{e}\rangle$. Tri-Bimaximal (TBM) mixing can be produced at both leading order (LO) and next to next to leading order (NNLO) in neutrino sector. All the masses of up-type quarks are obtained at LO. We also get the bottom-tau unification $m_{\tau}=m_{b}$ and the popular Georgi-Jarlskog relation $m_{\mu}=3m_{s}$ as well as a new mass relation $m_{e}=\frac{8}{27}m_{d}$ in which the novel Clebsch-Gordan (CG) factor arises from the adjoint field $H_{24}$. The GUT relation leads to a sizable mixing angle $\theta^{e}_{12} \sim \theta_{c}$ and the correct quark mixing matrix $V_{CKM}$ can also be realised in the model. The resulting CKM-like mixing matrix of charged leptons modifies the vanishing $\theta^{\nu}_{13}$ in TBM mixing to a large $\theta^{PMNS}_{13}\simeq\theta_{c}/\sqrt{2}$, in excellent agreement with experimental results. A Dirac CP violation phase $\phi_{12}\simeq\pm\pi/2$ is required to make the deviation from $\theta^{\nu}_{12}$ small. We also present some phenomenological numerical results predicted by the model.


Introduction
Solar and atmospherical neutrino oscillation experiments have measured leptonic mixing angles with great accuracy. The resulting lepton mixing Pontecorvo-Maki-Nakagawa-Sakata matrix U P M N S [1,2] within the errors of the experimental measured solar and atmospherical mixing angles, can be well compatible with the simple Tri-Bimaximal (TBM) mixing pattern, introduced by Harrison, Perkins and Scott [3]: which predicts that θ ν 23 = π 4 and θ ν 12 = arcsin 1 √ 3 , but θ ν 13 = 0. The two nonzero leptonic mixing angles θ ν 12 and θ ν 23 are predicted to be rather large by contrast the quark mixing angles, which are known to be very small [4]. The pattern suggests some kind of underlying non-Abelian discrete flavor symmetry G f would exist in the lepton sector at least. Indeed the Tri-Bimaximal neutrino mixing is successfully explained by a certain number of flavor models based on different discrete non-Abelian symmetry groups, such as tetrahedral group A 4 , which has only 12 elements, seems to be suitable and economical to reproduce the simple mixing structure as a good first order approximation after symmetry broken properly. An alternative approach is to enlarge the A 4 flavor symmetry group, such as the S 4 group is the most popular one. In most of models either Abelian discrete groups Z N or continuous U (1) groups are needed to eliminate unwanted couplings both in vacuum alignment and fermion masses. For the flavor models based on the typical discrete symmetries, please see Ref. [5][6][7] for a review. The other models based on continuous groups have been proposed [8][9][10][11][12][13][14]. Most of the models based on flavor symmetry groups can obtain a non-zero reactor angle θ 13 ∼ O(λ 2 C ), with λ c = sin θ c ≃ 0.22 being Wolfenstein parameter [15], where θ c is Cabibbo angle, by adding higher order corrections to both charged lepton sector and neutrino sector. The resulting small θ 13 was within the range of global neutrino data fits [16] before the determinate large θ 13 Daya Bay [17] neutrino experiment measured. The deviations from the TBM values of θ 12 and θ 23 are also at most O(λ 2 C ) from both charged lepton and neutrino corrections, in agreement at 3σ error range with the experimental data or say global fits. Even before the Daya Bay result, however, there have emerged direct evidence of large θ 13 while take account of the data from T2K [18], MINOS [19] and Double Chooz [20]. The evidence from the reactor neutrino experiments brings great challenge in flavor model building because subleading corrections are usually too small to generate the sizable θ 13 . Another problem is these models usually work well in lepton sector to produce suitable mass hierarchies and mixing for both charged leptons and neutrinos, but they are difficult to reproduce the correct mass hierarchies and/or mixing of quark sector without fine tuning [21] or accidental enhancement [22].
Besides the models based on Tri-Bimaximal mixing pattern ansatz, other models based on similar simple mixing patterns within the known solar and atmospherical neutrino data were proposed, such as Bimaximal [23][24][25], Golden-Ratio [26][27][28] and Democratic [29] mixing patterns are the most widely discussed in flavor models. Many of them, however, regardless of the specific seesaw type, have in common that vanishing θ ν 13 in lepton sector. Since the reactor neutrino experiments verified the relatively larger reactor angle θ 13 ∼ λ C / √ 2, most of the models based on the TBM ansatz and else would be ruled out, only a small part that explain the experimentally determined value would be maintained. The global fits of solar, reactor and atmospherical neutrino experimental data only have hinted a non-zero but small θ 13 years ago. The Daya Bay Collaboration [17] now has confirmed a relatively larger θ 13 with a significance of 7.7 (the first is 5.2) standard deviations from the reactor ν e → ν e oscillations. The best-fit result in 1σ range is sin 2 2θ 13 = 0.089 ± 0.010(stat) ± 0.005(syst) (2) which is equivalent to θ 13 ≃ 8.7 • ± 0.8 • . And the RENO [30] also reported that sin 2 2θ 13 = 0.113 ± 0.013(stat) ± 0.019(syst) The accurate nonzero reactor angle θ 13 implies Tri-Bimaximal lepton mixing pattern and else would be ruled out. However the opposite is true because there still exists the possibility of maintaining the mixing angles θ ν 12 = arcsin 1 √ 3 and θ ν 23 = π 4 which TBM predicted as the leading order (LO) result of a model. In a general scheme considerations for obtaining the large θ 13 from zero, the most popular and well motivated correction to TBM mixing is the contributions from charged lepton mixing, especially inducing a sizable θ e 12 ∼ θ c is viable in the Grand Unified Theories (GUT) flavor models, see [31,32] as example. In a series of models without GUT the correctional contributions to mixing angles from charged lepton sector are at most O(λ 2 c ), as demonstrated above. Since quarks and leptons can be unified in joint representations of GUT gauge symmetry group, the elements of mass matrices of quarks and leptons can be expressed by some simple single joint operators. In a generical GUT scheme the down-type quark and charged lepton Yukawa matrices are the crucial factor to create GUT relations that connect θ 13 with Cabibbo angle, see [33,34]. The gauge symmetry groups are usually chose as SU (5), SO (10) or Pati-Salam context SU(4) C × SU(2) L × SU(2) R . The simplest GUT gauge symmetry group is SU(5) [35], in which matter fields of standard model are assigned to be5 and 10 dimensional representations. In a family the right-handed down-type quark and left-handed lepton SU(2) doublet are unified in5(d c , ℓ e ) L , and the rest member of the family are collected in 10(q 1 , u c , e c ). We remark that a series of models based on discrete flavor symmetry group together with a GUT gauge group have been proposed, for example the SU(5) ×A 4 [36][37][38], SU(5)×S 4 [39][40][41][42][43] and SU(5)×T ′ [44], SO(10)×A 4 [45,46], SO(10)×S 4 [47][48][49][50], SO(10)×P SL 2 (7) [51,52] and SO(10)×∆ 27 [53].
In this paper we propose a SUSY SU(5) GUT flavor model, with S 4 × Z 4 × Z 6 as flavor symmetry group. The flavor symmetry S 4 can be spontaneously broken by Vacuum Expectation Values (VEV) of flavon fields, they are Φ e = (ϕ, η, χ, ξ, ρ) in charged fermion sector and Φ ν = (φ, ∆, ζ) in neutrino sector. The extra symmetry Z 4 × Z 6 guarantees both the desired vacuum alignment and observed fermion mass hierarchies together with mixings. The fundamental assumption we adopted is similar with Ref. [54], of which the VEVs of Φ e and Φ ν allow a moderate hierarchy: Φ ν ∼ λ c Φ e . The dynamical tricky assumption makes the θ P M N S 13 around O(λ c ) possible. In neutrino sector we introduce SU(5) singlets 1 as right-handed neutrinos which are triplet 3 1 of S 4 , the neutrino masses are simply generated through type-I see-saw mechanism. Tri-Bimaximal mixing pattern can be produced exactly at LO, and even still holds exactly at next to next to leading order (NNLO) in neutrino sector, it is a salient feature of our model. Charged fermion mass hierarchies are controlled by spontaneously broken of the flavor symmetry without introducing Froggatt-Nielsen mechanism [55], especially the mass of up quark is obtained at next-to-leading order (NLO) corrections since m u = 0 at LO in our model. Since the up-type quark masses are more hierarchical than the down-type quark ones, the mixing in downquark is larger than that in up-quark, the resulting θ d 12 ≃ θ c holds to the first order approximation. The correct quark mixing matrix V CKM can also be realised in the model. The same set of joint operators of down-type quarks and charged leptons imply the mass matrix of charged leptons is similar with that of down-type quarks regardless of different group-theoretical Clebsch coefficients and transposed relations. The resulting mixing among charged leptons implies a sizable mixing angle θ e 12 ∼ θ c and another two relative smaller angles θ e 13 ∼ θ 4 c , θ e 23 ∼ θ 4 c . The CKM-like mixing matrix of charged leptons would modify the vanishing θ ν 13 in TBM mixing to a large θ P M N S
The paper is organized as follows. In Section 2 we discuss the basic strategic considerations for obtaining the large θ e 12 ∼ λ c in the model, also we elucidate a simple assumptions about the hierarchy between VEVs of different flavon fields. In Section 3 we introduce all matter fields and flavon contents of our model, giving vacuum alignments of scalar fields. In Section 4 the LO neutrino masses are obtained and exact Tri-Bimaximal mixing matrix are reproduced. In Section 5 we obtain the masses of charged fermions (up-type quarks, down-type quarks and charged leptons) at LO. Also we can get the mixing matrices of quarks and leptons, V CKM and U P M N S , respectively. Section 6 is devoted to subleading corrections to VEVs of all flavons and LO masses of fermion masses and mixings. In Section 7 we show a bit of phenomenology of the model predicted in numerical results. Section 8 is our conclusion.

The strategy and assumptions
In a large classes of flavor models that give arise to TBM and else mixing patterns with or without GUT context, the order of θ 13 is usually about O(λ 2 c ) ∼ 3 • with higher dimensional corrections. In order to modify the Tri-Bimaximal mixing pattern to be compatible with sizable reactor neutrino mixing angle θ P M N S

13
, the most popular and well-motivated correctional approach is provided by the large charged lepton mixing contributions. The GUT flavor models can usually describe quark masses and mixing angles successfully, give arise Cabibbo angle θ q 12 ≃ θ c . The difficult is to generate a relative large mixing between the first and second generation charged leptons (electron and muon), namely that θ e 12 ≃ θ c , since the mixing is negligible small (about O(λ 2 c )) in most flavor models (for simplicity θ c and λ C are treated as the same below). In quark sector the mixing angles among up-type quarks are usually negligible small, the only larger Cabibbo angle of CKM matrix originates mainly from the mixing between first two generations of down-type quarks. In the context of unified theory the Yukawa matrices of charged leptons and down quarks are unified in single joint operators, which provide a possible approach to generating a larger mixing angle θ e 12 ≃ λ c . However traditional unified operators give rising different relations, such as Georgi-Jarlskog (GJ) relation [56], between quark and charged lepton will give θ e 12 ≃ λ c /3. The GJ factor of 3 leads to θ 13 ≃ λ c /3 √ 2 in large classes of GUT flavor models, which now contradicts with experimental data. In order to eliminate the factor of 3, we should consider a general form of Yukawa matrices for both down-type quarks and charged leptons, in which make sure the dominant 1-2 mixing angles are both about O(λ c ). Indeed the basic strategical considerations for achieving θ 13 ≃ λ c / √ 2 in GUT flavor model has been suggested in [33,34,57]. Considering the upper 2 × 2 part of Yukawa matrices and omitting the specific VEVs of Higgs and/or scalar fields, the desired Yukawa matrices of down quarks and charged leptons that lead to θ d,e 12 ≃ λ c are generically of the form where parameters c d,e ij (i, j=1,2) are complex numbers with modulus of order one. The desired mixing angle θ d 12 ≃ |(Y D ) 21 /(Y D ) 22 | ≃ λ c , the Cabibbo angle can be easily derived. Similarly in charged lepton sector the mixing angle θ e 12 ≃ |(Y e ) 21 /(Y e ) 22 | ≃ λ c . More precisely the formalism of Y D and Y e that satisfy Georgi-Jarlskog relation [56] in SU(5) GUT context are as following It is an important goal to obtain such form Y e above in the model. The VEVs alignment directions of flavons in different sectors are actually the key to approach the (12) and (22) relation in Y e . Since the matter fields are assigned to be different dimensional S 4 representations, their production together with flavon combinations will have different contributions to the elements of mass matrices, and finally affected the mixing matrices.
In order to produce the form of Y D,e in Eq. (5), the difficulties are to generate different vacuum alignment and guarantee the hierarchy between (21) and (22) elements. In the model we take VEVs of flavons in Φ e to be of same order of magnitude, i.e., Φ e /Λ ∼ λ 2 c , where Λ is the cutoff scale of the theory. We have assumed a moderate hierarchy between the VEVs of Φ e and Φ ν : Φ ν ∼ λ c Φ e , thus one can also obtain Φ ν /Λ ∼ λ 3 c . The hierarchical relation is allowed in principle, and the VEVs of Φ e and Φ ν are subject to some phenomenological constrains, especially the mass hierarchies of charged fermions and their mixings. We parameterize the VEVs/Λ by small expansion parameters ǫ for Φ e /Λ and δ for Φ ν /Λ, i.e.
The above hierarchy assumption is the fundamental ingredient of our model, since it provides a possibility to generate the desired correlation between (21) and (22) elements in Eq. (5). In fact the mixing angles θ d 12 ∼ λ c and θ e 12 ∼ λ c in the model are produced just by the relative hierarchy since The present of the specific dynamical tricks are in principle allowed by "separated" scalar potential, which is guaranteed by the auxiliary Abelian flavor symmetry group G A . The G A separates the Φ e and Φ ν generically as follows where the V (Φ ν , Φ e )| sub is the subleading scalar potential, at least at NLO, and is usually called "partially" separated. At LO the two sectors are naturally separated, but it is not the case for subleading corrections at NLO and/or NNLO. The scalar potential is "fully" separated while which makes a hierarchy between Φ e and Φ ν possible, see Lin's work in Ref. [54]. In our model the auxiliary flavor group is chosen as However in the present GUT flavor model the "fully" separated scalar potential is not exactly the same as that in Lin's proposal. In fact Φ e should be divided into two parts as well, the set Φ e 1 = (ϕ, η, χ) and the set Φ e 2 = (ξ, ρ). Flavons in Φ e 1 are charged under both Z 4 and Z 6 , while those in Φ e 2 are charged only under Z 4 . Note that the flavons in Φ ν are also charged only under Z 6 . Due to the different constrains on the flavons, the subleading corrections to Φ e 1 arising at NLO, while those to Φ e 2 arising at NNLO, and the subleading corrections to Φ ν are suppressed by 1/Λ 2 at NNLO but do not rely on Φ e . For clarity the subleading scalar potential is separated under G A = Z 4 × Z 6 as following form where the dots stand for higher order subleading terms. The scalar potential of Φ ν is separated at both LO and NNLO, thus the magnitude of Φ ν /Λ is not necessarily the same as that of Φ e /Λ. Actually it is feasible to build a model for TBM based on the S 4 symmetry with the allowed hierarchy in Eq. (6). The hierarchy assumption not only gives arise to the large θ e 12 ∼ λ c , but also produces the correct up quark mass at NLO. The mass hierarchies of up-type quarks are roughly as m u : m c : m t = λ 8 c : λ 4 c : 1, we note that λ 8 c is not only given by ǫ 4 , but also by ǫδ 2 . Indeed the mass of up quark is given by the flavon combinations of order ǫδ 2 at NLO corrections. To a certain extent the hierarchical VEVs are even necessary in the present model. More discussions about Φ e 1 and Φ e 2 : we remind the reader that the diagonal elements of charged fermions (both quarks and charged leptons) mass matrices are mainly determined by the VEVs of Φ e 1 but not Φ e 2 , and the off-diagonal elements (mainly mean (12) and/or (21) elements) should be comprised of Φ ν and Φ e . In the context the desired correlation between (21) and (22) elements of both downtype quarks and charged leptons mass matrices demands a bridge connects Φ e 1 and Φ ν , and the set Φ e 2 is in fact the bridge. The symmetry guarantees the matter fields couple to different flavons, then elements (21) and (22) contain VEVs of different sector, made the ratio (Y D,e ) 21 /(Y D,e ) 22 reduces to the ratio Φ ν / Φ e possible, finally the desired mixing angles θ d,e 12 ∼ δ/ǫ ∼ λ c would emerge.

The Construction of Model: Fields content and vacuum alignment
In this section we introduce the SUSY S 4 flavor model with auxiliary Z 4 and Z 6 symmetries. The flavor symmetry group S 4 is the permutation group of four objects, as well as the invariance group of octahedron and cube. It has 24 elements, which can be generated by two basic permutations S and T as the generators : In group theory the generic permutation can be expressed by (1, 2, 3, 4) → (n 1 , n 2 , n 3 , n 4 ) ≡ (n 1 n 2 n 3 n 4 ). The two basic permutations S = (2341) and T = (2314) is used to generate the elements of S 4 . The group has five inequivalent irreducible representations: two three-dimensional representations 3 1 and 3 2 , one 2-dimensional 2, and two one-dimensional 1 1 and 1 2 representations. The multiplication rules are presented as follows: For the detailed irreducible representation matrices and the matrices of generators S, T, one may see APPENDIX A. In the model matter fields, Higgs and flavon fields are assigned to be different representations of gauge group SU(5). All matter fields in SU(5) are unified into5 and 10 dimensional representations, denoted by F and T 1,2,3 , respectively. The Higgs fields include the SU(5) 5,5, 45 and 45-dimensional representations, right handed neutrinos and all flavon fields are SU(5) gauge singlets 1 in the model. All matter fields, Higgs and flavons are also assigned to be different representations of flavor symmetry group S 4 . The first and second generations of 10 dimensional representation T in SU(5) are assigned to be 2 of S 4 doublet, with the first family as second component of the doublet and the second family as first component of the doublet, respectively, i.e., T = (T 2 , T 1 ) T . This special assignment will avoid the down quark and strange quark masses to be the same order, see footnote in Ref. [42] for details. The third generation of 10 dimensional representation T 3 is assigned to be 1 1 of S 4 . Right-handed neutrinos of singlet 1 and three generations of5 are all assigned to be 3 1 of S 4 . All Higgs are assigned to be 1 1 of S 4 . The flavon fields, which include all different S 4 representations, are introduced to break the S 4 flavor symmetry spontaneously. To be specific, the left-handed down-type quarks in three colors and doublet leptons are collected in SU(5) representation5 as whereas the representation 10 contains SU(2) L doublet quarks as well as up-type quark singlet and charged leptons where i = 1, 2, 3 indicates the fermions family indices of standard model (SM), and R, B, G stand for color indices. The matter fields and flavons in the SU(5)×S 4 ×Z 4 ×Z 6 model, with their transformation properties under the flavor symmetry group, are listed in Table 1. The additional gauge singlets, the driving fields ϕ 0 , χ 0 , σ 0 , ξ 0 , φ 0 and ∆ 0 , are listed in Table 2. We also introduce a global U(1) R continuous symmetry which meant a R-parity discrete subgroup. The driving fields carry +2 U(1) R charge, which made them linearly appear in the superpotential. Matter fields and heavy right-handed neutrinos are charged with +1 U(1) R charge, while all Higgs fields and flavons are uncharged.
At first the vacuum alignment would be discussed in the model building. The problem can be solved by so-called supersymmetric driving field method introduced by Altarelli and Feruglio in Ref. [58]. At LO, the superpotential of driving fields, which is invariant under the flavor symmetry S 4 × Z 4 × Z 6 , is given by The vacuum alignments of all flavons in Φ e and Φ ν are determined by deriving w d with respect to each component of driving fields Φ 0 in SUSY limit. After minimized the derivative equations and solved each unknown component of all flavons, the VEVs structures of the flavons can be obtained.
Usually the solutions are not uniquely determined, we should choose one set by taking into account some constrained conditions. The detailed minimization equations of flavons ϕ, η and χ are given as There are two un-equivalent solutions for the set of equations, one gives another solution of the form v ϕ = (1, 1, 1) T v ϕ , χ = (1, 1, 1) T v χ , η = (1, −1) T v η is forbidden by the last term in Eq. (14). The three flavon fields will mainly determine the diagonal elements of mass matrices of charged fermions up to second order of 1/Λ at LO. Another two flavons ξ and ρ in Φ e also appear in second order of 1/Λ at LO, actually act as a bridge between (21) and (22) elements of mass matrices of down-type quarks and/or charged leptons, which would determine the mixing matrices in that sector. Their vacuum configurations are easily obtained by one minimization equation The equation above will give the following VEVs For neutrino sector the minimization equations for vacuum configuration of φ, ∆ and ζ are given by The solutions for above indeterminate equations are listed by three sets of vacuum alignment as below.
The first one is The second one is and the last one is The first solution (21) is used to produce the Tri-Bimaximal mixing pattern in the following sections. As indicated in section 2, the order of magnitude for VEVs of Φ e and Φ ν with respect to cutoff scale Λ The reason for restrict the magnitude of Φ e /Λ about O(λ 2 c ) is that Φ e should be responsible for the strong mass hierarchies of charged fermions, however it is unnecessary for weak neutrino mass hierarchy. It's natural to require the subleading corrections to Φ e should be smaller than m µ /m τ ∼ O(λ 2 c ), or even more strictly smaller than m e /m τ ∼ O(λ 3 c ). Because of the constrain of the auxiliary symmetry Z 4 × Z 6 , the subleading corrections to superpotential w e1 d are suppressed by 1/Λ and arising from invariant combinations of flavons in Φ ν . The situation for w e2 d is odd since the corrections are suppressed by 1/Λ 2 and depend on both Φ e and Φ ν . The w ν d receives subleading corrections which are suppressed by 1/Λ 2 and still arising from neutrino sector only, see APPENDIX B. We can realise the assumption that Φ ν ∼ λ 3 c Λ won't destroy the stability of TBM prediction in neutrino sector because the Φ ν receives extremely small shift along the same direction in the leading order alignment.

Neutrino
The right-handed neutrinos are SU(5) singlets 1 in the model, neutrino masses are generated through type-I seesaw mechanism only where the m D and M M are Dirac and Majorana mass matrices respectively. The two matrices are derived from the superpotential invariant under the flavor symmetry. Concretely the superpotential in neutrino sector is as follows The first three terms contribute to Dirac masses and the last one is Majorana righted-handed neutrinos mass. The Tri-Bimaximal mixing is reproduced by the vacuum alignments of scalar fields φ, ∆ and ζ, see Eq. (21). After Electroweak and flavor symmetry breaking as these flavon fields and Higgs developing their VEVs, neutrinos will gain masses. For Dirac neutrino mass matrix at LO we have and Majorana mass matrix is The eigenvalues of Majorana mass matrix M M can be diagonalized by unitary transformation where α is an phase parameter. All three right-handed neutrinos are degenerate with mass equal to M, the unitary transformation cannot be solely determined.
Using Eq. (28), the derived light neutrino mass matrix can be diagonalized by the transformation where the light neutrino masses m 1 , m 2 , m 3 are (33) is given by Thus the famous TBM mixing matrix is obtained at LO exactly. The phases ϑ i , i = 1, 2, 3 can be easily obtained and is as following Here we would estimate the scale of M. The mass square differences have been determined to high precision in several global neutrino data fits, such as [59] ∆m 2 sol = 7.50 +0.18 −0.19 × 10 −5 eV 2 , where NH (IH) stands for normal (inverted) hierarchy of mass spectrum, another work of global fit data can be seen in Ref. [60][61][62]. The magnitude of masses m i can be roughly estimated around 10 −2 eV − 10 −1 eV, and v 5 is electroweak sclale ∼ 10 2 GeV, then generally the scale of M will be Besides the operators in Eq. (29), it is worth to consider effective operators, i.e., the higher dimensional Weinberg operators [63] would also contribute to neutrino masses. In this model the effective operators for both Dirac and Majorana mass terms are 1 With the VEVs of Higgs fields, the operators generate the neutrino mass terms as follows One can realize the structure of m W is exactly the same as that of m D in Eq. (30), and it can be diagonalized by TBM mixing matrix where the light effective neutrino masses m W 1 , m W 2 , m W 3 come from the Weinberg operators. Obviously we can compare the relative magnitude between m diag ν and m diag W with the assumption that all couplings, y ν 1,2,3 , and y 1,2 , are of order one, then the ratio could be The contribution to neutrino mass would be larger than the seesaw one if the Weinberg operator has cutoff Λ W ∼ Λ GU T . In order to avoid the problem we will require Λ W ≫ Λ GU T , such that Λ W ∼ Λ P lanck , then the 5-dimensional effective operator can be neglected.

Charged Fermions
Different from neutrino sector, charged fermions have large mass hierarchies, especially in up-type quarks sector, and the mixing is much smaller compared with neutrino sector. In the following the masses and mixings of charged fermions in different sectors are analyzed in detail. The first part is devoted to the up-type quarks, whose mass hierarchies are extremely large while the mixings are small. The second part is the down-type quarks and charged leptons.

Up-type quarks
The masses of up-type quarks are generated by S 4 symmetry breaking in the invariant superpotential at LO. The enormous mass hierarchies between up-type quarks should be guaranteed by VEVs of scalar flavons. The invariant superpotential under the whole symmetry group SU (5) where Then after the scalar fields develop their VEVs in Eq. (15), the mass matrix of up-type quarks can be written as Then mass matrix can be diagonalized by bi-unitary transformation in which the mass eigenvalues are The mass of up quark vanishes at LO, and top quark mass is produced at tree level, the mass hierarchy between charm and top quark are obtained given that the v ϕ , v η of order λ 2 c Λ. After diagonalize the mass matrix (45), one can note that the mixing of up-type quarks only exists between charm quark and top quark, it is an experimental acceptable feature of our model. Thus the form of the resulting mixing matrix U L is rather simple as follows with the mixing angles where S u ij is the sine of mixing angles and C u ij ≃ 1 is assumed. The only mixing of order λ 2 c is between the second and third generations, and another two mixing angles vanish at LO. In section 6.3 up quark mass is generated by corrections from higher dimensional operators and shifted VEVs, and the other two vanishing mixing angles at LO would also be corrected to be nonvanishing values.

Down-type quarks and charged leptons
The LO superpotential giving rise to the masses of down-type quarks and charged leptons is 3 The operators involved O D in Eq. (50) should include all possible independent S 4 contractions. The down-type quarks and charged leptons in the superpotential Eq. (50) would obtain their masses as the scalar sector developing their VEVs. Substituting the VEVs, the mass matrix of down-type quark is and that of charged leptons is equivalent to the transposed M D as follows where the coefficients y d ij (i, j=1, 2, 3) and those with primes are linear combinations of LO coefficients. The Georgi-Jarlskog relation [56] produces the factor of 3 differences between elements of M D and M ℓ , which is induced by Higgs H 45 . The two mass matrices M D and M ℓ can be diagonalized by the similar bi-unitary transformations as in up quark sector and that of charged leptons From the mass expressions Eq. (55) and Eq. (56), one can easily find that bottom quark and tau lepton have the same mass, and the mass of muon is three times that of strange quark The bottom-tau unification and Georgi-Jarlskog relation [56] are produced in the model. The unitary transformation matrices D L and V ℓ L are approximately as The complete quark mixing matrix V CKM is composed of mixing matrices of both up-type quark and down-type quark sector then we can directly get all elements of CKM matrix It is an experimental constrains that V ub and V td are of order λ 3 c , which demands a fine tuning between v 45 and v5: v 45 ∼ λ c v5. Adding this condition, we can easily check that the quark CKM mixing matrix is produced correctly. The Cabibbo angle are determined by the mixing between the first and second family down-type quarks, the parameters y d 21 and y d 22 are of order one and the moderate hierarchy Φ ν ∼ λ c Φ e was assumed.
Similarly the resulting lepton mixing matrix U P M N S is written as The desired CKM-like mixing matrix V ℓ L in Eq. (59) implies a large mixing angle between the first and the second generation of charged leptons, and it will remarkably change the lepton mixing, although the TBM mixing is exactly produced in neutrino sector. The large θ 13 arises from the charged lepton sector, and θ 12 prominent deviates from its TBM value. Thus the three leptonic mixing angles θ P M N S ij would be The relations are compatible with leptonic mixing sum rules [64]. Here we should noted that experimental data shows that θ P M N S 12 is close to TBM value, while the sum rules above manifests a sizable deviation. The reason is that CKM-like V ℓ L brings a larger |U e3 | requires large CP violation. Only considering CP phase around π/2 or −π/2 the mixing angle can be shifted to its appropriate range. The detailed discussions about the problem shall be presented in Section 6.4.

Corrections
The subleading corrections to superpotentials above arise from the higher dimensional operators in powers of 1/Λ expansions, which are suppressed by at least one power of 1/Λ and constrained by flavor symmetry. In this work the modified superpotential of driving fields w d would give rise to the shifts of LO VEVs along different directions in Φ e but the same direction in Φ ν . Fermions' superpotentials w ν , w U and w D also give the correctional mass and mixing matrices by adding the subleading operators.

Correction to Vacuum Alignment
Here we just present the final results of shifted VEVs. The detailed calculation procedure is presented in APPENDIX B. The modified VEVs of all flavons are of the form as follows with v η and v ρ undetermined. We remark that the correctional results of each component of the flavons in Φ e are different. Similarly we have Where v ζ still undetermined. We also remark that in fact the correctional results of each components of the flavons in Φ ν are exactly the same, it is an important feature of the model which makes the Tri-Bimaximal mixing pattern still holds at subleading corrections. Take into account the hierarchical VEVs Φ e /Λ ∼ λ 2 c and Φ ν /Λ ∼ λ 3 c at LO, and the subleading operators linear in driving fields are suppressed by different power of Λ, the shifts of both sectors also have different order of magnitude as follows

Corrections to neutrino
Due to the auxiliary symmetry Z 4 × Z 6 , the subleading corrections to neutrino Majorana mass matrix only appear at next to next to leading order (NNLO). The higher dimensional operators arising from inserting bilinear invariant combinations of φ, ∆ and ζ into the superfield N c N c in all possible ways, the corresponding terms are expressed explicitly as follows: Inserting the LO VEVs of the flavons in Φ ν , the NNLO corrections to Majorana mass matrix will be We can note the form of δM M is still compatible with TBM mixing. Denote the three terms that following the matrices as S 1 , S 2 , S 3 respectively, the eigenvalues of the subleading correctional Majorana masses are easily obtained as follows The possible content that spoils the TBM mixing only arising from Dirac mass terms, whose corresponding subleading superpotential is comprised of the shifted VEVs at LO and higher dimensional operators as following Take into account these subleading contributions to both Dirac and Majorana mass matrices, we can easily rewrite the modified light neutrino masses as follows Note that a, b and c in above expressions are not exactly the same as in Eq. (34), because of the redefinition of LO parameters y ν i by absorbing the corrections to Dirac masses in Eq. (71). However the subleading corrections are too small to change the order of magnitudes, we can treat them unchanged.

Corrections to up-type quarks
There are two sources for the correctional contributions to charged fermions, one is the higherdimensional operators, another is the shifted VEVs at LO. The top quark obtains mass at tree level, it signifies that the LO coupling parameter y t of the T 3 T 3 H 5 can not be changed the order of magnitude by the corrections from the higher dimensional operators, even take into account of the corrections from the operators T 3 T 3 H 45 with certain flavons. For the sake of completeness, however, the contributions to (33) elements from the operators are still listed in the following. The correctional contributions to up-type quarks caused by higher-dimensional operators, which are invariant contractions as following where the operators O U are simply written as We will discuss the contributions to elements of quark mass matrix caused by these operators in detail. Collected all contributions from both higher dimensional operators allowed by flavor symmetry group and modified VEVs, the subleading correctional mass matrix of up-type quarks are as follows The complex coefficients y u,u ′ ,u ′′ ij (i, j = 1, 2, 3) are also linear combinations of subleading coefficients. The subleading mass matrix δM U in Eq. (75) together with LO M U in Eq. (45) lead to up-type quark masses as following eigenvalues Note that the coefficients y u ij (i, j = 2, 3) in m c and m t are redefined ones by absorbed the subleading coefficients into leading ones. We should keep in mind that the values of the coefficients are different from corresponding LO ones, however, whose order of magnitudes are not changed by subleading corrections. It is remarkable that the nonvanishing up quark mass m u is given by the contributions from both subleading operators O U 4 i and VEV shifts at LO, though the m u vanishes at LO. Take into account the hierarchical order of magnitudes of parameters ǫ ∼ λ 2 c and δ ∼ λ 3 c , the up-type quarks obtain their correct mass hierarchies, i.e., m u : m c : m t ≃ λ 8 c : λ 4 c : 1. At last the three mixing angles among up quarks are modified as follows The mixing angles among up-type quarks are all small, especially the vanishing θ u 12 and θ u 13 at LO are corrected to be only of order λ 4 c and λ 6 c , respectively. The largest θ u 23 is still of order λ 2 c , the subleading corrections have no effect on the angle.

Corrections to down-type quarks and charged leptons
The contributions of subleading corrections to mass matrices of down-type quarks and charged leptons still come from both the higher dimensional operators by replacing the flavons combinations O D at LO, and those shifted vacuum alignments of flavons in the LO operators. The higher dimensional operators are written as follows ϕχφ∆, ϕχφζ, ϕχ∆∆, ϕχ∆ζ, ϕχζζ, ηχφφ, ηχφ∆, ηχφζ, ηχ∆∆, ηχ∆ζ, ηχζζ, ξφ 3 , ξφ 2 ∆, ξφ 2 ζ, ξφ∆ 2 , ξφ∆ζ, ξφζζ, ρφ 3 , ρφ 2 ∆, ρφ 2 ζ, ρφ∆ 2 , ρφ∆ζ, ρφζζ} For sake of distinguish with the coefficient parameters at LO, the subleading coefficients are denoted by x d ij . The correctional mass matrix of down-type quarks, denoted by δM D , will be kept the entries up to order λ The contributions to the correctional mass matrix δM D induced by shifted vacua are explained as follows. The operator T 3 F ϕH 5 in LO superpotential (50) leads to corrections of order δ 2 v 5 to all elements in the third column of δM D by plugging the nontrivial shifted vacua of ϕ. The shifts δv ϕ 1,2,3 , δv η 1 and δv χ 1,2,3 in the LO operators T F ϕχH 45 and T F ηχH 45 introduce corrections of order ǫδ 2 v 45 to the entries in the first two columns of δM D . We remark that correctional contributions from the operators T F ξφH 45 and T F ρφH 45 are also neglected because of the relative smaller order compared with the contributions from the former two operators. The operators arise at 1/Λ 3 in superpotential w D are also ignored their correctional contributions.
The subleading corrections above include all possible contributions to all entries of mass matrix δM D up to order λ 8 c . Considering the redefinition of LO parameters with contributions of the subleading terms, the correctional down quark mass matrix are expressed as follows The parameters x d ij (i, j = 1, 2, 3) are linear combinations of all subleading coefficients. Consequently the down-type quark mass eigenvalues will be modified as following The elements of CKM matrix will be modified to the delicate form as below The elements of CKM matrix are mainly determined by LO predictions, only small contributions from subleading corrections are included. One can conclude that the LO results are precisely enough to obtain the applicable quark mixing angles, which are coincide with experimental data. Next it is easy to write down the corrections to mass matrix of charged leptons and the correctional masses of charged leptons the final modified lepton mixing angles are almost unchanged when compared with LO predictions, the expressions are as follows The deviational terms clearly show that the corrections to all three lepton mixing angles originate completely from charged lepton sector, except that the induced negligible coefficient is close to TBM value, the large departure is seemingly unsuitable. Setting the expansion parameters ǫ and δ to be positive real numbers, we note that the real part of the complex number y d 12 /y d 22 contains the module and the phase. Thus the problem could be settled by taking into account of CP violating phase, which is in fact the complex phase, denoted by φ 12 , of the ratio y d 12 /y d 22 . Ignoring the higher order deviations, the deviation of angle θ P M N S 12 from the TBM prediction is approximately expressed as following The correlation has been given in [64] as mentioned before, and similar result with minus sign difference (because of different diagonalization conventions of fermion mass matrix) has been obtained in Ref. [32]. The experimental data within measurement error have shown that θ P M N S 12 is rather close to the value arcsin(1/ √ 3) predicted by TBM mixing. In order to be consistent with the TBM value, the phase φ 12 should be around π/2 or −π/2 so that the sizable departure vanishes, or at least decreases to be of order λ 2 c . The detailed study about the link between CP violation and charged lepton corrections to mixing angles is beyond the scope of the present work, one may refer [65] for example.

Phenomenology
The model we have constructed has only analytic form despite of some parameters within it. In the next step we shall present the phenomenological numerical results of some observables we interested. We require the input oscillation parameters ∆m 2 sol , ∆m 2 atm , sin 2 θ 12 , sin 2 θ 13 , and sin 2 θ 23 to lie in their 3σ intervals which are taken from Ref. [62]. It is easy to express the light neutrino mass spectrum with well determined mass differences and the lightest neutrino mass m lightest = m 1 (m 3 ) in NH (IH) spectrum as follows

Mixing angles
In the numerical analysis all the LO and NLO coefficients in analytic expressions of mixing angles are taken to random complex numbers with their modulus within an interval  Fig. 1(a) and Fig. 1(b) respectively. The horizontal lines show the 3σ (green), 2σ (black) and 1σ (red) boundaries of the mixing angles θ P M N S  Fig. 1(a), and it is obviously that φ 12 ∼ π/2 or −π/2 are favoured so that θ P M N S

Sum of neutrino masses, Neutrinoless double beta decay
First we consider two simple arithmetical relationship between light neutrino masses: the ratio |m 3 /m 2 | and sum of all masses against the lightest neutrino mass. The plot of the ratio |m 3 /m 2 | and sum of light neutrino mass k m k as function of the lightest neutrino mass m lightest , which is m 1 (m 3 ) for NH (IH) mass spectrum, are shown in Fig. 3(a) and Fig. 3(b), respectively. Note that the horizontal lines in Fig. 3(b) represents the cosmological bound at 0.19 eV (black), corresponding to the combined observational data from [66], and the upper bounds 0.23 eV from Planck [67]. The ratio tends to a degenerate mass spectrum in both cases as the value of m lightest → 0.1eV which is disfavoured in the model. The masses sum k m k in the model is predicted too similar for both hierarchical mass spectrums to be distinguished using the cosmological bound on the sum of neutrino masses. The effective Majorana mass |m ee | determines the 0νββ decay amplitude, which is also the (11) element of neutrino mass matrix in flavor basis, and meant a real and diagonal matrix of charged lepton mass. It is defined as follows and the β decay effective mass which could measure the non-zero neutrino masses is In the numerical analysis the mass differences are well-known input parameters, as explained in the beginning of this section, and the analytic form of m lightest (m 1 or m 3 ) in Eq. (72) are the single variant of physical quantities |m ee | and m β . The predicted allowed region of the two effective masses against lightest mass, are shown in Fig. 4(a) and Fig. 4(b). The dashed line in Fig. 4(a) shows the future sensitivity of CUORE [68] experiment at 15 meV. In the model |m ee | is predicted to below the upper limit of |m ee |, which is constrained from the Heidelberg-Moscow experiment [69]. As the m lightest grows to be around 0.1 eV, |m ee | tends to degenerate in both NH and IH mass spectrum.
The β decay effective mass m β is predicted to be below the future sensitivity 0.2eV of the KATRIN [70] experiment, as showed in Fig. 4(b). The vertical line represents the sensitivity.

Conclusion
In this paper, we have proposed a flavor model in the framework of SUSY SU(5) GUT based on S 4 × Z 4 × Z 6 flavor symmetry. In the model the first and second generations of matter fields 10 dimensional representation in SU(5) are assigned to be 2 of S 4 , with first family as second component of the doublet and second family as first component of the doublet, respectively. The third generation of 10 dimensional is to be 1 1 of S 4 . Right-handed neutrinos of singlet 1 in SU(5) and three generations of5 in SU (5) are all assigned to be 3 1 of S 4 . The flavons are divided into charged fermion sector Φ e and neutrino sector Φ ν , whose VEVs are of different orders of magnitude: Φ e /Λ ∼ λ 2 c and Φ ν /Λ ∼ λ 3 c , with λ c being Cabibbo angle. Also the energy scale Λ is below the GUT one.
The three right-handed neutrinos are SU(5) singlets, the light neutrino masses are generated via type-I seesaw mechanism only. The diagonalization of neutrino mass matrix leads to Tri-Bimaximal mixing pattern at LO, both normal and inverted hierarchy mass spectrums are allowed. The subleading corrections to both Majorana and Dirac masses arise from the higher dimensional operators and shifted VEVs lead to the same mass structures as that at LO, give no change to the TBM mixing pattern, it is a salient feature of the present model that TBM mixing still holds exactly at NNLO. The only correction to TBM mixing comes from the mixing of left-handed charged leptons, which is a CKM-like mixing matrix. The resulting θ e 12 ≃ λ c tremendously modified lepton mixing matrix, give rise an sizable θ P M N S 13 ≃ λ c / √ 2, in well compatible with the determinate reactor neutrino experiment result. The solar mixing angle is also corrected with departure from its TBM value of order λ c , only take CP violating phase around ±π/2 into account the large deviation may vanish, or at least decreases to be of order λ 2 c .
The LO masses of charm and top quarks are derived from combinations of the tenplets 10 together with flavons ϕ and η, while the mass of up quark vanishes at LO. When taking into account of the subleading corrections at NLO, nonvanishing up quark mass is obtained, and thus arrive at the suitable mass hierarchies m u : m c : m t = λ 8 c : λ 4 c : 1. The mixing at LO only exists between the last two generations of up quarks, given that θ u 23 ∼ λ 2 c . The mixing angles θ u 12 ∼ λ 4 c and θ u 13 ∼ λ 6 c are achieved at NLO corrections. Note that the mass hierarchies of both quarks (up-and down-type) and charged leptons are generated only via spontaneously flavor symmetry broken rather than the Froggatt-Nielsen mechanism. Combining the mixings from the up-and down-type quark sectors, the quark mixing CKM matrix is correctly generated. The Cabibbo mixing angle between first two generations is guaranteed by the relative order of Φ ν and Φ e , while the mixing between the first and third generations demands a fine tuning v 45 /v5 ∼ λ c .
We also present the subleading corrections to flavon alignment in detail, we find the all VEVs Φ ν actually receive very small shifts along the same directions of the LO alignment even at NNLO corrections, but it is not the case for Φ e . The stable solutions of Φ ν also guarantee the stability of TBM mixing in neutrino sector. In the end we show the phenomenological numerical results predicted by the model. Future long base line neutrino experiment with higher precision is possible to verify or falsify the hints about leptonic CP violating phase predicted in the model. The neutrinoless double beta decay experiment is also a tool for testing the model, which can discriminate between the NH spectrum and the IH one.
The five conjugate classes of S 4 group can be written as The explicit expression of S 4 elements in different representations can be found in Ref. [21]. The character table of S 4 group are shown in Table 3. The multiplication rules between various irreducible representations are shown in Eq. (10), and when taking into account of the generators S and T in different representations, we can straightforwardly obtain the decomposition of the product representations. The product rules of S 4 group, with ψ i , ϕ i as the elements of the first and second representation of the product, respectively, are given in detail as follows ⋆ 1 1 ⊗ r = r ⊗ 1 1 = r with r being any representation (A7) The product rules with two-dimensional representation are as follows: and the product rules with three-dimensional representation are as follows

APPENDIX B: Corrections to vacuum alignment
The subleading corrections to LO vacuum alignment stem from higher dimensional operators which are suppressed by one or more power of 1/Λ. The modified driving superpotential will consist of leading order term w 0 d , which is just Eq. (13), and correctional term δw d , which comes from all invariant operators linear in the driving fields that suppressed by 1/Λ at least one power The correctional term δw d contains all contractional invariant subdominant operators under the flavor symmetry in which the complex coefficients a i , b i , c, k i , s i and t i are all of order one but cannot be specific determined according to the flavor symmetry. Operators O i represent all the subdominant invariant contractional operators under the symmetry group The sub-dominate term δw d induces shifted VEVs of all flavon fields, and we can rewrite the modified vacuum alignment as following where three shifts δv η 2 , δv ρ and δv ζ have been absorbed into the redefinition of the undetermined v η , v ρ and v ζ respectively. With only terms linear in the shift δv retained and ignoring the δv/Λ terms, the new minimization equations are still derived by the zeros of F-terms, i.e. the first derivative of new w d in Eq. (B1) with respect to all driving fields. First the minimization equations for the set Φ e are showed as follows A little discussions about these coefficients: Each terms in A i or B i are all of the same magnitude, but it was not the same case for K. Take into account the hierarchy between Φ e and Φ ν , the terms of coefficient K actually have different hierarchies too. Suppose all k i s are of the same order of magnitude, then only the first four terms (k 2 − k 6 ) determine the LO magnitude of K, while the next four (k 7 − k 10 ) and the last six (k 11 − k 16 ) terms are subleading terms compared with the first four ones, with the relative magnitude λ c and λ 2 c with respect to the leading terms, respectively. The subdominate terms will not affect the magnitude of K, thus we can even remove them safely.
The solutions for Eq. (B10) and Eq. (B11) can be easily obtained as follows The shifted VEVs δv Φ e are all of order λ 6 c Λ from Eq. (B14) and Eq. (B15), which imply the relative order of shifted VEVs with respect the LO VEVs of Φ e are λ 4 c , i.e., δv Φ e /v Φ e ∼ λ 4 c . As illuminated in the end of Section 3, the subleading corrections to Φ e should be smaller than the mass ratio mµ mτ or, more strictly me mτ . The results above have shown the corrections are suitable for the model.
At last for the Φ ν sector we have the minimization equations where coefficients C 1 and C 2 are C 1 = 6(s 4 + s 5 ) + 4(s 11 + s 12 ) + 3(s 14 + s 15 ) The solutions to the Eqs. (B16) are also easily obtained as follows Obviously all the shifts in three and/or two components of all scalar fields in Φ e are different but within the same order of magnitude. All the shifts, however, in three and/or two components of all scalar fields in Φ ν are exactly the same, means the small shifts are in the same direction of LO alignment. The result shows the stability of Φ ν , thus no soft terms are needed to drive the superpotential into desired minimum. The stable solutions of Φ ν also guarantee the stability of TBM mixing in neutrino sector. Take into account the conditions in Eq. (27), it is easy to check the relative order of shifted VEVs with respect to LO VEVs as following