Correlation of normal neutrino mass ordering with upper octant of θ23 and third quadrant of δ via RGE-induced μ-τ symmetry breaking

The recent global analysis of three-flavor neutrino oscillation data indicates that the normal neutrino mass ordering is favored over the inverted one at the 3σ level, and the best-fit values of the largest neutrino mixing angle θ23 and the Dirac CP-violating phase δ are located in the higher octant and third quadrant, respectively. We show that all these important issues can be naturally explained by the μ-τ reflection symmetry breaking of massive neutrinos from a superhigh energy scale down to the electroweak scale owing to the one-loop renormalization-group equations (RGEs) in the minimal supersymmetric standard model (MSSM). The complete parameter space is explored for the first time in both the Majorana and Dirac cases, by allowing the smallest neutrino mass m1 and the MSSM parameter tanβ to vary within their reasonable regions.


Introduction
The striking phenomena of solar, atmospheric, reactor and accelerator neutrino oscillations have all been observed in the past twenty years [1], demonstrating that the standard model (SM) of particle physics is by no means complete, and it must be extended to explain both the origin of tiny but finite neutrino masses and of large lepton flavor mixing effects. Qualitatively, the smallness of neutrino masses might be attributed to the existence of some heavy degrees of freedom at a superhigh energy scale: a popular idea known as the seesaw mechanism [2]. The largeness of neutrino mixing angles and CP-violating phases might originate from an underlying flavor symmetry [3,4], which should also manifest itself at a superhigh energy scale. Thus, a combination of the seesaw and flavor symmetry conjectures is the most likely phenomenological approach to understanding what is behind the observed spectrum of neutrino masses and the observed pattern of lepton flavor mixing. In this case, the renormalization-group equations (RGEs) are imperative to bridge the gap between the (theoretically suggestive) superhigh and (experimentally measurable) low energy scales [5,6].
Because the oscillation experiments are insensitive to the Dirac or Majorana nature of massive neutrinos, one may describe the link between the three known neutrinos (ν e , ν µ , and ν τ ) and their mass eigenstates (ν 1 , ν 2 , and ν 3 ) in terms of a 3 × 3 unitary matrix, namely the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix [7]  (1) with c ij ≡ cosθ ij and s ij ≡ sinθ ij (for ij = 12,13,23). A global analysis of currently available data on neutrino oscillations indicates that the normal neutrino mass ordering (m 1 < m 2 < m 3 ) is favored over the inverted one (m 3 < m 1 < m 2 ) at the 3σ level 4) . Furthermore, the best-fit value of the largest neutrino mixing angle θ 23 is slightly larger than 45 • (i.e., located in the higher octant), and the best-fit value of the Dirac phase δ is slightly smaller than 270 • (i.e., located in the third quadrant) [8,9]. Because θ 23 = 45 • and δ = 270 • can naturally be derived from the neutrino mass matrix M ν constrained by the µ-τ reflection symmetry 1) , a simple flavor symmetry that assures that M ν is invariant under proper charge-conjugation transformations of the left-and right-handed neutrino fields [4,12,13], the symmetry is expected to be the minimal flavor symmetry responsible for nearly maximal atmospheric neutrino mixing and potentially maximal CP violation in neutrino oscillations. If this simple discrete symmetry is realized at a superhigh energy scale, such as the seesaw scale with Λ µτ ∼10 14 GeV, then it will be broken at the electroweak scale Λ EW ∼ 10 2 GeV owing the RGE running effects, leading to deviations of (θ 23 ,δ) from (45 • ,270 • ). Whether such quantum corrections evolve in the right direction to fit the experimental results for θ 23 and δ at low energies depends on the neutrino mass ordering and the theoretical framework accommodating the RGEs [14][15][16]. Thus, this dependence provides an elegant possibility for correlating three major issues in today's neutrino physics: the neutrino mass ordering, the octant of θ 23 , and the strength of leptonic CP violation.
Given the fact that a credible global analysis of relevant experimental data often points to the truth in particle physics 2) , it is necessary for us to take the latest global-fit results of neutrino oscillations seriously, and to investigate their implications independently of any details of model building. In this study, we show that the normal neutrino mass ordering, the slightly higher octant of θ 23 , and the possible location of δ in the third quadrant can naturally be correlated and explained via RGE-induced µ-τ reflection symmetry breaking of massive neutrinos; that is, via the one-loop RGE evolution of neutrino masses and flavor mixing parameters from Λ µτ ∼10 14 GeV down to Λ EW ∼10 2 GeV in the minimal supersymmetric standard model (MSSM). This kind of correlation will soon be tested by more accurate experimental data.
Although the same topic has been discussed previously, the present study differs from previous works at least in the following aspects: • Based on more reliable experimental data, especially the 3σ indication of the normal neutrino mass ordering, our work is the first to numerically explore the almost complete parameter space in both the Majorana and Dirac cases, by allowing the smallest neutrino mass m 1 and the MSSM parameter tanβ to vary within their reasonable regions. In comparison, previous work such as Ref. [6] has only estimated the RGE correction to θ 23 by assuming that θ 13 = 0. This is an approximation that has now been excluded. Such work also did not consider all four distinct cases of the Majorana phases in the µ-τ reflection symmetry limit, nor the RGE correction to δ.
• While all previous studies have simply assumed some special values of m 1 and tanβ to "illustrate" a possible correlation of the neutrino mass ordering with the octant of θ 23 and the quadrant of δ based on the RGE-induced µ-τ reflection symmetry breaking, our present work has made remarkable progress in "figuring out" which part of the parameter space is favored by current neutrino oscillation data, and which part is unfavored or ruled out. The outcome of our in-depth analysis is therefore more up-to-date , suggestive, and useful for model building.
• Our statistical analysis shows that currently the best-fit points of θ 23 and δ [8,9] can be explained by the µ-τ reflection symmetry breaking induced by the RGE running from Λ µτ down to Λ EW , but this simple flavor symmetry itself is in slight disagreement (at the 1σ confidence level) with the low-energy data. Although the best-fit points may shift when new data becomes available, the analysis method itself will remain useful.
• The upcoming precision measurement of θ 23 and an experimental determination of δ will test the considered scenario and help to locate the correct regions in parameter space for the smallest neutrino mass m 1 , the MSSM parameter tanβ, and even the Majorana phases. No previous approaches were able to achieve this.
The rest of the paper is organized as follows. In Section 2 we present the main analytical results of the µ-τ reflection symmetry breaking caused by the RGE running effects. In Section 3 the parameter space of tanβ and m 1 is extensively explored and constrained with the help of current experimental data. Section 4 presents a summary and some concluding remarks.
1) As in most literature, in this case the so-called µ-τ reflection symmetry actually means the ν µ -ν τ reflection symmetry in the neutrino sector. One may refer to the discussions above Eq. (2) or Eq. (9) in Section 2, and to Refs. [10,11] for the building of specific models that realize this interesting discrete flavor symmetry without involving charged leptons.
2) One successful example of this kind was the global-fit "prediction" for an unsuppressed value of θ 13 , made in 2008 [18], which proved to be essentially true after the direct measurement of θ 13 was reported in 2012 [19]. This work has recently been recognized by the prestigious Bruno Pontecorvo Prize.

