RG Running of Multiple Neutrino Mixing Parameters at Oscillation Experiments

Abstract

If the new physics scale is within the energy scale of neutrino oscillation experiments, it may lead to a renormalization group (RG) running effect between the production and detection processes as well as between different experiments. It is then possible to use multiple neutrino oscillation experiments to disentangle the multiple RG running parameters. We investigate this effect in a general model-independent sense for a variety of flavor structures in the context of upcoming experiments DUNE-ND, JUNO-TAO, and FASER$\nu$2 that span a large range in neutrino energies and many different flavor combinations. We find strong sensitivity to the running effects of new physics with combination of these experiments, especially the possibility of addressing the non-trivial degeneracies.

I Introduction

The neutrino sector is one of the most promising places to look for new physics beyond the Standard Model of particle physics [11]. While considerable progress has already been made to measure the oscillation parameters and future progress is expected to measure all six parameters [22, 24], the nature and source of neutrino masses still remain unresolved. From the model building point of view, possible new physics may appear at a wide variety of energy scales. It is often simplest to consider the new physics scale to be much higher than all scales in an experiment, typically the weak scale, in an effective field theory context. However, this may miss important physics, especially if the new physics happens at a lower scale. For example, in a compelling new physics model with non-standard interactions (NSI) [41, 16], the new physics scale could be heavy $>100$ GeV and it could also be much lighter, lying within the energy scale of neutrino experiments [27, 26, 14, 25, 21, 15, 1].

If there is new physics that arises at intermediate scales, its effect will appear as not just effective operators but also renormalization group (RG) running of the effective parameters. Such an effect can significantly change the neutrino oscillation phenomena [12, 13]. The oscillation amplitude,

$$\mathcal{A}_{\beta\alpha}\equiv\sum_{i}U_{\beta i}(Q^{2}_{d})e^{-iLm_{i}^{2}/2E_{\nu}}U^{*}_{\alpha i}(Q^{2}_{p}),$$ (1)

contains three parts that corresponds to the neutrino production ($U^{*}_{\alpha i}$), propagation ($e^{-iLm_{i}^{2}/2E_{\nu}}$), and detection ($U_{\beta i}$) processes. Conventionally, the mixing matrix $U^{*}_{\alpha i}$ at production is treated as the same as its counterpart $U_{\beta i}$ at detection. However, this needs not to be true if new physics appears at a scale below the characteristic scale $Q^{2}_{p}$ for production or $Q^{2}_{d}$ for detection. The two mixing matrices $U(Q^{2}_{p})$ and $U(Q^{2}_{d})$ could be different from each other in the presence of RG running [28, 29].

For example, the neutrinos at accelerator experiments are produced dominantly from $\pi^{\pm}$ or $\mu^{\pm}$ decays whose momentum transfer $\sqrt{Q^{2}_{p}}$ is around 100 MeV but are detected by nuclei scatterings with momentum transfer $\sqrt{Q^{2}_{d}}$ up to several GeV. This then leads to a variation of the effective neutrino mixing matrix between the production and detection [18, 10, 9, 44, 39, 34, 37, 38, 43, 33, 12, 13, 28, 29, 36]. Such momentum transfer mismatch can affect the oscillation probabilities. Especially, the RG running effect on the Dirac CP phase can be tested at not just the long-baseline accelerator experiments such as DUNE [28], but also the short-baseline reactor experiment such as JUNO-TAO [29]. In this paper, we focus on the expected future measurements from the near detector (DUNE-ND) at DUNE [32], JUNO-TAO the short-baseline detector apart of JUNO [3], and FASER$\nu$2 at the LHC [6]. Comparing with long-baseline experiments, the near detectors with short baseline can uniquely probe the RG running parameters without complication from the oscillation phases.

In a concrete model, it is unlikely that the running affects only a single oscillation parameter. A more natural case is that multiple parameters can experience RG running simultaneously. In this paper, we consider not just the Dirac CP phase but also the three mixing angles. Especially, we show the correlation among the RG running parameters and possible ways of breaking the degeneracy with multiple experiments. In Sec. II, we first review the phenomenological consequences of RG running on the neutrino oscillations. We then discuss the important role of considering multiple experiments at multiple different energies with different flavor channels accessible to address degeneracies in Sec. III. We present our numerical results in Sec. IV and then discuss them and conclude in Sec. V.

II Neutrino Oscillation with Multiple RG Running Parameters

Quantum radiative corrections can typically introduce energy dependence which can be parameterized as RG running of the neutrino mixing parameters. To be more concrete, the running of the neutrino oscillation parameters are described by the $\beta$-functions $\beta_{X}$ and the renormalization equation [9, 39],

$$\frac{dX}{d\ln|Q^{2}|}=\beta_{X},$$ (2)

where $X$ denotes the neutrino mixing parameters including the three mixing angles $\theta_{ij}$ with $\beta_{ij}$ and the Dirac CP phase $\delta_{D}$ with $\beta_{\delta}$. Without loss of generality, we take the Lorentz-invariant momentum transfer $|Q^{2}|$ as the renormalization scale, which is known as the Gell-Mann-Low scheme [30, 42] and widely used in the literature [17, 12, 13, 28, 29].

The $\beta_{X}$ term controls the running of $X$ and is possibly a function of the mixing parameters themselves [9, 39]. Below the corresponding new physics scale $Q^{2}_{0}$, $\beta(Q^{2}<Q^{2}_{0})=0$ and there is no running. The RG running effect appears above $Q^{2}_{0}$ when $\beta_{X}$ becomes nonzero. For concreteness, the value of $Q^{2}_{0}$ can originate from a hidden light mediator mass. In this case we generally expect $Q^{2}_{0}\gtrsim 1$ MeV² to avoid cosmological bounds [40], which means that there is quite a bit of interesting parameter space probed by neutrino experiments which are the MeV, GeV, and TeV scales.

For $Q^{2}\geq Q^{2}_{0}$, the RG running effect of mixing parameters is typically small [29] and the induced variation of $\beta_{X}$ does not exceed 10% as a rough estimation. Thus, $\beta_{X}$ can be regarded as a constant and perturbative variable ($\beta_{X}\lesssim\mathcal{O}(10^{-1})$) for typical oscillation experiments. In the small $\beta_{X}$ regime and $Q^{2}\geq Q^{2}_{0}$, the problem becomes model-independent with the solution of Eq. (2) obtained by a linear approximation,

