Confronting the prediction of leptonic Dirac CP-violating phase with experiments

We update and improve past efforts to predict the leptonic Dirac CP-violating phase with models that predict perturbatively modified tribimaximal or bimaximal mixing. Simple perturbations are applied to both mixing patterns in the form of rotations between two sectors. By translating these perturbed mixing matrices to the standard parameterization for the neutrino mixing matrix we derive relations between the Dirac CP-phase and the oscillation angles. We use these relations together with current experimental results to constrain the allowed range for the CP-phase and determine its probability density. Furthermore, we elaborate on the prospects for future experiments probing on the perturbations considered in this work. We present a model with $A_4$ modular symmetry that is consistent with one of the described perturbed scenarios and successfully predicts current oscillation parameter data.


I. INTRODUCTION
With the discovery of the Higgs scalar by the LHC in 2012 the standard model (SM) of particle physics took the seat as the most predictive high energy theory so far. While further experimental efforts keep giving results that are mostly consistent with the SM, one of its sectors has, since long ago, given the best motivation for physics beyond the SM: Neutrinos. First proposed as a way to fix conservation laws in beta decays, they have had an eventful history, while they went from massless to having tiny masses and changing flavoroscillate-while travelling due to mixing between flavor states. The first evidence of neutrino oscillations was reported in 1998 [1]. Neutrino oscillations were firmly established in 2001 using solar neutrinos [2], and, since then, experiments regularly close in on their oscillation pattern and the mass differences responsible of these oscillations. Fast forwarding to 2011 and 2012, the Double Chooz [3] and Daya Bay [4] experiments measured θ 13 = 0 with enough precision to open the possibility of a Dirac type CP-violating phase in the mixing of the leptonic sector of the SM described through the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix [5,6].
The usual approach to extend the SM to include neutrino masses and mixing employs a discrete flavor symmetry at a very high energy. After the spontaneous breaking of this symmetry at lower energies, residual symmetries remain in the charged and neutral leptons mass masses, thus, resulting in particular mixing patterns in the PMNS matrix, U PMNS . Before the measurement of the reactor angle, θ 13 = 0, models that predicted no mixing between first and third family were popular, in particular models that predicted two maximal oscillations popularly known as bimaximal (BM) mixing [7][8][9][10][11][12][13] and, as more data accumulated, other works appeared suggesting maximal mixing of two and three families, known as tribimaximal (TBM) mixing [14][15][16][17]. Naturally, after the measurement of a non-zero reactor angle, the exact BM and TBM mixing patterns were ruled out. In more complicated formulations, these patterns can be considered the result of residual symmetries that need to be broken by perturbations that permit the appearance of a non-zero reactor angle. Interestingly, this type of formulations often result in relationships between oscillation parameters that allow an estimation of the size of the Dirac type leptonic CP-violating phase. This is the idea that was developed in Ref. [18] as well as in other several works [19][20][21][22][23][24][25][26][27][28][29][30][31][32]. In the present work we attempt to follow up on the scenarios explored in Ref. [18] and extend the analysis to probability densities for the CP-violating phase based on currently allowed ranges for oscillation parameters from experiments. Moreover, we simulate the effects of the constraints from these scenarios to estimate their chances of survival in three long-baseline experiments that may be operative in the near future. We complete this work by showing how one of these scenarios can be realized in a flavor model of neutrino masses and mixing.
The rest of the paper is laid out as follows: In Sec. II we introduce the perturbations to TBM mixing that will be used along the rest of this work and their constraints on oscillation parameters. In Sec. III we present probability densities related to the CP-violating phase considering constraints from the cases of Sec. II. In Sec. IV we describe the perturbed scenarios applied to BM mixing and comment on the effects of current experimental constraints. In Sec. V we present the prospects of future experiments expected to constrain the scenarios considered here. In Sec. VI we construct a model using A 4 modular symmetry that is consistent with one of the perturbed scenarios and expand on its properties. Finally, in Sec. VII we discuss the most relevant details of this work and conclude.