The Majorana case
Let us assume that the tiny masses of the three known neutrinos originate from a viable seesaw mechanism at a superhigh energy scale Λ µτ ∼ 10 14 GeV. Without loss of generality, we choose the basis in which the mass eigenstates of the three charged leptons are identical with their flavor eigenstates. In this case, only the neutrino sector is responsible for lepton flavor mixing and CP violation. If the effective Majorana neutrino mass term is invariant under the charge-conjugation transformations ν eL ↔ν c eR , ν µL ↔ ν c τ R and ν τ L ↔ ν c µR , then the corresponding mass matrix must take the form constrained by m ee = m * ee , m eµ = m * eτ , m µµ = m * τ τ , and m µτ = m * µτ [4]. One may diagonalize a generic 3×3 Majorana mass matrix through U † M ν U * = Diag{m 1 ,m 2 ,m 3 }, where the unitary matrix U can be decomposed as U =P l V P ν , with P l =Diag{e iφ e ,e iφ µ ,e iφ τ } and P ν = Diag{e iρ ,e iσ ,1} being the phase matrices, and V as shown in Eq. (1). Taking account of the µ-τ reflection symmetry of M ν , we immediately arrive at the following constraints on U : θ 23 = 45 • , δ = 90 • or 270 • , ρ = 0 • or 90 • , and σ = 0 • or 90 • for the four physical parameters, as well as φ e = 90 • and φ µ +φ τ = 0 • for the three unphysical phases. Because a clear preference for sinδ <0 has been obtained from a recent global analysis [8], it is fairly reasonable for us to focus only on δ=270 • in the µ-τ symmetry limit.
In the framework of the MSSM, the evolution of M ν from Λ µτ down to Λ EW through the one-loop RGE can be expressed as [20] with T l =Diag{I e ,I µ ,I τ }, in which Here, χ = ln µ/Λ µτ , where µ is an arbitrary renormalization scale between Λ EW and Λ µτ , g 1 and g 2 denote the gauge couplings, and y t and y α (for α = e,µ,τ ) denote the Yukawa coupling eigenvalues of the top quark and charged leptons, respectively. The smallness of y e and y µ assures that I e ≃I µ ≃1 holds to an excellent degree of accuracy. Note that is also a fairly small quantity, but may affect the running behaviors of the flavor mixing parameters in an appreciable manner [20]. To illustrate, Fig. 1 presents the twodimensional maps of ∆ τ (left panel) and I 0 (right panel) versus the variables Λ µτ and tanβ. One can observe that I 0 does not change significantly with different settings of Λ µτ and tanβ. In comparison, ∆ τ can change from 0.001 to 0.05. Note that shifting the energy scale is equivalent to altering tanβ, and the outputs at Λ µτ = 10 14 GeV and Λ µτ = 10 16 GeV are reasonably similar in magnitude. If we shift the energy scale from Λ µτ =10 14 GeV to Λ µτ =10 9 GeV, then ∆ τ will lie in the range (0.001,0.03) instead of (0.001,0.05). In the following numerical cal-culations, we shall fix Λ µτ ∼ 10 14 GeV as the µ-τ flavor symmetry scale. One may diagonalize the neutrino mass matrix at Λ EW , and then obtain the mass eigenvalues (m 1 ,m 2 ,m 3 ), flavor mixing angles (θ 12 ,θ 13 ,θ 23 ), and CP-violating phases (δ,ρ,σ). Here, we define ∆θ ij ≡ θ ij (Λ EW ) − θ ij Λ µτ (for ij = 12,13,23), ∆δ ≡ δ(Λ EW ) − δ Λ µτ , ∆ρ ≡ ρ(Λ EW )−ρ Λ µτ , and ∆σ ≡ σ(Λ EW )−σ Λ µτ to measure the strengths of the RGE-induced corrections to the parameters of U . As a good approximation, the three neutrino masses at Λ EW are found to be where θ 12 and θ 13 take their values at Λ EW . Unless otherwise specified, the nine physical flavor parameters appearing in the subsequent text and equations are all those at Λ EW . In a reasonable analytical approximation, we can also arrive at 1) , 1) Note that our analytical results are not exactly the same as those obtained in Ref. [15], where a different phase convention for the PMNS matrix has been adopted.