$$X(Q^{2})\approx X(Q^{2}_{0})+\beta_{X}\ln\left(\left|\frac{Q^{2}}{Q^{2}_{0}}\right|\right).$$ (3)

While the exact values of the new physics scale $Q^{2}_{0}$ and $\beta_{X}$ depend on the concrete model [39, 38, 12] in the linear regime, $Q^{2}_{0}$ and $\beta_{X}$ can be treated as phenomenological parameters to be constrained through experimental data with $\beta_{X}=0$ returning to the standard physics picture.

The oscillation probability from the $\alpha$ flavor to the $\beta$ flavor can be obtained as $P_{\alpha\beta}(L)\equiv|\mathcal{A}_{\beta\alpha}|^{2}$. Besides the neutrino mixing parameters with the RG running effect, the experimental setup characterized by the neutrino energy $E_{\nu}$ and the baseline $L$ can also affect the oscillation probability via the evolution phase $e^{-iLm^{2}_{i}/2E_{\nu}}$ as shown in Eq. (1). A prominent feature of quantum interference is that it is the phase differences $e^{-iL\Delta m^{2}_{ij}/2E_{\nu}}$ with the mass squared difference $\Delta m^{2}_{ij}\equiv m^{2}_{i}-m^{2}_{j}$ that finally determine the interference probability.

However, it is useful to look to environments where the phase accumulation does not affect the observables as this removes much of the complications arising from the standard three-flavor oscillation parameters not all of which are yet determined. Notably, the evolution phases can reduce to 1 in the zero-distance limit of $L\ll 2E_{\nu}/\Delta m^{2}_{ij}$ [12, 28, 29]. This limit can be achieved in the typical short-baseline neutrino experiments and at the near detectors of long-baseline experiments. The oscillation probability then reduces to,

$$P_{\alpha\beta}(L=0)=\left|\left[U_{d}U_{p}^{\dagger}\right]_{\beta\alpha}\right|^{2}.$$ (4)

The mismatch in the momentum transfers of the neutrino production and detection processes will manifest itself as mismatch, $\delta U\equiv U_{d}-U_{p}$, between the two mixing matrices. So Eq. (4) reduces to,

$$P_{\alpha\beta}(L=0)=\left|\left[\mathbb{1}+\delta UU_{p}^{\dagger}\right]_{\beta\alpha}\right|^{2},$$ (5)

where we have used the fact that $U_{p}$ is a unitary matrix, $U_{p}U^{\dagger}_{p}=\mathbb{1}$. As the zero-distance oscillation probability reduces back to the usual case, $P_{\alpha\beta}(L=0)\rightarrow\delta_{\alpha\beta}$ as $\delta U\to 0$, it would naively seem that any deviation could appear with both linear and quadratic terms.

Actually, we see that the linear term vanishes and only the quadratic order can survive in the zero-distance limit. For $\alpha\neq\beta$, the unit matrix $\mathbb{1}$ in Eq. (5) does not contribute and the transition probability,

$$P_{\alpha\beta}(L=0)=\left|\left[\delta UU_{p}^{\dagger}\right]_{\beta\alpha}\right|^{2},$$ (6)

is at the quadratic order only. For $\alpha=\beta$, the transition probability is

$$P_{\alpha\alpha}(L=0)=1+2\real[(\delta UU^{\dagger}_{p})_{\alpha\alpha}]+|(\delta UU^{\dagger}_{p})_{\alpha\alpha}|^{2},$$ (7)

with terms containing both the linear and quadratic orders of the deviations $\delta U$ induced by RG running. However, the unitarity condition, $U_{d}U_{d}^{\dagger}=\mathbb{1}$, implies that $2{\rm Re}[(\delta UU_{p}^{\dagger})_{\alpha\alpha}]=-(\delta U\delta U^{\dagger})_{\alpha\alpha}$. In other words, the linear term cancels and the correction is always quadratic.

Putting things together, the transition probability in the zero-distance limit takes the form as,

$$P_{\alpha\beta}(L=0)=\left[1-(\delta U\delta U^{\dagger})_{\alpha\alpha}\right]\delta_{\alpha\beta}+[|\delta UU^{\dagger}_{p}|_{\beta\alpha}]^{2},$$ (8)

from which one can clearly see that the deviation from the usual $\delta_{\alpha\beta}$ is always at the quadratic order.

Consequently, in the linear regime of Eq. (3), the zero-distance oscillation probability is quadratic in $\ln^{2}\left(|Q^{2}_{d}|/|Q^{2}_{p}|\right)$ and its scale dependence can be generally parametrized in terms of effective parameters $\mathcal{B}_{\alpha\beta}$,

$$P_{\alpha\beta}(L=0)\equiv\delta_{\alpha\beta}+\mathcal{B}_{\alpha\beta}\ln^{2}\left(\frac{|Q^{2}_{d}|}{|Q^{2}_{p}|}\right),$$ (9)

with $\delta_{\alpha\beta}$ being the Kronecker delta. The effective running parameters $\mathcal{B}_{\alpha\beta}$ for the oscillation probability $P_{\alpha\beta}(L=0)$ are functions of the mixing angles and CP phases calculated at production. Note that the running of $P_{\alpha\beta}(L=0)$ is governed by a double logarithm rather than the single one in Eq. (3).

While $P_{\alpha\beta}(L=0)$ has only quadratic dependence, the long-baseline oscillation can have linear dependence on the RG running parameters $\beta_{X}$ [12]. For a perturbative $\beta_{X}$, the linear dependence would have larger effect than its quadratic counterpart. However, longer baseline means significant dilution of the neutrino fluxes at detector. It is actually of great advantages to use the short-baseline experiment for probing the RG running parameters. For simplicity, we would omit $L=0$ for the zero-distance limit and $P_{\alpha\beta}$ as abbreviation of $P_{\alpha\beta}(L=0)$ in the later part of this paper.

As we have already noticed, there are four degrees of freedom in the mixing matrix relevant for oscillation experiments. We consider the usual parameterization with three mixing angles $\theta_{ij}$ and the Dirac CP phase $\delta_{D}$, see [23] for a discussion of the choice of parameterization. One can map the new physics into RG running of the four parameters, each with their own $\beta_{X}$. Once taking multiple RG running parameters into account, degeneracies can exist among them in the transition probabilities. For simplicity, we consider only two RG running parameters at a time in the following discussions which we find to be sufficient to understand the parameter space.

