Compact Perturbative Expressions For Neutrino Oscillations in Matter

We further develop and extend a recent perturbative framework for neutrino oscillations in uniform matter density so that the resulting oscillation probabilities are accurate for the complete matter potential versus baseline divided by neutrino energy plane. This extension also gives the exact oscillation probabilities in vacuum for all values of baseline divided by neutrino energy. The expansion parameter used is related to the ratio of the solar to the atmospheric $\Delta m^2$ scales but with a unique choice of the atmospheric $\Delta m^2$ such that certain first-order effects are taken into account in the zeroth-order Hamiltonian. Using a mixing matrix formulation, this framework has the exceptional feature that the neutrino oscillation probability in matter has the same structure as in vacuum, to all orders in the expansion parameter. It also contains all orders in the matter potential and $\sin\theta_{13}$. It facilitates immediate physical interpretation of the analytic results, and makes the expressions for the neutrino oscillation probabilities extremely compact and very accurate even at zeroth order in our perturbative expansion. The first and second order results are also given which improve the precision by approximately two or more orders of magnitude per perturbative order.


Introduction
Neutrino oscillation based on the standard three flavor scheme provides the best possible theoretical paradigm which can describe most of the experimental results obtained in the atmospheric, solar, reactor, and the accelerator neutrino experiments.In matter, the propogation of neutrinos is significantly modified by the Wolfenstein matter effect [1].The theoretical derivation and understanding of the neutrino oscillation probabilities in matter have been pursued by various means.The exact expressions of the eigenvalues, mixing angles, and the oscillation probabilities have been obtained [2][3][4], albeit under the assumption of uniform matter density.But, the resulting expressions of the oscillation probabilities are way too complex to facilitate understanding of the structure of the three flavor neutrino oscillations.For this reason, analytic approaches to the phenomena are mostly based on variety of perturbative frameworks.For a comprehensive treatment of neutrino oscillation in the matter, see ref. [5].
What is the appropriate expansion parameter in such a perturbative framework?We now know that sin θ 13 , once used as the expansion parameter (there are an enormous number of references, see e.g., [6]), is not so small, sin θ 13 ≃ 0.15.Moreover, expansion around sin θ 13 = 0 misses the physics of the resonance which exists at an energy around E ∼ 10 GeV for earth densities.Therefore, in the environments in which the matter effect is comparable to the vacuum mixing effect, the only available small expansion parameter known to us is the ratio of the solar-scale ∆m 2  ⊙ to the atmospheric-scale ∆m 2 ⊕ , ∆m 2 ⊙ /∆m 2 ⊕ ≃ 0.03.This framework was examined in the past, to our knowledge in refs.[6][7][8][9].
Recently, two of us, see [10], presented a new perturbative framework for neutrino oscillation in matter using a modified ∆m 2  ⊙ /∆m 2 ⊕ expansion.We identified a unique ∆m 2 ⊕ that absorb certain "first-order" terms into the "zeroth-order" Hamiltonian.The resulting expansion parameter, ǫ ≡ ∆m 2 21 /∆m 2 ee where ∆m 2 ee ≡ ∆m 2 31 − sin 2 θ 12 ∆m 2 21 , multiplies a particularly simple perturbing Hamiltonian with zero diagonal entries.This re-organization of the perturbation expansion lead to simple and compact oscillation probabilities in all channels.The ν e disappearance channel is particularly simple, being of a pure two flavor form.
As was noted in [10], this new perturbation expansion, while valid in most of the baseline, L, divided by neutrino energy, E, versus matter potential plane, has issues around vacuum values for the matter potential at large values of L/E.These issues are caused by the crossing of two of the eigenvalues of the new zeroth order Hamiltonian at the solar resonance.In this paper, we solve these issues by performing an additional rotation of the neutrino basis in matter by introducing an additional matter mixing angle which is identical to θ 12 in vacuum.With this extra rotation, the new eigenvalues of the unperturbed Hamiltonian do not cross and the perturbing Hamiltonian remains non-diagonal and is multiplied by an additional factor which is always less than unity and is zero in vacuum.With this additional rotation our perturbative expansion is valid in the full L/E versus matter potential plane and the zeroth order gives the exact result in vacuum.
The sectional plan of this paper is as follows: in section 2 we describe in detail the sequence of rotations of the neutrino basis that leads us to the simple Hamiltonian that will be used in the perturbative expansion.The zeroth order eigenvalues and mixing matrix are given in this section.Then, in section, 3 we explicitly calculate the first and second order corrections for both the eigenvalues and the mixing matrix.In section 4, we give compact analytic expressions for ν e and ν µ disappearance channels as well as ν µ → ν e appearance channel at both zeroth and first order in our perturbative expansion.All other channels can by obtained by unitarity.Here we discuss the precision of the perturbative treatment.Finally, in section 5 there is a conclusion.A number of technical details are contained in the appendices, see A. We have also published the new Nu-Pert code used in this paper online. 1