The Dirac case
Although the Majorana nature of massive neutrinos is well motivated from a theoretical viewpoint, it is also interesting to consider the possibility of a pure Dirac mass term for the three known neutrinos, and combine this with a certain flavor symmetry that can be realized at a superhigh energy scale Λ µτ [17]. In this case, the µ-τ reflection symmetry means that the Dirac neutrino mass term is invariant under the charge-conjugation transfor- and is constrained by m ee = m * ee , m eµ = m * eτ , m µe = m * τ e , m µτ = m * τ µ , and m µµ = m * τ τ [17]. One can observe that in general M ν is neither symmetric nor Hermitian, even though it possesses the µ-τ reflection symmetry. However, one may diagonalize the Hermitian combination M ν M † ν by means of the unitary transfor- e iφ µ ,e iφ τ } being an unphysical phase matrix and V taking the form given in Eq. (1). Then, the µ-τ reflection symmetry constraints on M ν naturally lead us to the constraints on the PMNS matrix U ; namely, Given the global-fit preference for sinδ < 0 [8], we focus only on the possibility that δ(Λ µτ )=270 • .
In the MSSM, the evolution of the Dirac neutrino mass matrix M ν from Λ µτ down to Λ EW via the one-loop RGE can be described as [17] where the definitions of I 0 and T l are the same as those in Eqs. (3) and (4), respectively. Diagonalizing M ν (Λ EW ) will yield the seven physical flavor parameters at Λ EW . Utilizing the same notations as in the Majorana case, we summarize our approximate analytical results as follows: for the masses of the three neutrinos, for the deviations of the three flavor mixing angles between Λ EW and Λ µτ , and for the deviation of δ between Λ EW and Λ µτ . In obtaining Eqs. (11)-(13), we have taken into account the µ-τ reflection symmetry conditions θ 23 (Λ µτ ) = 45 • and δ(Λ µτ )=270 • . The analytical approximations presented in Eqs. (7) and (8) for Majorana neutrinos or in Eqs. (12) and (13) for Dirac neutrinos are instructive, and helpful for understanding the RGE corrections to relevant flavor mixing and CP-violating parameters. However, the accuracy will be rather poor if the neutrino masses are strongly degenerate. In the subsequent section, we shall numerically evaluate the effects of µ-τ reflection symmetry breaking, and explore the allowed parameter space to fit current experimental data.