For a typical short-baseline experiment like the DUNE near detector (DUNE-ND), the relevant channels are $\nu_{\mu}\rightarrow\nu_{e}$ and $\nu_{\mu}\rightarrow\nu_{\mu}$ as well as their anti-neutrino counterparts. In the following, unless otherwise specified, we focus on the discussion about the neutrino mode, while the discussion about the anti-neutrino mode can be derived in a similar way. When expanded to the quadratic order of $\beta_{X}$, the transition probability effective parameters that are involved in the DUNE-ND are,

$$\mathcal{B}_{\mu e} =\beta_{12}^{2}c_{13}^{2}c_{23}^{2}+\beta_{12}\beta_{13}c_{\delta}c_{13}s_{2\theta_{23}}+\beta_{13}^{2}s^{2}_{23}$$

$$+\beta_{\delta}^{2}c_{13}^{2}s_{13}^{2}s_{23}^{2}-\beta_{\delta}\beta_{12}c_{13}^{2}s_{\delta}s_{13}s_{2\theta_{23}}+\mathcal{O}(\beta_{ij}^{3})$$ (10a)

$$\mathcal{B}_{\mu\mu} =-\beta_{12}^{2}\left(1-s_{23}^{2}c_{13}^{2}-s_{13}^{2}s^{2}_{\delta}s^{2}_{2\theta_{23}}\right)-\beta_{13}^{2}s^{2}_{23}$$

$$\quad-\beta_{23}^{2}-\beta_{\delta}^{2}s^{2}_{13}s^{2}_{23}\left(1-s^{2}_{13}s^{2}_{23}\right)$$

$$\quad-\beta_{12}\beta_{13}c_{\delta}c_{13}c_{2\theta_{23}}-2\beta_{12}\beta_{23}c_{\delta}s_{13}$$

$$\quad+\beta_{\delta}\beta_{12}s_{\delta}s_{13}s_{2\theta_{23}}\left(1-2s^{2}_{13}s^{2}_{23}\right)+\mathcal{O}(\beta_{ij}^{3})\,,$$ (10b)

where we denote $(c_{ij},s_{ij})\equiv(\cos\theta_{ij},\sin\theta_{ij})$ and $(c_{\delta},s_{\delta})\equiv(\cos\delta_{D},\sin\delta_{D})$. Here we only keep the RG running parameters up to the second order, while the higher-order terms are neglected. As discussed above, the presence of the RG running effect can induce non-trivial transition probabilities. From the above formula, different transition channels are determined by different RG running parameters. For example, there is no $\beta_{23}$ term in the expression of $\mathcal{B}_{\mu e}$ while all four $\beta_{ij}$ and $\beta_{\delta}$ are involved in $\mathcal{B}_{\mu\mu}$.

Nevertheless, the oscillation channels $\nu_{\mu}\rightarrow\nu_{e}$ and $\nu_{\mu}\rightarrow\nu_{\mu}$ involves all the four RG running beta functions $\beta_{ij}$ and $\beta_{\delta}$. Using only these two channels is not enough to break the degeneracy. It is then necessary to search for other transition channels to break the degeneracy.

However, not all oscillation probabilities are independent. Although the unitarity property for a single transition probability has been lost [29], the sum of transition probabilities still preserve some unitarity properties,

$$\sum_{\alpha}P_{\alpha\beta}(E_{\nu};Q^{2}_{p},Q^{2}_{d})=\sum_{\beta}P_{\alpha\beta}(E_{\nu};Q^{2}_{p},Q^{2}_{d})=1.$$ (11)

In the zero-distance limit, the transition probabilities only involve the mixing matrices as shown in Eq. (4) and the summation above can then be expressed as,

$$\sum_{\alpha}P_{\alpha\beta}=\left[U_{p}U^{\dagger}_{d}U_{d}U^{\dagger}_{p}\right]_{\beta\beta}=\delta_{\beta\beta}.$$ (12)

Even with nonzero baseline, the evolution phase sandwiched between the two mixing matrices would not change the final conclusion. This unitarity relation shows that among the total nine transition channels, only four transition channels are independent.

Taking the unitarity condition in Eq. (11) for the probability sum into account, only 4 of the 9 oscillation probabilities could be independent. Besides the aforementioned $\nu_{\mu}\rightarrow\nu_{e}$ and $\nu_{\mu}\rightarrow\nu_{\mu}$ at the near detector of the DUNE accelerator experiment, one may also consider the JUNO-TAO near detector at the reactor experiment whose transition probability is determined by

$$\mathcal{B}_{ee}=-\beta_{13}^{2}-\beta_{12}^{2}c^{2}_{13}-\beta_{\delta}^{2}c_{13}^{2}s_{13}^{2}.$$ (13)

In addition, we take the FASER$\nu$2 experiment that consists of $\mathcal{O}$(TeV) neutrinos into account. Due to its high-energy and multiple neutrino sources with three flavors produced, FASER$\nu$2 can detect all three flavor neutrinos. Particularly interesting is the $P_{e\tau}$ transition probability with,

$$\mathcal{B}_{e\tau} =\beta_{12}^{2}c_{13}^{2}s_{23}^{2}-\beta_{12}\beta_{13}c_{\delta}c_{13}s_{2\theta_{23}}+\beta_{13}^{2}c_{23}^{2}.$$ (14)

Putting things together, we may see that $\mathcal{B}_{\mu e}$ in Eq. (10a), $\mathcal{B}_{ee}$ in Eq. (13), and $\mathcal{B}_{e\tau}$ above share the same set of beta functions ($\beta_{12}$, $\beta_{13}$, and $\beta_{\delta}$). Such combination of three independent oscillation channels can really constrain the three beta functions. Then the remaining $\mathcal{B}_{\mu\mu}$ can constrain the fourth beta function $\beta_{23}$.