Rotations of the neutrino basis and the Hamiltonian
In this section we perform a sequence of rotations on the neutrino basis and the corresponding Hamiltonian such that the following conditions are satisfied: • The diagonal elements of the rotated Hamiltonian are excellent approximations to the eigenvalues of the exact Hamiltonian and do not cross for any values of the matter potential.These diagonal elements will form our H 0 .
• The size of non-diagonal elements are controlled by our small parameter, ǫ ′ , which vanishes in vacuum.The non-diagonal elements will form our perturbing Hamiltonian, H 1 .
The first two of these rotations are identical to the rotations performed in [10], while the last rotation is needed to deal with the remaining eigenvalue crossing at the solar resonance.With these three rotations the resulting Hamiltonian satisfies the conditions above and leads us to a rapidly converging perturbative expansion for the oscillation probabilities that covers all of the L/E versus matter potential plane.

Overview
Neutrino evolution in matter is governed by a Schrödinger like equation where in the flavor basis U MNS is the lepton mixing matrix in vacuum, given by and the matter potential, assumed to be constant, is given by (2.1.5) We will perform a sequence of rotations on the flavor basis by multiplying the left and right hand side of eq.2.1.1 by an appropriate unitary matrix, U † and inserting unity (U U † ) between H and |ν .These rotations are chosen such that the final resulting Hamiltonian satisfies the following properties: the diagonal elements are an excellent approximations to the exact eigenvalues and the size of off-diagonal elements are controlled by a small parameter (ratio of the ∆m2 's) and are identically zero in vacuum.
The sequence of rotations applied to the eigenstates is performed in the following order with the corresponding Hamiltonians The first rotation undoes the θ 23 − δ rotation, whereas the φ followed by ψ rotations are matter analogues to the vacuum θ 13 and θ 12 rotations, respectively.In vacuum, the final Schrödinger equation is just the trivial mass eigenstate evolution equation.as defined in [11,12], and the ratio of the ∆m 2 's In terms of the |a| → ∞ eigenvalues the exact Hamiltonian is simple given by3 Note that H is real and does not depend on θ 23 or δ.

U 13 (φ) rotation
Since s 13 ∼ O( √ ǫ), it is natural to diagonalize the (1-3) sector next, using U 13 (φ), again see [10].After this rotation the neutrino basis is and the Hamiltonian is given by where which is identical to eq. 3.1 of [10].The angle, φ, that achieves this diagonalization of the (1-3) sub-matrix (see appendix A.1), satisfies from which it is easy to derive ) The Hamiltonian given in eq.2.3.2 was used to derive simple, compact and accurate oscillation probabilities for a wide range of the L/E versus ρE plane, see [10].However, as was noted in that paper, there is a region of this plane for which a perturbation theory based on Ĥ is insufficient to describe the physics accurately.This region is small ρE and large L/E given by To address this region of the L/E versus ρE plane, we perform one further rotation on the Hamiltonian.This rotation removes the degeneracy of the zeroth order eigenvalues at the solar resonance when λ − = λ 0 .This is performed in the next subsection.