II. PERTURBATIVE MODIFICATIONS TO TRIBIMAXIMAL MIXING
Let us begin by recalling the form of the exact TBM mixing matrix [14] As mentioned before, this mixing matrix form has been the motivation for a great number of models that attempt to predict the neutrino oscillation parameters employing discrete symmetries. It is this sort of pattern with a vanishing 1-3 matrix element that were ruled out by the measurement of the non-zero reactor angle θ 13 . In this paper we will consider the following minimal perturbations to the TBM mixing matrix where the U ij (θ, φ) matrices are given by Finding the equivalence between the mixing matrix V of each case and the U PMNS can be The exact TBM pattern of Eq. (1) can be regarded as result of residual symmetries in the charged lepton and neutrino sectors from a flavor model defined at a higher energy. In this case, the mixing matrices of cases A and B in Eq. (2) can be considered the consequence of additional effects that break these residual symmetries on the planes (2,3) and (1,3) in the side of the neutrino sector, respectively. Similarly, cases C and D would break the residual symmetries on the side of the charged leptons on the planes (1,2) and (1,3), respectively.
Note that cases A and C were also studied in Refs. [33][34][35]. To simplify the notation, we will be using the shorthand s ij = sin θ ij and c ij = cos θ ij in the rest of the paper.

III. CP-VIOLATING PHASE FROM PERTURBATIVE SCENARIOS
One of the most relevant points of enabling a non-zero θ 13 is opening up the possibility of having a Dirac-type CP-violating phase in the PMNS mixing matrix. Due to the features of the cases mentioned in Eq. (2) it is possible to relate either θ 12 (cases A and B) or θ 23 (cases C and D) with the reactor angle θ 13 and, lastly, to relate the δ CP phase to the pair of free mixing angles. This is achieved by identifying the parameterizations that result from Eq. (2) with the standard PDG parameterization of the PMNS matrix. In this way, in Ref. [ following relations between oscillation parameters were worked out: A: , B: , C: , D: , where η 23 = 2 tan 2θ 23 and ξ = sin 2θ 12 . Considering the form of the matrices of Eq. (3) we can write the following expressions for the other oscillation parameters in terms of the angle θ and the phase φ: C, D: Note that, for every case, there is a relationship between θ 13 and θ, consistent with the idea that the matrices in Eq. (3) are perturbations that deviate θ 13 from zero. Using Eqs. (6) to (12), other noteworthy consequences of these perturbations include that for case A s 2 12 < 1/3, while for case B s 2 12 > 1/3, resulting in case B not being able to reproduce the current best fit value for this oscillation parameter. For cases C and D we obtain s 2 23 < 1/2 and s 2 23 > 1/2, respectively, meaning that whenever the octant of θ 23 is resolved at least one of these two cases will be ruled out.
A. Probability densities of cos δ CP Using the expressions in Eqs. (6) to (9) and the measured oscillation parameters from NuFIT 5.1 global fit [36], we can calculate probability densities for the predictions of the δ CP phase in every scenario. The process for calculating these densities follows Ref. [37].
There are three facts that simplify the process in the present case: 1. Eqs. (10) to (12) imply an upper bound on s 2 13 that is well above the acceptable experimental range and, thus, has no relevant effect in this analysis.
2. With the same equations, the values we can get for s 2 23 in cases A and B are not particularly limited by specific values of s 2 13 in the range of interest from the global fit, therefore, we can consider s 2 13 independent of s 2 23 .
3. Using input values around 3σ range for s 2 13 and s 2 23 in cases A and B predicts only physical values for cos δ CP . Points 2 and 3 above are also true for cases C and D replacing s 2 23 by s 2 12 . Considering these details, we can calculate the probability density for cos δ CP directly using the probability densities of s 2 13 , s 2 23 and s 2 12 . The integral that we have to perform to calculate the probability density at some particular value z of cos δ CP is given by where w, x, y represent values of s 2 12 , s 2 13 , s 2 23 , respectively, that we have to integrate over. The functions f j , with j ∈ {A,B,C,D}, represent cos δ CP for each case and the delta function ensures that the integration is performed over a line where cos δ CP = z. Independently of the three points enumerated before, the probability densities P s 2 ij can be any normalized function where the values of z are well defined in the integration intervals. Note that, in general, one of the two probability densities in each integral should be a conditional probability distribution dependent on the input of the other, however, given point 2 above, we are considering both distributions in each integral as independent.
For this work, we are interested in using Eqs. (13) and (14) to calculate cos δ CP probability densities from currently observed oscillation parameters. For this purpose we use the χ 2 tables provided by the NuFIT collaboration available on their website [38]. The data used corresponds to the normal (NO) and inverted (IO) ordering results that include Super-Kamiokande's tabulated χ 2 data (lower part of Table 3 in Ref. [36]), these have been collected in Table I  The changes in the distributions of cos δ CP from considering NO or IO data from Table I is mostly related to changes in central values and χ 2 projections. Cases C and D, which depend on s 2 13 and s 2 12 , do not change notably between using NO or IO data. However, in the case of A and B, due to the significant change in the χ 2 projection of s 2 23 , the distribution of cos δ CP changes to display two more leveled peaks with the higher peak changing side in both cases. This can be see in detail in Fig. 2

IV. PERTURBATIVE MODIFICATIONS TO BIMAXIMAL MIXING
In the same way we modified the TBM mixing case in Sec. II, we can apply perturbations to the very well known BM mixing [7][8][9][10][11][12][13]. The exact form of the BM mixing matrix is given by From here, perturbations proceed identically as for the TBM case. We can define the following scenarios with the U ij (θ, φ) matrices given in Eq. (3). Cases G and H were considered ruled out by experimental data above 3σ when they were studied on Ref. [18]. For cases E and F the expressions for s 2 23 are identical to those of cases C and D, respectively. The relationships between mixing angles and CP-violating phase are given by E: where η 12 = 2 tan 2θ 12 .
To provide an update for cases E and F, we find that they cannot predict physical values for cos δ CP within the 3σ ranges using current results from Ref. [36]. In Fig. 3 the 3σ rectangle for s 2 13 and s 2 12 is shown together with the closer physical boundary (colored contours) for the predicted cos δ CP for both cases E and F. Interestingly, in both panes of Fig. 3, the boundary of the physical predictions for cos δ CP is barely outside the 3σ rectangle, almost touching the upper right corner, indicating that this level of exclusion must be quite recent.
With these results, all the cases with U BM 0 considered in Ref. [18] can be considered ruled out at 3σ or above. Considering this, we will not follow the detailed analysis of the previous section on the CP-violating Dirac phase for the cases of this section and the rest of this work will be focused on cases A, B, C and D.

V. PROSPECTS AT FUTURE EXPERIMENTS
To analyze the prospects for the four cases considered in this work, we will employ simulations using the GLoBES software package [41,42]. We will consider three longbaseline experiments: DUNE, T2HK and ESSnuSB. For the DUNE experiment we consider the configuration detailed in their technical design report [43]. According to Ref. [43], the DUNE experiment is planned to have a long-baseline of 1300 km, with a 1. we assume that the 10 years run time is distributed with a ratio of 3:1 for antineutrino to neutrino modes. For ESSnuSB we consider the experimental setup outlined in Ref. [45].
The neutrino beam would be produced at the European Spallation Source with a power of 5 MW. Neutrinos would be received at a MEMPHYS-like [46] water Cherenkov detector with a (fiducial) mass of (507 kt) 1 Mt, at a distance of 540 km. With this configuration, ESSnuSB will reach 2.7 × 10 23 POT per year. In our analysis, we assume a run of 10 years with a ratio of 8:2 for antineutrino to neutrino modes as mentioned in the "Nominal value" column of Table 1.1 of Ref. [45].
The statistical analysis follows the methodology described in Sec. III of Ref. [47]. To summarize the steps in this methodology, we start with GLoBES χ 2 G function comparing the N obs events observed in the simulation of the experiment against N th events expected from theory. GLoBES χ 2 G function can be written as where (θ, φ) refers to a set of parameters in the theory and the summation run over bins.
Additionally, we include two Gaussian prior contributions to the total χ 2 using the reported central values, s 2 12,obs and s 2 13,obs , and their corresponding errors, σ 12 and σ 23 , given in Table I [36]. Considering that currently the octant of s 2 23 is not known, for its prior we use an interpolation of the χ 2 table provided in NuFIT's website [38]. The full χ 2 pr is given by The total χ 2 to be minimized is given by Following Ref. [47], our results will be presented for ∆χ 2 = χ 2 mod − χ 2 free , where χ 2 mod is the result of minimizing Eq. (21) over the model parameters θ and φ, while χ 2 free is the minimization over oscillation parameters ignoring constraints from the scenarios of Sec. II.
The results of our simulations are presented in Fig. 4  The experimental results used in the simulation correspond to normal ordering, however, there is no significant change for inverted ordering other than the current best fit point, indicated as a red thick × for NO and a blue-white + for IO. C and D, the compatible s 2 23 is more strongly constrained by DUNE and T2HK, while for all cases ESSnuSB reduces the δ CP phase range. Assuming that future experimental results will be close to the current best fit point (NO: red thick ×, IO: blue-white +), we can see that T2HK (red regions) alone could exclude cases C and D for both NO and IO at 5σ or more, while DUNE (blue regions) could exclude C and D only for the NO result, with the IO result remaining within 3 to 5σ. Under the same assumed future results, ESSnuSB could not exclude any case above 5σ. However, the combination of the three experiments (black contours) has the capacity of excluding cases B, C and D for both orderings to 5σ or more.
Case A has the best chances of survival, with a NO result disfavoured only between 3 to σ and IO staying well below 3σ.

VI. A MODEL WITH A 4 MODULAR SYMMETRY
In this section we will construct a model that predicts the neutrino masses and mixing within the measured limits, and we will show that symmetry breaking in this model results in a mixing pattern that is consistent with case A studied in previous sections.
The properties of modular forms are described in detail in Ref. [48]. To summarize the modular approach to flavor models, consider the group Γ(N ) defined by where SL(2, Z) is the special linear group of 2×2 matrices with integer elements and determinant equal to 1. The elements of the group Γ(N ) transform a complex variable τ , constrained by Im(τ ) > 0, according to  relation to the permutation groups mentioned, the interested reader may check Refs. [48][49][50].
The model that we develop here is based on modular forms of level N = 3 which have a quotient group, Γ 3 , isomorphic to A 4 , the symmetry group of the tetrahedron. For A 4 the generators have the properties The modular forms of level 3 were constructed on Appendix C of Ref. [48] and correspond to where η is the Dedekind eta function defined by and

A. Lepton masses
We will consider Majorana neutrinos that acquire small masses via the see-saw mechanism. This model is based on an extension by the symmetry group A 4 × U (1) X with the right handed neutrinos in a triplet of chiral supermultiplets N c . The field content of the model, A 4 representation, U (1) X charges and modular weights k I are collected in Table II.
With those assignments for the fields, the superpotential is given by where α i , β i and γ i are dimensionless couplings. In the case of α i and β i they can be made real by field redefinitions and are taken as O(1) coefficients. The γ i couplings are complex and we take their modulus as O(1). The weight 4 modular forms are given by At energies below the electroweak scale, the scalars χ, H u and H d acquire vacuum expectation value (VEV) giving masses to the fields in the superpotential of Eq. (28). The charged lepton masses can be extracted from the first line of the superpotential and correspond to the diagonal matrix where a ≡ Y 2 /Y 1 and b ≡ Y 3 /Y 1 . One can choose integers f e,µ,τ and χ /Λ in order for ( χ /Λ) fe−fτ = 0.0003 and ( χ /Λ) fµ−fτ = 0.06 to satisfy the empirical results of charged lepton masses. From the second and third lines in the superpotential we can read off the following mass matrices Note that the light mass matrix is proportional to H u 2 / χ , therefore, neutrino masses are expected to be small for very large χ .
Since the A 4 flavor models of leptons give large flavor mixing angles clearly [51,52], several A 4 modular invariant models have been proposed [53][54][55][56][57]. It may be useful to comment on the distinctive features of our model. Our charged lepton mass matrix is diagonal, in contrast to previous models, by assignment of A 4 singlets for both left-handed and right-handed leptons apart from the right-handed neutrinos. Then, the lepton mixing angles come from flavor structure of the neutrino mass matrices. Therefore, our model is advantageous for discussing the TBM and the case A in the context of A 4 flavor symmetry.
It is emphasized that the Dirac phase δ CP ≈ 3/2π could be reproduced around the fixed point τ = i as seen in the next subsection.

B. Perturbative modifications to TBM mixing
As mentioned in Sec. II, in the cases where the perturbation to TBM mixing is due to the breaking of residual symmetries in the neutrinos sector, such as the model presented here,  we can expect perturbations of the form of cases A or B.
Considering the constraints imposed on s 2 12 and mentioned before the start of Sec. III A, case B is unable to reproduce the current best fit value for s 2 12 . This leaves case A as the most appropriate candidate for realistic phenomenological studies. Here we will attempt to show that the A 4 model presented above can predict oscillation parameters that are consistent with case A and are in complete agreement with the current best fit limits summarized in Table I.
As a first step, we find a few parameter choices that give predictions with good agreement with current experimental values and are consistent with case A, characterized by Eqs. (6).
We provide a few benchmark points in Table III labeled as Aj as well as their predictions in Table IV. The next step is finding a neighboring point that reproduces TBM with good accuracy. Such point should be considered only illustrative, since TBM is in disagreement with current bounds, namely with the measured range for s 2 13 . The TBM points neighboring the Aj points are given in Table III with the label TBMj. Finally, we can compare these two types of points to assess how much each parameter changes between TBM and the perturbed case A.
We found that TBM mixing, particularly s 2 13 = 0, would require β 1 ≈ 0. Getting the other s 2 ij requires τ ≈ −0.16 + i1.0 and β 2 /β 3 ≈ 0.9. Note that the points TBMj in Table III do not use the value β 1 = 0 since it would make the first column of m D exactly zero. The overall factor of M ν is taken as H u 2 / χ |γ 1 | = 1.0092 × 10 −11 GeV. When comparing the points in Table III we see that β 1 and β 3 change the least remaining identical to 2 significant digits. Expectedly, β 1 changes the most since this parameter is related to the appearance of non-zero θ 13 . The parameter τ is dominated by its imaginary part which changes roughly 3%. The parameter γ is dominated by its real part and is the one that changes the most, and is mostly related to the requirement that the Aj points predict appropriate mass differences which we do not require from the TBM points for simplicity.
To conclude this section with a comment, while the A 4 model presented above permits mixing patterns far more complicated, the study of this section illustrates how case A may arise in a realistic model. Moreover, the relation that exists between case A and TBM mixing is made explicit in the comparison between model parameter values. This analysis is independent of the model and could be an starting point for a detailed study of the effects of breaking the residual symmetries that led to TBM mixing in the first place.

VII. CONCLUSION
In this work we revisit the perturbed mixing patterns that were considered in Ref. [18] for the popular BM and TBM mixings. Using current best fit values and 3σ ranges for the oscillation parameters we found that the considered perturbations to BM mixing, labeled E and F, cannot predict physical values for cos δ CP with s 2 12 and s 2 13 inside their 3σ ranges, while the four cases that consider perturbations of TBM mixing survived. We extended on previous efforts to predict the leptonic CP-violating Dirac phase by calculating distributions for its allowed values in light of the relations between oscillations parameters. For cases A, B and C we found that the preferred δ CP phase is located around 270 • while for case B the most favoured values spanned a range roughly from 200 • to 320 • . These values consider that, according to Ref. [36], the observed preferred range for δ CP is between 144 • and 350 • . Interestingly, planned experiments will have the power to constrain these simple perturbations, particularly cases B, C and D, which have the most constraining conditions.
In the case of B, s 2 12 > 1/3 is in tension with the currently measured value, and if future experiments keep this tendency we will see the tension increased. For cases C and D, due to each case predicting s 2 23 in different octants, one of them will be excluded when the octant problem is resolved. Nonetheless, both cases, C and D, predict s 2 23 quite close to 1/2, and if s 2 23 stays in close proximity to its current central value both cases could eventually be ruled out. The simulations performed and described in Sec. V show that DUNE2, T2HK and ESSnuSB experiments have the combined capacity to rule out cases B, C and D by more than 5σ, while case A could be left disfavoured by more than 3σ. We finalize by showing the emergence of case A from an A 4 modular symmetry flavor model. This model is capable of predicting currently measured oscillation parameters within their acceptable ranges. Moreover, we showed the existence of nearby points that predict TBM mixing to illustrate the degree of perturbation in the parameters required to obtain the mixing pattern of case A. The results of this study can be applied to any model that results in a mixing pattern consistent with the list in Eq. (2). Furthermore, any of the steps performed in this study could be applied to different neutrino masses and mixing models for which one can obtain relations like those in Eqs. (6) to (9), and may help reveal details brought about by the existence of such constraints.