As multiple RG running parameters are simultaneously involved, certain degeneracy among these RG running parameters exist in any single channel. Explicitly, $\mathcal{B}_{\mu e}$ has a cross term of $\beta_{12}\beta_{13}$ which indicates a correlation between them. Notice that the correlation depends on the sign of $c_{\delta}$. Namely, since both prefactors of $\beta_{12}^{2}$ and $\beta_{13}^{2}$ are positive in $\mathcal{B}_{\mu e}$, there exists a positive (negative) correlation between $\beta_{12}$ and $\beta_{13}$ for $c_{\delta}<0$ ($c_{\delta}>0$). In the case that $c_{\delta}=0$, there is no correlation between $\beta_{12}$ and $\beta_{13}$. Similarly, although the cross term of $\beta_{12}$ and $\beta_{13}$ in $\mathcal{B}_{\mu\mu}$ has an opposite sign to its counterpart in $\mathcal{B}_{\mu e}$, the prefactors of $\beta_{12}^{2}$ and $\beta_{13}^{2}$ are also opposite. Consequently, the correlation behavior between $\beta_{12}$ and $\beta_{13}$ in $\mathcal{B}_{\mu\mu}$ is the same as the one in $\mathcal{B}_{\mu e}$. Below we will show with more details how the three experiments (DUNE-ND, JUNO-TAO, and FASER$\nu$2) can break the degeneracy among the RG running beta functions.

Figure 1: The momentum transfer averaged transition probabilities $\overline{P}^{0}_{\mu\mu}=0.99$ (blue solid) and $\overline{P}^{0}_{\mu e}=0.01$ (red dashed) on the $\beta_{12}$–$\beta_{13}$ space. The neutrino energy $E_{\nu}$ is set as 2.5 GeV and the new physical scale $Q^{2}_{0}$ is set as 1  MeV². The production momentum transfer is fixed at $m_{\pi}^{2}$. We take the global-fit result of $\delta_{D}=194^{\circ}$ from [20].

III Breaking Degeneracies with Low- and High-energy SBL Experiments

As shown above, degeneracies exist among multiple RG running parameters especially when considering only a single experiment at a time. This will lead to a correlation between $\beta_{12}$ and $\beta_{13}$ in both detection channels, namely, the $\nu_{\mu}$ disappearance channel and $\nu_{e}$ appearance channel. Such a degeneracy can significantly reduce the experimental sensitivity to the multiple RG running parameters when only one experiment is considered at a time. To numerically calculate the transition probability and show the correlation behavior for a realistic experiment, we take a momentum transfer averaged transition probability,

$$\overline{P^{0}_{\alpha\beta}}\equiv\int_{Q^{2}_{d,\rm min}}^{Q^{2}_{d,\rm max}}P^{0}_{\alpha\beta}(Q^{2}_{d})\left[\frac{1}{\sigma}\frac{d\sigma}{dQ^{2}_{d}}\right]dQ^{2}_{d}\,.$$ (15)

where $d\sigma/dQ^{2}_{d}$ is the differential detection cross section and $\sigma$ is the total cross section. We use this approach since in most neutrino experiments the detection momentum transfer is not fixed, but follows a distribution determined by the neutrino energy and detection cross-section. However, neutrinos are produced from particle decays via the charged current weak interaction and hence the production momentum transfer is constant and equal to the parent particle mass squared.

Fig. 1 shows the averaged probabilities $\overline{P^{0}_{\mu\mu}}$ (blue) and $\overline{P^{0}_{\mu e}}$ (red) that give a 1% correction to the standard picture, on the $\beta_{12}$–$\beta_{13}$ space. For illustration, we take the neutrino peak energy $E_{\nu}=2.5$ GeV of the DUNE experiment and $Q^{2}_{p}=m_{\pi}^{2}$ [2]. We also take the global-fit result $\delta_{D}=194^{\circ}$ [20] and the detection momentum distribution was obtained with GENIE [7]. In both curves, there is a positive correlation between $\beta_{12}$ and $\beta_{13}$, which is consistent with the analytical formula of Eq. (10). The blue curve has an ellipse-like shape with the major axis directed towards the 1st quadrant direction. The red curve matches the blue curve for $|\beta_{ij}|\lesssim 0.05$. While for $|\beta_{ij}|>0.05$ the red curve has an elongated tail extending above $|\beta_{ij}|>0.2$ which comes from the $\beta_{X}^{4}$ terms omitted in Eq. (10a) and Eq. (10b). The existence of the positive correlation generates the degeneracy between the RG running parameters which can significantly reduce the experimental sensitivity.

Figure 2: The averaged transition probability $\overline{P}^{0}_{ee}=0.9$ at JUNO (orange solid) and $\overline{P}^{0}_{e\tau}=0.01$ at FASER$\nu$2 (green dashed) on the $\beta_{12}$–$\beta_{13}$ space. The neutrino energy $E_{\nu}$ is set as 5 MeV and the new physical scale $Q^{2}_{0}$ is set as 1 MeV². The production momentum transfer is fixed as $(m_{n}-m_{p})^{2}\approx 1.3\,$MeV².

In general, there are two possible ways to look for complementary experiments involved with different transition channels from DUNE-ND. One is using different neutrino sources. As mentioned above, the major neutrino flux in the production at DUNE-ND is muon flavor. Therefore, the possible different neutrino sources are electron and tau flavors. However, due to experimental challenges, it is difficult to produce a sufficient tau-flavor neutrino flux to have a large statistics. Therefore, it is promising to focus on the electron-flavor flux.

Indeed, the $\bar{\nu}_{e}$ flux can be sufficiently produced from a reactor and can be detected by the near detector in a typical reactor experiment. In this paper, we consider JUNO-TAO [3] as a benchmark experiment. Especially, the corresponding transition probability for the reactor neutrino experiment has been summarized in Eq. (13). Being different from Eq. (10), there is no cross term among the RG running parameters. Therefore, no correlation between $\beta_{12}$ and $\beta_{13}$ appears in the electron disappearance channel.

The orange contour in Fig. 2 on the $\beta_{12}$–$\beta_{13}$ space shows the transition probability $\bar{P}_{ee}=0.9$ with typical neutrino energy $E_{\nu}=5\,$MeV and production momentum transfer $Q^{2}_{p}=1.3\,$MeV² at JUNO-TAO [29]. The contour is a circle, indicating that there is no correlation between $\beta_{12}$ and $\beta_{13}$ which is consistent with Eq. (13). This feature provides JUNO-TAO as a powerful platform to constrain on $\beta_{12}$ and $\beta_{13}$ without any contamination from the degeneracy.

