Nonstandard interaction effects on neutrino parameters at medium-baseline reactor antineutrino experiments

Precision measurements of leptonic mixing parameters and the determination of the neutrino mass hierarchy are the primary goals of the forthcoming medium-baseline reactor antineutrino experiments, such as JUNO and RENO-50. In this work, we investigate the impact of nonstandard neutrino interactions (NSIs) on the measurements of {sin^2 theta_{12}, Delta m^2_{21}} and {sin^2 theta_{13}, Delta m^2_{31}}, and on the sensitivity to the neutrino mass hierarchy, at the medium-baseline reactor experiments by assuming a typical experimental setup. It turns out that the true mixing parameter sin^2 theta_{12} can be excluded at a more than 3 sigma level if the NSI parameter varepsilon_{e mu} or varepsilon_{e tau} is as large as 2% in the most optimistic case. However, the discovery reach of NSI effects has been found to be small, and depends crucially on the CP-violating phases. Finally, we show that NSI effects could enhance or reduce the discrimination power of the JUNO and RENO-50 experiments between the normal and inverted neutrino mass hierarchies.


I. INTRODUCTION
The experimental endeavor in the past decades has established that the phenomenon of neutrino flavor transitions is described by neutrino oscillations at leading order. Now that the smallest leptonic mixing angle θ 13 has been measured very accurately at reactor [1][2][3][4] and accelerator [5] neutrino experiments, one of the next major and open problems is the discrimination between normal and inverted neutrino mass hierarchies. Towards this challenge, different types of experiments have been suggested using long-baseline acceleratorbased neutrinos, atmospheric neutrinos, supernova neutrinos, or reactor neutrinos. Recent studies indicate that a dedicated medium-baseline reactor antineutrino experiment with sufficient statistics and unprecedented detector performance provides an opportunity to determine the neutrino mass hierarchy and probe with high precision the other neutrino parameters (see, e.g., Ref. [6] and references therein). In the near future, potential projects along this direction include the JUNO [7] and RENO-50 [8] experiments.
Beyond the standard oscillation picture, new physics may appear in future neutrino experiments in the form of unknown couplings involving neutrinos, which are usually referred to as nonstandard neutrino interactions [9][10][11][12][13] (NSIs). In fact, NSIs are predicted as dimension-six and higher-order operators in many interesting extensions of the Standard Model, e.g., R-parity violating supersymmetric theories, left-right symmetric models, grand unification theories, extra dimensions, and various seesaw models (see, e.g., the recent review on NSIs [14] and references therein). Basically, all modern extensions could give rise to NSIs, and the investigation of NSIs could be very helpful in revealing additional new physics behind neutrino flavor transitions. In addition, it plays an important complementary role to direct searches of new physics at colliders. NSI effects at reactor antineutrino experiments have been studied with emphasis on mimicking effects [15], the discrepancy between source and detector NSIs [16], comparisons between reactor and accelerator experiments [17,18], and modifications of event rates and the impact on the measurements of standard oscillation parameters [19]. On the other hand, the model-independent bounds on production and detection NSIs have been derived in Ref. [20], and all bounds are basically of order 10 −2 .
In this work, we investigate NSI effects at a typical medium-baseline reactor antineutrino experiment, in particular, the mimicking effects in the precision measurements of neutrino parameters, the discovery reach of NSI effects, and the distortion of the neutrino mass hierarchy determination. Specifically, we concentrate on an experimental setup similar to the JUNO experiment and show that, for five years of running, a clear hint for nonstandard neutrino physics can be provided. The remainder of this Letter is organized as follows.
In Section II, we present the general formalism and derive relevant antineutrino survival probability formulas used in the subsequent analysis. Illustrations of NSI effects on the corresponding oscillation probabilities and energy distributions of neutrino events are also given. In Section III, a detailed numerical analysis of the observability of the NSI effects at JUNO is performed. Finally, in Section IV, we summarize and state our conclusions.