U 12 (ψ) rotation
Since λ − and λ 0 cross at the solar resonance, a ≈ ǫ∆m 2 ee cos 2θ 12 / cos 2 θ 13 , to describe the physics near this degeneracy we need to diagonalize the (1-2) submatrix of Ĥ, using U 12 (ψ).The new neutrino basis is The resulting Hamiltonian, split into a zeroth order Hamiltonian and a perturbing Hamiltonian, is given by where (2.4.4) The diagonal elements of the zeroth order Hamiltonian are (2.4.5) The angle, ψ, that achieves this diagonalization of the (1-2) sub-matrix of Ĥ (see appendix A.1), satisfies ) where we introduce the useful shorthand notation, It is easy to derive that4 and (2.4.10) Figure 1 shows φ and ψ as functions of the matter potential as well as the eigenvalues of Ȟ for both the normal ordering (NO) and the inverted ordering (IO).Several additional useful identities used in the calculations throughout this paper are listed in appendix A.2.

Remarks
A number of summarizing and useful comments are warranted at this point.
• The neutrino basis that will be used in our perturbation theory, |ν is related to the flavor basis, |ν by where • The Hamiltonian, eqs.2.4.3 and 2.4.4,that will used as the basis for our perturbation theory is given by Ȟ with the diagonal elements the zeroth order Hamiltonian and the off-diagonal elements the perturbing Hamiltonian.While the λ a,b,c eigenvalues cross twice and the λ −,0,+ eigenvalues cross once, the new λ 1,2,3 eigenvalues do not cross, see figure 1, which allows for the perturbation theory to be well defined everywhere.• The size of the perturbing Hamiltonian, Ȟ1 , is controlled by the parameter which is never larger than 1.4%.
• In vacuum,  • Since perturbing Hamiltonian, Ȟ1 , has only non-diagonal entries the first order correction to the eigenvalues are zero.The diagonal elements multiplied by 2E are, to an excellent approximation, the mass squares of the neutrinos in matter.
• There is a very useful interchange symmetry involving λ 1,2 and ψ.The Hamiltonian is invariant under the pair of transformations λ 1 ↔ λ 2 and ψ → ψ ± π/2.Our expressions for s ψ and c ψ , see eq. 2.4.10,satisfy this interchange symmetry with the + in front of the π/2.Since the transition probabilities always have an even number of ψ trig functions, this interchange symmetry can be simply expressed as In the rest of this paper we call this the λ 1,2 − ψ interchange symmetry.
• An antineutrino with energy E is equivalent to a neutrino with energy −E.
• The values of all of the eigenvalues in vacuum and for a → ±∞ are shown in appendix A.3.

Perturbation expansion
To calculate the neutrino oscillation probabilities at zeroth order, all that is needed is eigenvalues and mixing matrix, given by eq.2.4.5 and eq.2.5.2 respectively.For higher order calculations we need not only the corrections to the eigenvalues but also the corrections to the mixing matrix.In this section we first given the corrections to the eigenvalues at both first and second order in our expansion parameter, ǫ ′ .This is followed by the corrections to the same order for the mixing matrix.Note that all corrections to both the eigenvalues and the mixing matrix vanish in vacuum as our expansion parameter is zero in vacuum, i.e. the zero order oscillation probabilities are exact in vacuum.

Corrections to the eigenvalues
Since the diagonal terms of Ȟ1 = 0 by construction, the first order corrections to the eigenvalues are exactly zero, since The second order corrections to the eigenvalues are given by5 λ (2) Using Ȟ1 from eq. 2.4.4,we see that the corrections are We verified that the eigenvalues satisfy the characteristic equation to second order, see appendix A.4.The eigenvalues are correct at zeroth order to a fractional precision of about 10 −4 or better, and through second order to a precision of 10 −8 or better.In fact, the precision of λ 1 + λ 1 for sign(∆m 2 ee )Y e ρE < 0 is completely saturated by the limits of double precision computer calculations.