Besides the different neutrino flavor source from production, the other practical way for breaking degeneracies is expanding the detection possibility. Typically, the detection channels are limited by the energy threshold. Explicitly, the JUNO-TAO can only detect electron flavor since its $\mathcal{O}$(MeV) neutrino energies are below the muon or tau mass. Similarly, with $\mathcal{O}$(GeV) energies, DUNE-ND can only efficiently detect electron and muon flavors, although a small number of tau events can be measured with high uncertainties. Therefore, in order to expand the detection possibility, we need to consider other higher-energy experiments.

As illustrated in Fig. 2, we have also shown the transition probability $\overline{P^{0}_{e\tau}}=0.01$ for the electron flavor transition into the tau flavor. Here we take $500\,$GeV neutrino energy and $m_{\pi}^{2}$ production momentum transfer as representative values. Unlike the positive correlation observed in Fig. 1, $\overline{P^{0}_{e\tau}}$ has a negative correlation. The combination of measurements with oppositely signed correlations allows for the removal of the degeneracy. More intuitively, although the overlapping area of the contours in Fig. 2 and Fig. 1 remains, it is much smaller when compared with each region separately.

To further illustrate the transition probabilities in the three experiments mentioned above, we define the average of the logarithm term $\overline{\ln^{2}|Q^{2}_{d}|/|Q^{2}_{p}|}$ to obtain the detection-averaged probability,

$$\overline{P^{0}_{\alpha\beta}}\equiv\delta_{\alpha\beta}+\mathcal{B}_{\alpha\beta}~~\overline{\ln^{2}\left(\frac{|Q^{2}_{d}|}{|Q^{2}_{p}|}\right)}.$$ (16)

The transition channels and values of the average logarithm squared for JUNO-TAO [3], DUNE-ND[32], and FASER$\nu$2 [35] are summarized in Table  1. The $\ln^{2}$ terms are $\mathcal{O}(10)$, meaning that a $\beta_{X}\sim\mathcal{O}(0.1)$ is close to the non-linear regime. In the non-linear regime, large changes to oscillation parameters take place and the oscillation data needs to be carefully re-analyzed. Moreover, due to the correlations described above, constraining the $\beta_{X}$ parameters within the linear regime is a difficult task for any experiment alone. Fortunately, as we will show, the combination of the three future experiments JUNO-TAO, DUNE-ND and FASER$\nu$2 is ideal to bring down the bounds on the RG running of neutrino mixing parameters to the linear regime.

Exp.	Oscillation Channels	$\overline{\ln^{2}|Q^{2}_{d}|/|Q^{2}_{p}|}$
JUNO-TAO	$\overline{\nu}_{e}\rightarrow\overline{\nu}_{e}$	5 ($E_{\nu}=3.6$ MeV)
DUNE-ND	$\nu_{\mu}\rightarrow\nu_{e,\mu}$  &  $\overline{\nu}_{\mu}\rightarrow\overline{\nu}_{e,\mu}$ $\nu_{e}\rightarrow\nu_{e,\mu}$  &  $\overline{\nu}_{e}\rightarrow\overline{\nu}_{e,\mu}$	12 ($E_{\nu}=2.5$ GeV)
FASER$\nu$2	$\nu_{\mu}\rightarrow\nu_{e,\mu,\tau}$  &  $\overline{\nu}_{\mu}\rightarrow\overline{\nu}_{e,\mu,\tau}$ $\nu_{e}\rightarrow\nu_{e,\mu,\tau}$  &  $\overline{\nu}_{e}\rightarrow\overline{\nu}_{e,\mu,\tau}$ $\nu_{\tau}\rightarrow\nu_{e,\mu,\tau}$  &  $\overline{\nu}_{\tau}\rightarrow\overline{\nu}_{e,\mu,\tau}$	82 ($E_{\nu}=0.5$ TeV)

Table 1: Oscillation channels at JUNO-TAO, DUNE-ND, and FASER$\nu$2 and the detection-averaged squared logarithm for a representative neutrino energy $E_{\nu}$.

IV Sensitivities to RG Running Parameters

To optimize the sensitivity of the RG parameters, we combine three SBL experiments: JUNO-TAO, DUNE-ND, and FASER$\nu$2, which cover a large energy range from MeV up to TeV and provide independent measurement of different transition channels: $P^{0}_{ee}$, $P^{0}_{\mu e}$, and $P^{0}_{e\tau}$, respectively. As shown above, these transition channels have totally different correlation behavior. For completeness, the transition channels of these three experiments are summarized in Table  1.

As the neutrino production and detection processes are different among these three experiments, we need different considerations to calculate the predicted event rate. Since JUNO-TAO is already taking data [4] and its systematics is established [3], we use the complete event distribution over the reconstructed energy and include the systematic in our simulations. For the future DUNE-ND and FASER$\nu$2 detectors, we take a more conservative approach and consider only the total number of events expected to be observed in each experiment. Even with such crude observable, the sensitivity obtained reaches competitive levels, especially when combining different experiments.

For the JUNO-TAO experiment, the initial neutrino flux contains only electron antineutrinos. Also, because the neutrino energies are below the muon or tau mass, only $\bar{\nu}_{e}$ can be measured at the detector. The predicted true neutrino event rate contains only a single channel that can be calculated,

$$\frac{\differential N^{\rm JUNO-TAO}_{\bar{\nu}_{e}\rightarrow\bar{\nu}_{e}}}{\differential E_{\nu}}=N_{t}\Phi_{\bar{\nu}_{e}}\overline{P^{0}_{ee}}\sigma_{\bar{\nu}_{e}},$$ (17)

where $N_{t}$ is the normalization factor, $\Phi_{\bar{\nu}_{e}}$ the initial $\bar{\nu}_{e}$ flux, which is a combination of several reactors each containing different fractions of isotope. Also, the detection cross section $\sigma_{\bar{\nu}_{e}}$ is taken to be the inverse beta decay cross section. In order to simulate the energy resolution of the detector, we implement a Gaussian smearing function with energy resolution of $1.5\%/\sqrt{E/{\rm MeV}}$ [3] where $E\equiv E_{\bar{\nu}_{e}}-0.8\,$MeV is the deposit energy.