II. BASIC FORMALISM
In this section, we present the general formulas for reactor antineutrino oscillations with NSIs. For a realistic experiment, NSIs may appear both in the production and detection processes, and the neutrino states produced in the source and observed at the detector can be treated as superpositions of pure orthonormal flavor states: where the superscripts 's' and 'd' denote the source and the detector, respectively, with the normalization factors 1 In general, the NSI parameter matrices ε s and ε d are arbitrary and nonunitary, indicating that neither |ν s α nor ν d β | are orthonormal states. As first observed in Ref. [22], heavy neutrino states responsible for neutrino mass generation decouple from oscillation processes, so that the unitarity of the leptonic mixing matrix is slightly violated. The nonunitary 1 In the calculation of the number of events in detectors, the normalization factors are canceled with the NSI factors in charged-current cross sections. However, for reactor antineutrino experiments, the neutrino fluxes are extracted from the measurement instead of a Monte Carlo simulation. Hence, the normalization factors should be taken into account. See also Ref. [21] for a detailed discussion. effects, which could be significant in low-scale seesaw models [23,24], can be regarded as one type of NSIs with the requirement ε s = ε d † [25]. For reactor antineutrinos, the leadingorder NSIs are of the V ± A type as long as CPT is conserved, and it is very common to assume ε s eα = ε d * αe [17]. Thus, we will take ε s eα = ε d * αe = ε eα e iφeα (with ε eα being the modulus of ε s eα ) in the current consideration and neglect the superscripts 's' and 'd' throughout the following parts of this work. Furthermore, the typical energy of antineutrinos produced in nuclear reactors is around a few MeV, which indicates that Earth matter effects are very small and can safely be neglected in most situations. Hence, the propagation of neutrino flavor states is governed by the vacuum Hamiltonian where the leptonic mixing matrix U is usually parametrized in the standard form by using three mixing angles and one CP-violating phase with c ij ≡ cos θ ij and s ij ≡ sin θ ij (for ij = 12, 13 and 23).
Including the NSI effects, we arrive at the amplitude for theν e →ν e oscillation channel where L is the propagation distance and the explicit form of A is a coherent sum over the contributions of all mass eigenstates ν i , i.e., A αβ = i U * αi U βi exp (−im 2 i L/2E). Inserting A αβ into Eq. (7), one can rewrite the amplitude in a compact form with It is a general feature from the first row of Eq. (9) that R i is real for the survival probabilities, which is actually a reflection of CPT conservation. Therefore, the amplitudes are identical for neutrino and antineutrino oscillations. However, this is no longer true when matter effects are non-negligible, since the ordinary Earth environment is CP asymmetric [26,27].
Furthermore, in the limit ε eα → 0, we have R i → |U ei | 2 , and hence, Eq. (8) is reduced to the standard oscillation amplitude.
In practice, it is more useful to express the R i parameters in terms of the leptonic mixing parameters, i.e., where the auxiliary parameters ε φ and ε δ are defined as One observes that the ε ee -dependent terms disappear from the survival probability, as a consequence of our chosen normalization. Moreover, if all CP-violating phases are vanishing, the auxiliary NSI parameters {ε φ , ε δ } are related to the original ones {ε eµ , ε eτ } by a rotation with the angle θ 23 . Now, it is straightforward to derive the survival probability for electron antineutrinos where ∆m 2 31 ≡ m 2 3 − m 2 1 > 0 for normal neutrino mass hierarchy (NH) and ∆m 2 31 < 0 for inverted neutrino mass hierarchy (IH). Here and henceforth, the hat in the antineutrino state |ν e will be omitted, hopefully without causing any confusion.
Then, considering a medium-baseline neutrino experiment, such as JUNO or RENO-50, we find that Eq. where the higher-order terms O(ε 2 ) have been neglected. Note that Eq. (16) reproduces the standard survival probability in the limit ε φ , ε δ → 0. It is worthwhile to note that ε φ and ε δ can be vanishingly small if, e.g., φ eµ = φ eτ = ±π/2 and δ = 0. In this case, only the higherorder terms O(ε 2 ) appear in the survival probability. In addition, the difference between the oscillation probability in the NH case and that in the IH case, i.e., ∆P ≡ P NH (ν e → ν e ) − P IH (ν e →ν e ), can be written as ∆P = 4s 13 c 13 s 13 c 13 s 2 12 + sin 2θ 12 s 13 ε φ + 2s 2 12 cos 2θ 13 ε δ sin indicating the impact of the NSI parameters on the discrimination between NH and IH.
For illustration, we show the survival probability P (ν e →ν e ) with and without NSIs in Fig. 1. The NSI polluted survival probability deviates from the standard one especially at energies around 3 MeV. For ε eτ = 0.03, the minimum of the probability is slightly larger than in the standard case, indicating that the measured value of θ 12 should be smaller than its true value. Such a feature can be understood from the oscillation probability in Eq. (16) that a negative ε φ is expected for ε eτ = 0.03 and ε eµ = φ eτ = δ = 0, reducing the amplitude of the ∆m 2 21 -driven oscillation. Consequently, one requires a smaller θ 12 compared to the true input parameter in order to reduce the effective amplitude so as to gain a better fit.
On the other hand, a positive ε δ enhances the amplitudes of the ∆m 2 31 and ∆m 2 32 -driven oscillation modes, suggesting that a larger θ 13 can be extracted from the experimental data.
These observations are also confirmed later by our numerical analysis in Section III B.
If ε eµ = 0.03 is assumed instead, the NSI-polluted survival probability becomes smaller than the standard one around the oscillation maximum associated with ∆m 2 21 , which can be understood from the opposite signs in front of ε eµ and ε eτ in Eq. (13). In an analogous way, one arrives at the conclusion that a larger value of θ 12 or θ 13 than its true value is expected from the fit to experimental data. However, at the probability level, there will be no visible differences between the NH and IH cases, which applies to the case with either ε eµ = 0.03 or ε eτ = 0.03. Finally, we remark that the oscillation frequencies are not affected by the NSIs, which can be observed from Eq. (16). Therefore, the NSIs affect the determination of the neutrino mass hierarchy at a medium-baseline reactor experiment by modifying the oscillation amplitudes instead of the oscillation frequencies.