Corrections to the eigenvectors
Here we present the corrections to the eigenvectors which allows us to calculate the transition probabilities to arbitrary order.This was called the V -matrix approach in [13].
First, we relate the flavor eigenvectors to the zeroth order eigenvectors (no subscript) using U m MNS , as in eq.2.1.4, Next, the exact eigenvectors of Ȟ, labeled with subscript (ex), are related to the eigenvectors of Ȟ0 (the zeroth order eigenvectors) by a unitary matrix, which we call Combining the above gives, where 3) The exact V matrix transforms the exact eigenvectors of Ȟ to the flavor basis.In vacuum (a = 0), U m MNS = U MNS and W = ½, so V = U MNS as expected.
Standard perturbation theory in Ȟ1 , which contains the small parameter ǫ ′ , can be used to calculate W † .Here we use a slightly modified perturbation theory to calculate W directly.Expanding W as a power series in ǫ ′ , we define It is clear from eq. 3.2.2 that W 0 = ½.
The first order correction to the W matrix is given by The second order correction, after using the facts that Ȟ1 is symmetric and has no diagonal elements, eq.2.4.4, is This series can be continued to reach arbitrary precision.However, we have found that second order provides more than sufficient precision.In summary the matrix relating the zeroth order eigenvalues of Ȟ0 to the flavor basis is given by to second order in ǫ ′ .Demonstration of the unitary nature of V , to the appropriate order, is given in appendix A.5.With the eigenvalues and eigenvectors determined to second order we can now calculate the neutrino oscillation probabilities.

Oscillation probabilities
In vacuum and in matter with constant density, it is well known that the neutrino oscillation probabilities for ν α → ν β for three-flavor mixing (i, j = 1, 2, 3) can be written in the following form6   .Both V and λ (ex) i s depend on the energy of the neutrino E, and the matter density ρ but the baseline L, dependence only appears in ∆ ij .
By unitarity and using the fact that the sin 2 functions and the triple sine function are linearly independent functions of L, as determined by their non-zero Wronskian, we have to the following powerful statements, Since D αα = 0, we also note that D αβ = −D αγ for α, β, γ all different.So, up to one overall sign, there is only one D term for all channels.
To determine the oscillation probability to n-th order in our perturbative expansion we must evaluate C, D, and ∆λ (ex) ij to the n-th order.We denote this perturbative expansion as follows ∆λ

The first order probabilities
At first order the ∆λ's are again given by eq.2.4.5, since λ (1) i = 0, see eq. 2.4.4,because the diagonal elements of Ȟ1 are zero.The first order corrections to C, D only have terms proportional to ∆λ −1 31 , ∆λ −1 32 .This comes from the form of W 1 , eq. 3.2.5, which follows from the position of the non-zero elements in Ȟ1 .In fact, all of the coefficients can be written in the following general form, where the F 1,2 , G 1,2 and K 1,2 are related by λ 1,2 , ψ interchange previously discussed.Thus only three modest expressions are required to describe the C's and D coefficients to first order for each channel.The F, G, K terms can be calculated from U m MNS by F and G are even under the interchange of α and β whereas K is odd.Their explicit values are given in table 2.
In the appearance channels the CP violating term must be of the following form where in the denominator one needs the exact eigenvalues in matter.This is the Naumov-Harrison-Scott identity, see refs.[14,15].We have checked this identity to the appropriate order, see appendix A.7.
The P (ν α → β) and P (ν α → νβ ) probabilities are related by δ → −δ and the P (ν α → ν β ) and P (ν β → ν α ) transition probabilities are related by L → −L.From eq. 4.0.1,we see that the D term is the only term odd in L. From tables 1 and 2, we see that the D term is also the only one odd in δ, confirming the CPT invariance of these equations.Moreover, all of the D αβ terms are the same order by order up to a coefficient of −1, 0, 1.

The second order probabilities
Although we have not expanded the second order oscillation probabilities analytically, the second order corrections to the eigenvalues, λ i , as well as the second order corrections to  The fractional uncertainties at zeroth and first order are plotted using the analytic formulas in tables 1 and 2 respectively.The probability to second order is calculated by using λ's and W through second order, see eqs.3.1.3and 3.2.6 .
the mixing matrix, W 2 , have been used to calculate the oscillation probabilities to second order.The resulting oscillation probabilities are more than two orders of magnitude closer to the exact values than the first order probabilities.