Numerical exploration of the parameter space
In the framework of the MSSM, we numerically run the RGEs from Λ µτ ∼10 14 GeV down to Λ EW ∼10 2 GeV by taking into account the initial conditions θ 23 =45 • and δ = 270 • , as well as the initial values of ρ and σ in the following four different cases: Case A ρ=σ =0 • , Case B ρ=σ=90 • , Case C ρ=0 • and σ=90 • , and Case D ρ=90 • and σ = 0 • . For any given values of the MSSM parameter tanβ and the smallest neutrino mass m 1 at Λ EW , the other relevant neutrino oscillation parameters, such as {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm } at Λ µτ are scanned over suitably wide ranges with the help of the Multi-Nest program [21]. Here, we have adopted the notations ∆m 2 sol ≡m 2 2 −m 2 1 and ∆m 2 atm ≡m 2 3 −(m 2 1 +m 2 2 )/2, as defined in Ref. [8] 1) . For each scan, the neutrino flavor parameters at Λ EW are obtained, and these are immediately compared with their global-fit values by minimizing where the ξ i terms denote the parameters at Λ EW produced from the RGE evolution, the ξ i represent the bestfit values of the global analysis, and the σ i are the corresponding symmetrized 1σ errors (i.e., σ i = (σ + i +σ − i )/2). In our numerical calculations, we concentrated on the normal neutrino mass ordering 2) and adopt the best-fit values and 1σ-level deviations of the six neutrino oscillation parameters [8]: In this case, the smallest neutrino mass m 1 is allowed to take values in the range [0,0.1] eV, and the MSSM parameter tanβ may vary from 10 to 50 based on a reasonable phenomenological argument 3) . It should be pointed out that our numerical results are independent of the analytical approximations obtained in the previous section, but the latter will be helpful for understanding some salient features of the former.