III. NUMERICAL ANALYSIS
We proceed to perform a numerical analysis of NSI effects on the measurements of neutrino parameters. To this end, we employ the GLoBES software [28,29] and consider the configuration of the JUNO experiment described in detail in Ref. [7]. Explicitly, we take the reactor thermal power to be P = 35.8 GW and the baseline length L = 52.5 km. The detector mass is assumed to be 20 kt together with the energy resolution 3%/ √ E as a benchmark, where the energy is given in units of MeV. As discussed in Ref. [7], such an experimental setup allows one to determine the neutrino mass hierarchy at a confidence level of 4σ for a six-year running (with 300 effective days per year). We further make use of a single overall factor f together with the error δf = 3% to parametrize the uncertainties of the total reactor antineutrino flux, the cross section of inverse beta decay, the fiducial mass, and the weight fraction of free protons.
The original Abstract Experiment Definition Language (AEDL) file for reactor antineutrino experiments in GLoBES [30] is modified for our purpose. The parametrization of reactor antineutrino flux is taken from Ref. [31], in which the same fuel composition as in Ref. [32] is assumed, and the cross section for the inverse beta decayν e + p → e + + n from Ref. [33] is adopted. The neutrino events are simulated by using the following true values for the oscillation parameters together with the NSI parameters. The running time is taken to be five years, equivalent to the JUNO setup in Ref. [7]. The total number of non-oscillated neutrino events is around 3.9 × 10 5 , and that of the oscillated ones is about 1.2 × 10 5 for the central values of standard mixing parameters in Eq. (18). Note that no backgrounds are assumed in our simulations.
The simulated data are then processed using the standard χ 2 analysis through the following where N i denotes the number of events in the i-th energy bin, the parameter vector p contains the standard oscillation parameters and the systematical errors, and ε represents the NSI parameters. In the fit, all the standard oscillation parameters for theν e →ν e channel are marginalized over, but the NSI parameters are fixed to zero. The reason is that we are interested in how the standard oscillation fit is modified when there are NSIs involved. The prior terms implement external input from other experiments and have the form (p − p true ) 2 /σ 2 p with σ p being the corresponding externally given uncertainty.