The quantity $\overline{P^{0}_{ee}}$ is the zero-distance momentum transfer averaged transition probability for $\bar{\nu}_{e}\rightarrow\bar{\nu}_{e}$. In order to calculate the transition probability, one needs three momentum transfer scales, $Q^{2}_{0}$, $Q^{2}_{p}$, and $Q^{2}_{d}$. For illustration, we take $Q^{2}_{0}=1\,$MeV² throughout this paper. We also take $Q^{2}_{p}=(m_{n}-m_{p})^{2}$ corresponding to the beta decay momentum transfer and an analytical distribution function for $Q^{2}_{d}$ [31]. The backgrounds come from accidental, fast neutron, and ⁹Li/⁸He backgrounds and their spectra are extracted from [5]. The sensitivity is calculated by a sum over all the bins through the function $\chi^{2}_{\rm JUNO-TAO}$,

$$\hskip-8.5359pt\chi^{2}_{\rm JUNO-TAO} \equiv\sum_{i}^{\rm bins}\left(\frac{N_{i}^{\rm true}-N_{i}^{\rm test}}{\sqrt{N_{i}^{\rm true}}}\right)^{2}+\chi_{\rm JUNO-SYS}^{2}.$$ (18)

The systematic uncertainties in $\chi_{\rm JUNO-SYS}^{2}$ includes several normalization and tilt parameters [3]. For a more detailed description of our JUNO-TAO simulation refer to [29]. To derive the sensitivity, we assume 6.5 years of data-taking time.

For DUNE-ND, we take only the ND-LAr module [32] into account. There are two running modes at DUNE, namely, the neutrino mode and anti-neutrino mode. In each mode, there are four neutrino components in the initial neutrino flux: $\nu_{\mu},\bar{\nu}_{\mu},\nu_{e}$, and $\bar{\nu}_{e}$ [2]. As for neutrino detection, we follow the DUNE near detector CDR [32] and classify them into two categories. For the neutrino mode, they are $\nu_{\mu}$ and $\nu_{e}+\bar{\nu}_{e}$, while for the anti-neutrino mode, they are $\bar{\nu}_{\mu}$ and $\nu_{e}+\bar{\nu}_{e}$. Therefore, the predicted event rates in the neutrino mode are,

$$\frac{\differential N^{\rm DUNE-ND}_{\nu_{\alpha}\rightarrow\nu_{\mu}}}{\differential E_{\nu}} \equiv N_{t}\Phi_{\nu_{\alpha}}\overline{P^{0}_{\alpha\mu}}\sigma_{\nu_{\mu}}\,,$$ (19a)

$$\frac{\differential N^{\rm DUNE-ND}_{\nu_{\alpha}\rightarrow\nu_{e}}}{\differential E_{\nu}} \equiv N_{t}\left[\Phi_{\nu_{\alpha}}\overline{P^{0}_{\alpha e}}\sigma_{\nu_{e}}+\Phi_{\bar{\nu}_{\alpha}}\overline{P^{0}_{\bar{\alpha}\bar{e}}}\sigma_{\bar{\nu}_{e}}\right]\,,$$ (19b)

and the predicted event rates in the anti-neutrino mode can be calculated in a similar way.

Since most DUNE neutrinos are produced from the decay of pions and muons which have a similar mass, we fix the neutrino production momentum transfer at $Q^{2}_{p}=m_{\pi}^{2}$ for simplicity. For the detection momentum transfer, we take the distribution from the GENIE [7] simulation results. To be conservative, we don’t include the energy spectrum information for the DUNE-ND, but only the total number of events for each mode and channel. Moreover, we include an overall nuisance parameter for each flux source when calculating the experimental sensitivity. For example, the muon-neutrino detection channel in the neutrino mode has,

$$\hskip-8.53581pt\chi^{2}_{\nu\text{-mode},\nu_{\mu}}\equiv\left(\frac{\sum_{\alpha}(1+a_{\alpha}^{\nu\text{-mode}})N^{\rm RG}_{\nu_{\alpha}\rightarrow\nu_{\mu}}-N^{\rm std}_{\nu_{\mu}}}{\sqrt{N^{\rm std}_{\nu_{\mu}}}}\right)^{2},$$ (20)

where $N^{\rm RG}_{\nu_{\alpha}\rightarrow\nu_{\mu}}$ is the bin-summed prediction according to Eq. (19) with the RG running effect and $N^{\rm std}_{\nu_{\mu}}$ is the total predicted number of muon neutrino events in the standard oscillation scenario. Here $a_{\alpha}^{\nu\text{-mode}}$ is the overall nuisance parameter of $\alpha$ in the neutrino mode. Summing up two neutrino modes and all relevant channels, the total DUNE-ND $\chi^{2}$ function is,

$$\chi^{2}_{\rm DUNE\!-\!ND}=\sum_{\begin{subarray}{c}\rm mode,\end{subarray}\beta}\chi^{2}_{\rm mode,\beta}+\sum_{\begin{subarray}{c}\rm mode,\end{subarray}\alpha}\left(\frac{a^{\rm mode}_{\alpha}}{\sigma_{a^{\rm mode}_{\alpha}}}\right)^{2}\,.$$ (21)

Here $\beta\in\{\nu_{\mu},\nu_{e}+\bar{\nu}_{e}\}/\{\bar{\nu}_{\mu},\nu_{e}+\bar{\nu}_{e}\}$ denotes the observable in the neutrino/anti-neutrino mode and $\alpha\in\{\nu_{\mu},\bar{\nu}_{\mu},\nu_{e},\bar{\nu}_{e}\}$ denotes the neutrino source. Each nuisance parameter $a_{\alpha}^{\rm mode}$ is minimized over with a Gaussian prior with uncertainty $\sigma_{\alpha}^{\rm mode}=10\%$.

For FASER$\nu$2, the neutrinos in both production and detection processes can be classified into three categories: $\nu_{e}+\bar{\nu}_{e}$, $\nu_{\mu}+\bar{\nu}_{\mu}$, and $\nu_{\tau}+\bar{\nu}_{\tau}$ [35]. Compared with the above two experiments, the neutrino production process is more complicated at the FASER$\nu$2 experiment. Explicitly, neutrinos are produced from the decay of multiple particles which have different masses with comparable contributions. Recall that the neutrino production momentum transfer is equal to the mass squared of the parent particle. Consequently, the calculation of event rates at FASER$\nu$2 should be divided into contributions from various parent particles and then summed up,