The Majorana case
The strategy of our numerical analysis is rather straightforward. Let us first examine how significantly θ 23 and δ at Λ EW can deviate from their initial values at Λ µτ , incorporating the recent global-fit results. To this end, we only need to take into account the 1) In the normal neutrino mass ordering case, one may therefore express m 2 and m 3 in terms of m 1 as m 2 = m 2 1 +∆m 2 sol and m 3 = m 2 1 +0.5∆m 2 sol +∆m 2 atm . 2) This is not only our phenomenological preference but also a recent 3σ global-fit indication that the true neutrino mass spectrum should exhibit a normal ordering m 1 <m 2 <m 3 .
3) Note that tanβ >50 is unfavored, because the heavy-quark Yukawa couplings would fall into the non-perturbative region, while for tanβ <10 the RGE-induced corrections to the relevant neutrino parameters are negligibly small, and thus less interesting for our purpose.
global-fit information concerning the parameter set ξ = {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm }. One can observe that the RGE running effects always push θ 23 to the higher octant, and in most cases lead δ to the third quadrant. This is the right direction as indicated by the best-fit values of these two quantities [8]. After this preliminary diagnosis is reached, the experimental information on θ 23 and δ will be included in order to evaluate the statistical compatibility between the RGE-triggered µ-τ reflection symmetry breaking and the global-fit values of θ 23 and δ. For each given value of m 1 or tanβ at Λ EW , we obtain the associated χ 2 min from Eq. (14), which is minimized over the chosen parameter set ξ. In Fig. 2, we illustrate χ 2 min with respect to m 1 for Case C as an example. It is obvious that χ 2 min can reach 0 (i.e., the best-fit point) in most cases, but for m 1 ∈[0.05,0.08] eV the value of χ 2 min increases to nearly 100. Here, χ 2 min > 0 means that the µ-τ reflection symmetry limit at Λ µτ cannot be reached if one runs the RGEs inversely, starting from the best-fit points of the six neutrino oscillation parameters at Λ EW . This observation will also be true even if one allows θ 23 and δ to take arbitrary values at Λ EW . The reason for this should be ascribed to the nontrivial differential structures of the RGEs [5,6], especially when the evolution becomes wild in a narrow parameter space. In such a case, the Majorana phases ρ and σ play a potentially significant role. We have demonstrated that χ 2 min in the bumped region is actually the local minimum. The reliability of this result is also ver-ified by running the system from Λ EW up to Λ µτ . The squared distance to the µ-τ reflection symmetry limit η≡{θ 23 =45 • ,δ=270 • ,ρ=0 • ,σ=90 • } at Λ µτ is defined as where the relevant parameters η i are obtained by running the RGEs from Λ EW up to Λ µτ with {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm } all taking their bestfit values at Λ EW . As shown in Fig. 2, D 2 min exhibits a similar bump, the range of which coincides with that of χ 2 min . In Eq. (16), the normalization factor 1/1000 is introduced simply to make the amplitude of the D 2 bump the same order of magnitude as that of χ 2 min . This check confirms the observation we obtained by running the RGEs from Λ µτ down to Λ EW . Figure 3 depicts the RGE-corrected result of θ 23 at Λ EW for different values of m 1 and tanβ with χ 2 = 0, where the four possible initial options for (ρ,σ) at Λ µτ have been considered. Note that the boundary conditions for the RGEs include both the initial values of {θ 23 ,δ,ρ,σ} at Λ µτ and the experimental constraints on {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm } at Λ EW , which are all specified in our numerical calculations. In this case, one should keep in mind that the low-and high-scale boundary requirements are likely to be sufficiently strong that the RGEs do not yield a realistic solution, illustrated by the gray-gap regions in Fig. 3 (Cases C and D), corresponding to the χ 2 min > 0 bump region in Fig. 2. Some additional comments on Fig. 3 are in order: • The gray-gap regions in Cases C and D are a result of the χ 2 min -bump, as shown in Fig. 2, but there is no such gap for Cases A and B. In Case A, the RGE running effect is quite significant, and θ 23 may run to almost 75 • for m 1 ≃ 0.1 eV and tanβ ≃ 50. Of course, such an evolution will be strongly constrained by the experimental information on θ 23 , which has not yet been included in our analysis. In contrast, θ 23 is not sensitive to the RGE corrections in Case B, and it maximally changes by only approximately 1 • . For Cases C and D, if one conservatively requires that m 1 < 0.07 eV from the present cosmological bound [22], then it will be possible for θ 23 to run to 46.6 • and 50 • , respectively. The best-fit value θ 23 ≃48 • [8] can easily be reached in Cases A and D.
• The RGE correction to θ 23 illustrated in Fig. 3 can be effectively understood with the help of the analytical approximations obtained in Eq. (7) if the neutrino masses are not too degenerate and the RGE evolution is not too strong. Eq. (7) tells us that the sign of ∆θ 23 is positive, because ∆ τ and ζ 31 ≃ ζ 32 are all positive for the normal neutrino mass ordering. The factor ∆ τ is essentially proportional to tan 2 β, because y 2 τ ∝(1+tan 2 β)≃tan 2 β for tanβ 10, and therefore ∆θ 23 always increases with tan 2 β. On the other hand, the dependence of θ 23 on the neutrino mass m 1 is different for the four options of ρ and σ at Λ µτ . For example, ζ −ηρ 31 is proportional to m 1 for ρ(Λ µτ ) = 0 • , but is proportional to 1/m 1 when ρ(Λ µτ )=90 • . In the region of small m 1 and tanβ, the radiative correction to θ 23 is proportional to m 1 for Cases A, C, and D, but is inversely proportional to m 1 in Case B with η ρ =η σ =−1.
In Fig. 4, we plot the allowed region of δ at Λ EW . Note that for each point in the m 1 -tanβ plane, δ and θ 23 are determined concurrently. Some remarks are in order.
• The RGE-induced corrections to δ in Cases A and B are very weak, being only approximately 0.5 • . Even though the higher octant of θ 23 (including its best-fit value) can easily be reached in Case A, it is impossible to approach the best-fit value of δ (i.e., δ ≃ 238 • ). However, the best-fit value of δ can be reached in both Cases C and D. Case D warrants highlighting, where the best-fit point (θ 23 ,δ)≃(48 • ,238 • ) is reachable from the same settings of m 1 and tanβ.
• Similar to the case of θ 23 , the radiative correction to δ is also proportional to tan 2 β regardless of its sign. However, the dependence of δ on m 1 in Eq. (8) is not as straightforward as that of θ 23 in Eq. (7). There are two terms in ∆δ, with one enhanced by 1/sinθ 13 and the other suppressed by sinθ 13 , but the latter can become dominant in some cases. In Case A, the first term ∝ 1/sinθ 13 is positive and dominant when the neutrino mass m 1 is relatively small, while the second term ∝ sinθ 13 is negative, and will gradually dominate when the value of m 1 increases. These analytical features can explain the numerical evolution behavior of δ for Case A shown in Fig. 4. In Case B with η ρ = η σ = −1, both terms of ∆δ are negative and inversely proportional to m 1 . Note that the first term of ∆δ in either Case A or B is suppressed, owing to the cancellation between ζ −η σ 32 and ζ −η ρ 31 , and this largely explains the smallness of ∆δ. In Case C (or Case D), the first term of ∆δ is negative (or positive) and initially dominant, but the second term containing ζ ηρησ 21 = ζ −1 21 will eventually dominate for relatively large values of m 1 . Hence, the RGE-induced corrections to δ in these two cases can be considerably enhanced by 1/∆m 2 sol .
To see the correlation between θ 23 and δ at Λ EW , let us marginalize m 1 and tanβ over the reasonable ranges m 1 ∈ [0,0.1] eV and tanβ ∈ [10,50]. Our numerical outputs are summarized in Fig. 5, where the recent global-fit results [8] are plotted as black lines with the 1σ (dashed), 2σ (solid), and 3σ (dotted) contours. The corresponding best-fit point of each analysis is marked as a black star. We note that the µ-τ reflection symmetry point (θ 23 ,δ) = (45 • ,270 • ) at Λ EW , which is marked by the red circled cross in the plot, lies on the dashed contour. This means that θ 23 (Λ EW )=45 • and δ(Λ EW )=270 • are statistically unfavored at the 1σ level [8]. The θ 23δ correlation at low energies, which arises from RGEinduced µ-τ symmetry breaking, is described by the pink or green region. In the pink region, the best-fit values of {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm } can be simultaneously reached (i.e., χ 2 min = 0). If the value of χ 2 min is relaxed to 11.83 (i.e., the 3σ confidence level for two degrees of freedom), then the wider green region for θ 23 and δ will be allowed. In the two upper panels of Fig. 5, which correspond to Cases A and B, the allowed range of δ is very narrow. This feature is compatible with the two upper panels of Fig. 4, where δ varies less than 1 • . In these . The pink region is allowed for θ 23 and δ when {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm } at ΛEW take their best-fit values, and the green region is allowed when these four observables deviate from their best-fit values by a 3σ level (i.e., χ 2 =11.83 for two degrees of freedom). two cases, the green region almost overlaps with the pink region. In the two lower panels of Fig. 5, corresponding to Cases C and D, the RGE-induced corrections are significant. Note that there is a separate shaded region around θ 23 ≃50 • in Fig. 5(C), and this is associated with the small upper-right corner of the parameter space in Fig. 3(C) or Fig. 4(C). There is a similar separate shaded region in Case D, but this is outside the chosen ranges of θ 23 and δ in plotting Fig. 5(D), and its confidence level is significantly weaker, being outside the 3σ region of the global analysis.
To numerically verify the compatibility between our µ-τ symmetry breaking scenario and the current experimental data, or how well the RGE-triggered µ-τ reflection symmetry breaking effect can fit the complete set of data at Λ EW , now we include the globalfit information on θ 23 and δ. Namely, we defineχ 2 ≡ χ 2 +χ 2 θ 23 +χ 2 δ , where χ 2 is formed with the parameter set {sin 2 θ 12 , sin 2 θ 13 , ∆m 2 sol , ∆m 2 atm } as before, and χ 2 θ 23 and χ 2 δ are the contributions from θ 23 and δ, respectively. In this case, the minimum ofχ 2 can be calculated for each of the four cases of ρ and σ at Λ µτ by marginalizing the relevant quantities over tanβ ∈ [10,50], and the result is plotted in Fig. 6 as a function of m 1 . For the special point (θ 23 ,δ)=(45 • ,270 • ) at Λ EW , the corresponding χ 2 min value is 2.08. Fig. 6 demonstrates the very effective reduction ofχ 2 min by incorporating the RGE running effect in the framework of the MSSM. Among the four cases under discussion, the red curve for Case D with ρ = 90 • and σ = 0 • at Λ µτ is most favorable. Even given the Planck limit on the sum of the neutrino masses ≡ m 1 +m 2 +m 3 < 0.23 eV at the 95% confidence level [22], it is still possible to reduce the value ofχ 2 min to almost 0.05.