A. Energy distribution
First, in Fig. 2 distinction between the black (thick) and red (thin) curves. However, such a spectrum shift can be absorbed into the mixing angle θ 12 , which does not lead to significant contributions to the χ 2 function but may change the fitted value of sin 2 θ 12 dramatically. When ε eµ is positive, the event rate decreases compared to the standard case, implying that the extracted value of θ 12 is slightly larger than its true value. Similarly, a nonvanishing ε eτ would result in a smaller fitted value for θ 12 . In both cases, the amplitude of high-frequency oscillations is enhanced, so the a larger fitted value of θ 13 is expected. These observations are in accordance with our discussions on the survival probability in Section II, and in agreement with our numerical fit in the following subsection. The impact of NSIs on the event rates and measurements of standard oscillation parameters has also been discussed in Ref. [19], where a different experimental setup was considered.

B. NSI effects in the parameter fit
Second, in order to examine the NSI-induced offsets in the standard oscillation fit, we use the true parameters from Eq. (18)   the neutrino events. Then, in the fit, all the standard oscillation parameters are marginalized over, while the NSI parameters are fixed to zero. Assuming two degrees of freedom, we show the effects of NSIs in the fit in Figs. 3 and 4 for NH and some benchmark NSI parameters.
Using Fig. 3, one finds that the deviation of the best-fit sin 2 θ 12 from its true value is remarkable. Even for a relatively small values of the NSI parameters, the true θ 12 may be ruled out erroneously. Consistent with our analytical studies in the previous section, a positive ε eµ leads to a fake θ 12 larger than its true value, whereas a positive ε eτ signifies a positive ε φ , and hence, a smaller fitted θ 12 is favored.
As for θ 13 , one observes from Fig. 4 that the best-fit value only mildly deviates from its true value, which can be viewed as a result of the moderate sensitivity for the JUNO experiment to θ 13 . The result for ε eτ = 0.02 is quite similar to that for ε eµ = 0.02 in the left panel. For comparison, we show in the right panel that the best-fit value of θ 13 becomes smaller than the true value if ε eµ = −0.02. As mentioned before, the determination of the neutrino mass-squared differences are essentially not spoiled by NSI effects, and both ∆m 2 21 and ∆m 2 31 can be fixed at a very-high confidence level as shown in the plots. It is worthwhile to note that the precisions of ∆m 2 21 and ∆m 2 31 are slightly better for larger values of θ 12 and θ 13 , respectively, as observed from Figs. 3 and 4.

C. Discovery reach
Third, we continue to discuss the experimental prospect of detecting NSI effects. For this purpose, we define the discovery reach as the ε ranges where the quality of a standard oscillation fit is below a given confidence level [17]. The numerical results are shown in Figs. 5 and 6. In the left panel of Fig. 5, only the NSI parameters ε eµ and φ eµ are switched on, while in the right panel only ε eτ and φ eτ . One can observe that the best sensitivity to ε eµ or ε eτ appears at φ ∼ 0 or π. This is in agreement with Eqs. (13) and (14), in which the maximal values of |ε φ | and |ε δ | are reached for cos φ eµ = 1 or cos φ eτ = 1. It is worth mentioning that a nonzero value of the Dirac CP-violating phase δ results in only a global shift of the contour lines in Fig. 5. Current bounds on ε eµ and ε eτ at the 90% C.L., i.e., ε eµ < 0.025 and ε eτ < 0.041, are already very stringent [20]. The medium-baseline reactor experiments can improve the bound on ε eτ to below 0.03, but they can hardly set any constraints on ε eµ .
In Fig. 6, we switch on both ε eµ and ε eτ , but set the corresponding phases to zero. In this case, the input value of δ changes the sensitivity dramatically. As shown in the right panel of Fig. 6, a maximal CP-violating phase δ = π/2 leads to a vanishing value of ε δ and a suppressed value of ε φ due to the approximate µ-τ symmetry. Compared to the current bounds on the NSI parameters, a medium-baseline reactor experiment only moderately improves the constraints for δ = 0. In this case, a 3σ hint can be obtained if both ε eµ and ε eτ are larger than 0.02. We remark that the discovery reach on NSIs depends strongly on the assumed uncertainty of θ 13 . In our numerical analysis, a relatively optimistic uncertainty on θ 13 is assumed, which actually relies on the precision of the measurement of θ 13 in ongoing and future non-reactor experiments (e.g., the long-baseline accelerator experiments T2K and NOνA). In general, NSIs induce conflicts between measurement obtained in different type of experiments, and therefore, a combined analysis of both short-, medium-baseline reactor experiments, and accelerator experiments will be more powerful in constraining the NSI parameters. However, such an analysis is beyond the scope of this work and will be performed elsewhere.