$$\frac{\differential N^{{\rm FASER}\nu 2}_{\nu_{\alpha}\rightarrow\nu_{\beta}}}{\differential E_{\nu}} =N_{t}\sum_{\rm parents}\left[f_{\rm p}\Phi_{\nu_{\alpha}}\overline{P^{0,~p}_{\alpha\beta}}\sigma_{\nu_{\beta}}\right.$$

$$\qquad\qquad\qquad\left.+f_{\bar{\rm p}}\Phi_{\bar{\nu}_{\alpha}}\overline{P^{0,~p}_{\bar{\alpha}\bar{\beta}}}\sigma_{\bar{\nu}_{\beta}}\right],$$ (22)

where $f_{\rm p}$ is the fraction of a given parent particle among all parents and $\overline{P^{0,~p}_{\bar{\alpha}\bar{\beta}}}$ is the momentum transfer averaged probability calculated with the production momentum transfer at the parent particle mass, $Q^{2}_{p}=m_{p}^{2}$. The main neutrino parent particles in the electron/muon/tau flavor production are (kaon, Dmeson, hyperon)/(pion, kaon, Dmeson)/(Dmeson) [35], respectively. Moreover, the fraction of each source depends on the neutrino energies [35]. For the detection momentum transfer, we still take its distribution from the GENIE [7] simulation.

The $\chi^{2}$ for the detection of $\nu_{\beta}$ at FASER$\nu$2 is constructed similarly to the DUNE-ND, by summing over all the energy bins and introducing an overall nuisance parameter for each neutrino source $\alpha$,

$$\chi^{2}_{\beta} \equiv\left(\frac{\sum_{\alpha}(1+a_{\alpha})N^{\rm RG}_{\alpha\beta}-N^{\rm std}_{\beta}}{\sqrt{N^{\rm std}_{\beta}}}\right)^{2}.$$ (23)

Similarly, $N^{\rm RG}_{\alpha\beta}$ ($N^{\rm std}_{\nu_{\mu}}$ ) is the bin-summed prediction with RG running effect (standard oscillation) using Eq. (22). Summing up the relevant channels, the total $\chi^{2}$ for FASER$\nu$2 is,

$$\chi^{2}_{\rm FASER\nu 2}=\sum_{\begin{subarray}{c}\beta\end{subarray}}\chi^{2}_{\beta}+\sum_{\alpha}\left(\frac{a_{\alpha}}{\sigma_{a_{\alpha}}}\right)^{2}.$$ (24)

There are three nuisance parameters $a_{\alpha}$, $\alpha=e,\mu$ and $\tau$ corresponding to a Gaussian prior with uncertainty $\sigma_{e}=\sigma_{\mu}=50\%$, and $\sigma_{\tau}=100\%$, which is more conservative than the values used in the literature [8, 19].

The predicted event rates in Eq. (17), Eq. (19), and Eq. (22) give the event distribution as a function of neutrino energy. In the presence of the RG running effect, a non-zero transition probability will lead to a different event rate from the standard case even at the total number of event level. Therefore, the short-baseline neutrino experiments provide a unique platform to test the RG running effect by comparing the measured event rates with the predicted one. The sensitivity should improve with a full energy spectrum analysis and also if the detection momentum transfer can be measured [28]. However, such analysis requires complicated treatment of the systematics, for the bin-to-bin uncertainties, which is beyond the scope of this work. Therefore, to simplify the calculation of the experimental sensitivity, we focus on the total event numbers for DUNE-ND and FASER$\nu$2 and include the energy spectrum information only for JUNO-TAO.

Figure 3: The total event difference fraction of flavor $\nu_{\alpha}$ between the standard ($N_{\alpha,\rm std}$) oscillation and with RG-running effect ($N_{\alpha,\rm RG}$) for JUNO-TAO, DUNE-ND and FASER$\nu$2 experiments. For illustration, we take the only non-zero RG-running parameters to be $\beta_{12}=0.03$ and $\beta_{13}=-0.03$. The error in the $x$-direction is estimated by taking the full energy range of each experiment, while the $y$-direction error is estimated using Eq. (25).

Figure 4: The JUNO-TAO (green), DUNE-ND (red), and FASER$\nu$2 (blue) expected sensitivity to the $\beta_{12}$ and $\beta_{13}$ RG running parameters at 68% C. L. The solid lines represents fixed at the best-fit value of [20] ($\delta_{D}=1.08\pi$) while the dashed lines minimize over the value of $\delta_{D}$ with prior knowledge on its uncertainty. The black curve represents the combination of all the three experiments sensitivities and represents a significant improvement on the sensitivity when compared with the separated cases

For illustrating the power of total event number analysis, in Fig. 3 we show the expected event difference fraction, $\Delta N_{\beta}/N_{\beta}^{\rm std}\equiv(N_{\beta}^{\rm RG}-N_{\beta}^{\rm std})/N_{\beta}^{\rm std}$ of $\beta$ flavor detection for the experiments JUNO, DUNE and FASER$\nu$2. When displaying the uncertainty in the y-axis direction, we combine the statistical and systematic ones. For example, FASER$\nu$2 has $\mathcal{O}(10^{6})$ $\nu_{\mu}$ events, which makes the statistical error tiny $\Delta N^{\rm error}_{\mu}\approx\sqrt{N_{\mu}^{\rm std}}/N_{\mu}^{\rm std}\sim 0.05\%$. On the other side, the systematical error for each flux is of the order of $50\%$ to $100\%$, however, the overall error is reduced due to the correlation between the different measuring channels. In order to take into account some of these features, we make an estimate about the error using the diagonal of a covariance matrix $\Delta N_{\beta}^{\rm error}\approx\sqrt{(\Sigma^{-1})_{\beta\beta}}$ such that

$$\Sigma_{\beta\alpha}\equiv\frac{\delta_{\alpha\beta}}{N_{\beta}^{\rm std}}-\frac{({\bf N}_{\beta}^{\rm RG}/N_{\beta}^{\rm std})({\bf N}_{\alpha}^{\rm RG}/N_{\alpha}^{\rm std})}{\sum_{\gamma}({\bf N}_{\gamma}^{\rm RG})^{2}(N_{\gamma}^{\rm std})^{-1}+\sigma^{-2}}\,,$$ (25)

where $\sigma$ is the normalization error. Notice that the actual error and correlations included in our sensitive analysis is more complicated than the one in the equation above and provides better sensitivity. In any case, the error bars displayed in Fig. 3 provide a good estimate on the expected results. In fact we clearly see that the $\nu_{e}$ channel (red) for DUNE-ND and the $\tau$ channel (green) for FASER$\nu$2 can provide a powerful sensitivity to the RG effects at the total event rate level since the gray-dashed line representing zero difference from the standard case is below the error bars.