The Dirac case
Because there is only a single CP-violating phase in the Dirac case, it is considerably easier to perform a numerical analysis of the parameter space, which is constrained by both the RGE-induced µ-τ reflection symmetry breaking effect and the recent global-fit results of neutrino oscillation data. Fig. 7 illustrates the allowed regions of θ 23 and δ at Λ EW , and their intimate correlation is illustrated in Fig. 8. Here, the evolution behaviors of these two parameters can be understood in a more straightforward manner than in the Majorana case, simply because of the absence of the two Majorana phases. In the leading-order approximation, the analytical expressions of ∆θ 23 and ∆δ in Eqs. (12) and (13) It becomes obvious that (θ 23 ,δ) may have larger deviations from (45 • ,270 • ) for bigger values of m 1 and tanβ. In particular, θ 23 and δ are always located in the upper octant and third quadrant, respectively. The gray regions in Fig. 7 exist for the same reason as in the Majorana case, as discussed above. Similar to Case D in the Majorana scenario, the RGE-induced corrections can take (θ 23 ,δ) very close to their best-fit point, and the correspondingχ 2 min is represented in Fig. 6 by the gray dashed line.

Concluding remarks
In neutrino physics, it is usually necessary (and popular) to introduce some heavy degrees of freedom and certain flavor symmetries at a superhigh energy scale in order to explain the tiny masses of the three known neutrinos and the striking pattern of lepton flavor mixing observed at low energies. In this case, it is also necessary to employ the RGEs as a powerful tool to bridge the gap between these two considerably different energy scales. Such RGE-induced quantum corrections may naturally break the given flavor symmetry, thus leading to some phenomenologically interesting consequences, including a possible correlation between the neutrino mass ordering and flavor mixing parameters.
In this work, we have considered the intriguing µτ reflection symmetry and its RGE-induced breaking as an instructive playground to realize the aforementioned idea, especially in view of the fact that the recent global analysis of neutrino oscillation data has indicated a remarkable preference for the normal neutrino mass ordering at the 3σ level, together with a slightly higher octant of θ 23 and the possible location of δ in the third quadrant. We have shown that all these important issues can be naturally correlated and explained by the RGE-triggered µ-τ reflection symmetry breaking from a superhigh energy scale Λ µτ ∼10 14 GeV down to the electroweak scale Λ EW ∼10 2 GeV in the MSSM. Unlike previous attempts along these lines, our study represents the first numerical exploration of the complete parameter space in both the Majorana and Dirac cases, by allowing the smallest neutrino mass m 1 and the MSSM parameter tanβ to vary within their respective reasonable regions [0,0.1] eV and [10,50]. We believe that direct measurements of the neutrino mass ordering and precision measurements of θ 23 and δ in the near future will test the simple but suggestive scenario under consideration.
Of course, some of our main observations are subject to the MSSM itself, and the current best-fit values of θ 23 and δ. The reason why we have chosen the MSSM instead of the SM is three-fold: (a) in the SM it is ex-tremely difficult to generate an appreciable value of ∆θ 23 via the RGE-induced µ-τ symmetry breaking effect, no matter which neutrino mass ordering is considered; (b) in the SM the running direction of θ 23 from Λ µτ down to Λ EW seems to be "wrong" if one takes today's bestfit result θ 23 >45 • seriously in the normal mass ordering case; and (c) the SM itself may suffer from the vacuumstability problem, as the energy scale is above 10 10 GeV [23]. When the two-Higgs doublet models (2HDMs) are considered [24], one may perform a similar analysis to reveal the correlation between the neutrino mass ordering, the octant of θ 23 , and the quadrant of δ via the RGEinduced breaking effects. The deviations of θ 23 and δ from their values in the µ-τ reflection symmetry limit can be quite different from those in the MSSM case, depending largely on which 2HDM scenario is taken into account. Ref. [17] has provided an explicit example of this kind in the type-II 2HDM scenario.
On the other hand, we admit that the best-fit values of θ 23 and δ will unavoidably "fluctuate" in the coming years, when more accurate experimental data are accumulated and incorporated into the global analysis framework. That is why we have numerically explored the complete parameter space to illustrate the tolerable ranges of m 1 and tanβ that allow us to correlate the normal neutrino mass ordering with the higher octant of θ 23 and the third quadrant of δ. If the inverted neutrino mass ordering, the lower octant of θ 23 , and (or) another quadrant of δ turned out to be favored by future precision measurements, it would be straightforward to consider a different correlation scenario for them, either within or beyond the MSSM. In the same spirit, one may study other interesting flavor symmetries and their RGE-induced breaking, in order to effectively link model building at high-energy scales to neutrino oscillation experiments at low energies.