D. Impact on the neutrino mass hierarchy determination
Fourth, the NSI parameters modify the oscillation probability, and hence may diminish or improve the experimental sensitivity to the neutrino mass hierarchy. In order to see the NSI effects in the neutrino mass hierarchy fit, we simulate the data in the NH case with nonzero NSI parameters, and then perform a standard oscillation fit with either NH or IH.
The dependence of the minimal χ 2 on the NSI parameters is shown in Fig. 7, where one can observe that the χ 2 can be slightly diminished compared to the standard case for specific choices of the NSI parameters. As in the upper-left plot, the χ 2 value is reduced approximately from 20 to 16 for ε eµ ≃ −0.025 and δ = 0, and the difference between the NH and IH fit decreases to 12. On the other hand, a positive ε eµ or ε eτ leads to a larger χ 2 in the wrong hierarchy fit. The deterioration effects in the neutrino mass hierarchy fit are mainly due to the mimicking effects on θ 13 . Namely, for some choices of NSI parameters, the probability difference between the NH and IH cases defined in Eq. (17) can be smaller than the standard situation, and thus the χ 2 fit becomes a little worse.
However, if current bounds on the NSI parameters are taken into account, the impact of NSIs on the determination of neutrino mass hierarchy is not important. As shown in the right column of Fig. 7, if δ = π/2 is assumed, the NSI effects will be even smaller.
In a different way, one can generate the neutrino events by assuming IH, and fit the data with either NH or IH. However, we have confirmed that the difference between these two approaches is insignificant.
Compared to the χ 2 analysis of the sensitivity to the neutrino mass hierarchy at JUNO in Ref. [7], our minimal χ 2 for vanishing NSI parameters for the wrong hierarchy fit χ 2 ∼ 20 is slightly larger. This difference can be ascribed to our optimistic treatment of the systematic uncertainties and the backgrounds, and to a larger number of simulated neutrino events.

IV. SUMMARY AND CONCLUSIONS
A reactor antineutrino experiment with a medium-baseline of around 50 km has been proposed to pin down the neutrino mass hierarchy, and to precisely measure leptonic mixing parameters and neutrino mass-squared differences. With the expected high-precision measurements of neutrino parameters in this kind of experiment, one can probe new physics beyond the standard-oscillation paradigm as sub-leading effects in neutrino flavor conversion.
In this Letter, we have investigated the impact of NSIs on the medium-baseline reactor antineutrino experiments. For reactor antineutrinos, only the NSI parameters ε eµ and ε eτ , and the corresponding phases φ eµ and φ eτ , are relevant. First of all, we demonstrate that the true value of the leptonic mixing angle θ 12 can be erroneously ruled out at more than 3σ C.L. if the NSI parameter ε eµ or ε eτ is as large as 0.02. However, the extraction of the leptonic mixing angle θ 13 from the experimental data is rarely affected, since the experimental sensitivity to θ 13 is not as high as the one to θ 12 . Then, we show the discovery reach of NSI effects at the medium-baseline experiment. It turns out that the CP-violating phases play an important role in constraining ε eµ and ε eτ . In the most optimistic case, a 3σ hint for NSI effects can be obtained when both ε eµ and ε eτ are larger than 0.02. Finally, the impact of NSIs on the determination of the neutrino mass hierarchy is considered. In principle, NSIs can diminish or enhance significantly the experimental power in discriminating neutrino mass hierarchies. However, when the current bounds on the NSI parameters are taken into account, such effects become insignificant.
It is worthwhile to emphasize that the precision measurements of neutrino parameters in the ongoing and forthcoming neutrino oscillation experiments provide a good opportunity to probe NSIs and the underlying new physics, which is complementary to the direct searches at collider experiments. A complete analysis of the NSI effects should combine all relevant experiments, and obviously deserves further detailed studies.