Figure 5: The JUNO-TAO (dashed-green), DUNE-ND (red), and FASER$\nu$2 (blue) expected sensitivity to the $\beta_{23}$ and $\beta_{\delta}$ renormalization parameters at 68% C. L. The JUNO-TAO experiment is insensitive to the $\beta_{23}$ parameter, so the region is obtained by the 1D interval obtained by the condition $\Delta\chi^{2}=0.05$, while the DUNE-ND and FASER$\nu$2 are sensitive to both parameters and the curves are obtained with the 2D confidence region obtained from $\Delta\chi^{2}=2.3$. In all cases we fix $\delta_{D}=1.08\pi$. The combined sensitivity (black) represents a significant improvement on the sensitivity when compared with the separated cases.

Our sensitivity analysis provides the separate sensitivity for JUNO-TAO, DUNE-ND, and FASER$\nu$2 according to the $\chi^{2}$ described by Eq. (18), Eq. (21), and Eq. (24), respectively. For the combined analysis, we adopt the total $\chi^{2}$ accordingly,

$$\chi^{2}\equiv\chi^{2}_{\rm TAO}+\chi^{2}_{\rm DUNE-ND}+\chi^{2}_{{\rm FASER}\nu 2}+\chi^{2}_{\rm prior}.$$ (26)

The last term $\chi^{2}_{\rm prior}$ includes a Gaussian prior for the Dirac CP phase $\delta_{\rm D}$ phase whose uncertainty is taken from the global-fit result [20]. All the other neutrino oscillation parameters are taken at their best-fit values [20].

The sensitivity result for the combination ($\beta_{12}$, $\beta_{13}$) is depicted in Fig. 4. As we can see, the JUNO-TAO (green) sensitivity has no correlation and goes from $\sim-0.03$ to $\sim 0.03$ for both parameters. At the same time, the DUNE-ND (red) sensitivity has a positive correlation that drastically reduces the sensitivity of both parameters creating a blind spot that can reach ($\beta_{12}$, $\beta_{13}$) = $\pm$(0.14, 0.14) when the $\delta_{D}$ is minimized over (dashed-red line). The $|\beta_{ij}|\sim\mathcal{O}(0.15)$ in combination with DUNE-ND’s average momentum transfer with $\ln^{2}(|Q_{p}^{2}|/|Q_{d}^{2})|\approx 10$ is essentially outside the linear approximation regime of Eq. (3). The presence of extra transition channels, especially for the muon/electron transition into the tau flavor, at FASER$\nu$2 (blue) generates a negative correlation and reduces the impact to the sensitivity when compared to DUNE-ND. However, the correlations are still important and reduces the FASER$\nu$2 capability to constrain the RG running parameters. Notably, the combined result (black curve) takes advantage of the different correlations among all experiments and greatly improves the sensitivity such that we obtain a projected sensitivity of $|\beta_{12}|,|\beta_{13}|\lesssim 0.015$ at 68% C.L., well within the linear regime.

Besides the degeneracy between $\beta_{12}$ and $\beta_{13}$, there are other combinations of parameters that benefit from the combined analysis. For illustration, in Fig. 5 we show the expected sensitivity of the combination ($\beta_{23}$, $\beta_{\delta}$). The JUNO-TAO experiment (dashed-green) is insensitive to $\beta_{23}$ as the $P_{ee}$ channel does not depend on the parameter as seen in Eq. (13). For that reason, we display the JUNO-TAO contour based on the $\Delta\chi^{2}=0.05$ cut. Since both DUNE-ND (red) and FASER$\nu$2 (blue) depend on both parameters we display the contours with the 68% C.L. cut ($\Delta\chi^{2}=2.3$). The shape of the DUNE-ND curve is a combination of the measurement of the survival probability $P_{\mu\mu}$ and the transition probability $P_{\mu e}$. For the two vertical curves, such behavior is determined by $P_{\mu\mu}$ as well as $P_{\mu e}$, which have no correlation between $\beta_{23}$ and $\beta_{\delta}$ as shown in Eq. (10). The two curved branches are affected by the $P_{\mu e}$ transition probability whose exact expression is,

$$P_{\mu e}=\sin^{2}\left(\frac{\Delta\delta}{2}\right)\sin^{2}2\theta_{13}\sin^{2}\theta_{23}^{p},$$ (27)

where $\Delta\delta\equiv\delta(Q^{2}_{d})-\delta(Q^{2}_{p})$ and $\theta_{23}^{p}\equiv\theta_{23}(Q^{2}_{p})$. For $\beta_{23}\approx-0.1$, the linear regime is not valid anymore. Together with $\ln(Q^{2}_{p}/Q^{2}_{0})\approx 10$, a large $\beta_{23}$ leads to $\sin^{2}(\theta_{23}+\Delta\theta_{23})\approx 0$, which drastically reduces the sensitivity to $\beta_{\delta}$. At the same time, FASER$\nu$2 is also not capable of constraining $\beta_{\delta}$ inside the linear regime. Fortunately, the combined result (black curve) provides a strong constraint on the parameters such that $|\beta_{23}|\lesssim 0.01$ and $|\beta_{\delta}|<0.045$ at 68% C. L. For completeness, we checked the impact on the sensitivity to the variation of $\delta_{\rm D}$, we find that the variation of the Dirac CP-phase has negligible effect on the sensitivity curves. For this reason, we only show the result with $\delta_{D}$ prior in Fig. 5.

V Discussion and Conclusions

New physics related to neutrinos may arise at intermediate scales and can generate a RG running of the neutrino mixing angles. In the linear regime, running can be parametrized by model-independent $\beta_{X}$ parameters. Interestingly, a zero-distance transition probability appears due to the mismatch of production and detection momentum transfer and allows for testing new physics in short-baseline detectors.

We show that the presence of several correlation among the $\beta_{X}$ parameters makes it a difficult task for any short-baseline experiment alone to model-independently constrain the RG running parameters. However, combining experiments with several detection channels and energies such as JUNO-TAO, DUNE-ND, and FASER$\nu$2 makes it possible to greatly reduce such degeneracies and provide constrains on $|\beta_{12}|,|\beta_{13}|\lesssim 0.015$, $|\beta_{23}|\lesssim 0.01$, and $|\beta_{\delta}|<0.045$ at 68% C.L.