Precision of the perturbation expansion
The oscillation probabilities that were perturbatively calculated in this section are only useful if they are more precise than the experimental uncertainties.In figure 3, we have plotted the fractional uncertainties8 at each order of our perturbative expansion for the ν µ → ν e channel at the DUNE [16], baseline of 1300 km.The precision at the first oscillation maximum and minimum for DUNE are shown in table 3. We note that the precision improves at lower energies, such as for NOνA [17] and T2K/T2HK [18,19].
The results are comparable for different values of δ, for the inverted ordering, for other channels, and for antineutrino mode.Therefore, even at zeroth order, the precision exceeds the precision of the expected experimental results.

Conclusions
In this paper we have further developed and expanded upon the recent perturbative framework for neutrino oscillations in uniform matter, introduced in [10].The new oscillation probabilities are of the same simple, compact functional form with slightly more complicated coefficients, yet, the range of applicability now includes the whole L/E versus matter potential, a, plane, i.e. the restriction that L/E be small, (L/E ≪ 1/∆m 2 21 ) around the vacuum values of the matter potential has been completely removed.In fact, with these new improvements, the oscillation probabilities in vacuum are exact at zeroth order in our perturbative expansion.This occurs because the expansion parameter s 12 c 12 ∆m 2 21 /∆m 2 ee = 0.014 is further multiplied by s (φ−θ 13 ) , where φ is the mixing angle θ 13 in matter.In vacuum, φ = θ 13 and therefore all corrections to zeroth order vanish.
To achieve this extended range of applicability, an additional rotation of the Hamiltonian is performed over that in [10].The third angle ψ is the mixing angle θ 12 in matter.In the resulting Hamiltonian, the diagonal elements are the eigenvalues of the zeroth order Hamiltonian and do not cross for any values of the matter potential, especially near the solar resonance (this occurred in [10]).The non-diagonal elements of the new Hamiltonian are the perturbing Hamiltonian for our perturbative expansion and their size is controlled by the small parameter s (φ−θ 13 ) s 12 c 12 ∆m 2 21 /∆m 2 ee , mentioned in the previous paragraph.The new perturbative expansion is now well defined for all values of the matter potential and gives very accurate oscillation probabilities.We have performed many cross checks on the perturbative expansion, e.g.we have checked the CP violating term recovers, order by order, the known form.We have calculated the oscillation probabilities for zeroth, first, and second order in our expansion parameter.For most practical applications related to experiments, the zeroth order oscillation probabilities are sufficiently accurate with a typical fractional uncertainty of better than 10 −3 .Including the first and second order corrections the accuracy improves that to better than 10 −6 and 10 −9 , respectively.

3 Figure 1 :
Figure1: The upper figure shows the angles, φ and ψ, as a function of the matter potential for both NO and IO.φ and ψ are the mixing angles θ 13 and θ 12 in matter respectively.For ψ, the curves for the two mass ordering are nearly identical.The two lower figures show the eigenvalues to zeroth order, λ 1,2,3 , in matter as a function of the matter potential for NO and for IO.For all our figures, Y e ρE ≥ 0 is for neutrinos and Y e ρE ≤ 0 for antineutrinos.

Figure 3 :
Figure 3: The ν µ → ν e oscillation probability is plotted in the upper part of the figure for DUNE parameters; a 1300 km baseline and Y e ρ = 1.4 g•cm −3 .The fractional uncertainties at zeroth and first order are plotted using the analytic formulas in tables 1 and 2 respectively.The probability to second order is calculated by using λ's and W through second order, see eqs.3.1.3and 3.2.6 .

Table 3 :
The transition probabilities, energies, and fractional uncertainties at zeroth, first, and second order.Values are calculated at DUNE for ν µ → ν e with the NO and δ = 3π/2.At higher maxima and minima the fractional uncertainties are even smaller.