Collective neutrino oscillations and neutrino wave packets

TL;DR

This paper investigates how wave packet separation causes decoherence in collective neutrino oscillations within dense environments, analyzing the effects in vacuum, matter, and non-adiabatic regimes using density matrix formalism.

Contribution

It provides a detailed analysis of decoherence effects on collective neutrino oscillations, including the derivation of damping terms in different regimes and models of adiabaticity violation.

Findings

Decoherence introduces damping in vacuum and adiabatic matter oscillations.

No damping term generally appears in non-adiabatic regimes.

Estimated neutrino coherence lengths relevant for supernovae environments.

Abstract

Effects of decoherence by wave packet separation on collective neutrino oscillations in dense neutrino gases are considered. We estimate the length of the wave packets of neutrinos produced in core collapse supernovae and the expected neutrino coherence length, and then proceed to consider the decoherence effects within the density matrix formalism of neutrino flavour transitions. First, we demonstrate that for neutrino oscillations in vacuum the decoherence effects are described by a damping term in the equation of motion of the density matrix of a neutrino as a whole (as contrasted to that of the fixed-momentum components of the neutrino density matrix). Next, we consider neutrino oscillations in ordinary matter and dense neutrino backgrounds, both in the adiabatic and non-adiabatic regimes. In the latter case we study two specific models of adiabaticity violation -- one with short-term and another with extended non-adiabaticity. It is demonstrated that, while in the adiabatic case a damping term is present in the equation of motion of the neutrino density matrix (just like in the vacuum oscillation case), no such term in general appears in the non-adiabatic regime.

Equations

$${\sigma }_t\sim \frac{{\sigma }_{x2}}{|{\vec{v}}_2-{\vec{v}}_1|},$$

$$\frac{m\bar{v}^{2}}{2}=\frac{3}{2}T\,,\qquad n_{bar}\approx\frac{\rho}{m}\approx\big{(}\tfrac{4}{3}\pi r_{0}^{3}\big{)}^{-1},$$

$$\bar{v}\approx 0.179\,\Big{(}\frac{T}{10\;\rm MeV}\Big{)}^{1/2}\,,\qquad r_{0}\approx 7.36\times 10^{-13}~{}{\rm cm}\,\Big{(}\frac{\rho}{10^{12}\;\rm g/cm^{3}}\Big{)}^{-1/3}\,.$$

$$\sigma_{x}\,\simeq\,v_{g}\sigma_{t}\,\lesssim\,4.1\times 10^{-12}~{}\rm{cm}\,\Big{(}\frac{\rho}{10^{12}\;\rm g/cm^{3}}\Big{)}^{-1/3}\Big{(}\frac{T}{10\;\rm MeV}\Big{)}^{-1/2}\,,$$

$${\sigma }_x\lesssim 4.1\times 10^{-12}\mathrm{cm}.$$

$$\frac{\sigma_{xe}}{r_{0}/|\vec{v}_{2}-\vec{v}_{1}|}\sim\frac{(4\pi\alpha^{2}n_{e})^{-1/3}}{[(4\pi/3)n_{bar})]^{-1/3}}\Big{(}\frac{3T}{m}\Big{)}^{1/2}\approx 4.2\times\Big{(}\frac{1}{2Y_{e}}\Big{)}^{1/3}\Big{(}\frac{T}{10\;{\rm MeV}}\Big{)}^{1/2}.$$

$$L_{\mathrm{coh}}\simeq {\sigma }_x\frac{v_g}{\Delta v_g},$$

$$\Delta v_g=\frac{\partial }{\partial p}\Delta E,$$

$$L_{\mathrm{coh}}\simeq \frac{2p^2}{\Delta m^2}{\sigma }_x.$$

$${\cal H}=\left(\begin{array}[]{cc}{\cal H}_{11}&{\cal H}_{12}\\ {\cal H}_{12}^{*}&{\cal H}_{22}\end{array}\right).$$

$$E_{1,2}=\frac{{\mathcal{H}}_{11}+{\mathcal{H}}_{22}}{2}\pm \sqrt{\frac{({\mathcal{H}}_{11}-{\mathcal{H}}_{22}{)}^2}{4}+|{\mathcal{H}}_{12}{|}^2}.$$

$$\Delta v_{g}\,=\,\frac{\partial}{\partial p}\bigg{(}2\sqrt{\frac{({\cal H}_{11}-{\cal H}_{22})^{2}}{4}+|{\cal H}_{12}|^{2}}\,\bigg{)}.$$

$$\frac{\Delta v_{g}}{v_{g}}\simeq\frac{\Delta m^{2}}{2p^{2}}\cdot\frac{\frac{\Delta m^{2}}{2p}-\sqrt{2}G_{F}n_{e}c_{20}}{\sqrt{\big{(}\frac{\Delta m^{2}}{2p}c_{20}-\sqrt{2}G_{F}n_{e}\big{)}^{2}+\big{(}\frac{\Delta m^{2}}{2p}s_{20}\big{)}^{2}}}\,,$$

$$\Delta v_g=\frac{\partial }{\partial p}|{\vec{H}}_{\omega }|.$$

$$\Delta v_{g}\,=\,\frac{\partial\omega}{\partial p}\Big{(}\frac{\partial}{\partial\omega}|\vec{H}_{\omega}|\Big{)}\,=\;\frac{\Delta m^{2}}{2p^{2}}\cdot\frac{\omega-\mu P_{0}c_{20}}{\sqrt{(\omega-\mu P_{0}c_{20})^{2}+\mu^{2}(P^{2}-c_{20}^{2}P_{0}^{2})}}\,,$$

$$|\nu (t,\vec{x})⟩=\int \frac{d^{3}p}{(2\pi {)}^3}{e}^{i\vec{p}\vec{x}}|\nu (t,\vec{p})⟩,$$

$$i\frac{d}{dt}|\nu (t,\vec{p})⟩={\mathcal{H}}_{\vec{p}}(t)|\nu (t,\vec{p})⟩,$$

$${\mathcal{H}}_{\vec{p}}(t)=U_0\frac{M^2}{2p}{U}_0^{†}+\hat{V}(t)+\sqrt{2}{G}_F\int \frac{d^{3}q}{(2\pi {)}^3}[\rho (t,\vec{q})-\bar{\rho }(t,\vec{q})](1-{\vec{v}}_{\vec{p}}{\vec{v}}_{\vec{q}}).$$

$$\rho (t)=\int \frac{d^{3}q}{(2\pi {)}^3}\rho (t,\vec{q})=\int \frac{4\pi q^{2}dq}{(2\pi {)}^3}\rho (t,|\vec{q}|),$$

$${\mathcal{H}}_{\vec{p}}(t)=U(t){\mathcal{H}}_{\vec{p}}^d(t)U^{†}(t),$$

$${\cal H}_{\vec{p}}^{d}(t)=\operatorname{diag}\big{(}E_{1}(t,\vec{p}),E_{2}(t,\vec{p}),...\big{)}\,.$$

$$|{\nu }_{\alpha }(\vec{p})⟩=\sum _j{U}_{\alpha j}^{\ast }(t)|{\nu }_j(t,\vec{p})⟩.$$

$$|\nu (t_0,\vec{p})⟩=|{\nu }_{\alpha }(\vec{p})⟩=\sum _j{U}_{\alpha j}^{\ast }(t_0)|{\nu }_j(t_0,\vec{p})⟩.\text{vspace}\ast -2.5mm$$

$$|{\nu }_j(t_0,\vec{p})⟩=f_{{\vec{p}}_j}(\vec{p})|{\nu }_j^{(0)}(t_0,\vec{p})⟩,$$

$$⟨{\nu }_j^{(0)}(t,\vec{p})|{\nu }_k^{(0)}(t,\vec{p}{}^{\prime })⟩=(2\pi {)}^3{\delta }^3(\vec{p}-\vec{p}{}^{\prime }){\delta }_{jk}.$$

$$\int \frac{d^{3}p}{(2\pi {)}^3}|f_{{\vec{p}}_j}(\vec{p}){|}^2=1.$$

$$f_{\vec{p}_{j}}(\vec{p})=(2\pi/\sigma_{p}^{2})^{3/4}\,\exp\Big{[}-\frac{(\vec{p}-\vec{p}_{j})^{2}}{4\sigma_{p}^{2}}\Big{]}.\vspace*{1.0mm}$$

$$|\nu (t_0,\vec{p})⟩=\sum _j{U}_{\alpha j}^{\ast }(t_0)f_{{\vec{p}}_j}(\vec{p})|{\nu }_j^{(0)}{t}_0,\vec{p}⟩.\text{vspace}\ast -2.5mm$$

$${\hat{\rho }}^a(t,\vec{x})=|{\nu }^a(t,\vec{x})⟩⟨{\nu }^a(t,\vec{x})|,$$

$$\hat{\rho }(t,\vec{x})=\sum _a{\hat{\rho }}^a(t,\vec{x})=\sum _a|{\nu }^a(t,\vec{x})⟩⟨{\nu }^a(t,\vec{x})|.$$

Full text

1 Introduction

Neutrino oscillations in dense neutrino environments existing in the early Universe and at certain stages of supernova explosions may differ significantly from oscillations in vacuum or in ordinary matter. In dense neutrino backgrounds, coherent neutrino–neutrino interactions may strongly affect the flavour evolution of the neutrino system, leading to a variety of collective oscillation phenomena. These include synchronized oscillations [1, 2, 3, 4, 5, 6, 7, 8], bipolar oscillations [9, 5, 10, 11, 12, 13, 14], spectral splits and swaps [15, 16, 17, 18] and multiple spectral splits [19]. These phenomena have attracted a great deal of attention recently, see Refs. [20, 21] for reviews and extensive lists of literature.

Since in supernovae and in the early Universe neutrinos are produced at very high densities, their production processes are characterized by a high degree of space-time localization and, as a result, the wave packets (WPs) describing the states of produced neutrinos are very short in configuration space. In particular, in [22, 23] the estimate $\sigma_{x}\sim 10^{-11}$ cm was obtained for the spatial length of supernova neutrino WPs by considering the properties of the WPs of the electrons involved in the neutrino production. In Section 2.1 of this paper we arrive at a similar result, basing our estimate on different considerations. The WPs of different propagation eigenstates222 By propagation eigenstates we mean the eigenstates of the effective neutrino Hamiltonian, which in general includes neutrino potentials coming from neutrino coherent forward scattering on the particles of the usual matter and on background neutrinos. composing a given produced flavour-eigenstate neutrino in general propagate with different group velocities. Therefore, very short neutrino WPs suggest that the packets of the individual propagation eigenstates would quickly separate in space, leading to decoherence and hence suppression of neutrino oscillations. Simple estimates of the coherence length of supernova neutrinos, i.e. the distance over which the WPs of their propagation eigenstate components fully separate, yield $L_{\rm coh}\sim 10$ km [22]. One might therefore expect that decoherence by WP separation would significantly affect flavour transitions of supernova neutrinos, and, in particular, collective neutrino oscillations, which are expected to occur at a distance $\sim 100$ km from the supernova’s center. Numerical simulations show, however, that in reality the situation is more complex.

Decoherence effects on collective neutrino oscillations are very little studied, with most existing work focusing on understanding the evolution of neutrinos in the coherent regime. To our knowledge, the only papers that directly address the conditions under which coherence is lost are Refs. [8] and [24].333Decoherence effects on flavour transformations of supernova neutrinos were also studied in Ref. [23]. However, the influence of the coherence loss on collective neutrino oscillations has not be considered there. In [8] it has been demonstrated that in a dense uniform and isotropic neutrino gas decoherence due to late-time dephasing of different neutrino momentum modes does not occur provided that the neutrino number density is high enough. As a result, synchronized neutrino oscillations that can take place in such a system are not destroyed by decoherence. In [8] an interpretation of this fact in terms of mode locking due to neutrino self-interaction was suggested; however, the WP picture of neutrino flavour transitions was not explicitly discussed. In Ref. [24] possible decoherence effects in dense neutrino gases (or their absence) were studied from the point of view of spatial separation of neutrino WPs, both in the adiabatic and non-adiabatic regimes of neutrino flavour evolution. The results of [8] were confirmed from this standpoint; in addition, some ad hoc assumptions used in the analytic approach developed in [8] were shown to be redundant.

In the present paper we continue studies of the effects of decoherence by WP separation on collective neutrino oscillations in dense neutrino gases. Our goal is in particular to develop a space-time picture of such decoherence effects, complementary to the more widely used momentum space approach, in which the neutrino ensemble is split into individual modes of fixed momentum. We address the question of whether WP separation can be described by effective evolution equations for the density matrix of a neutrino WP as a whole (rather than of its components for a given momentum or energy). We derive such effective equations and discuss in detail their domain of applicability.444The present paper supersedes and replaces our previous paper on this subject [25] (currently withdrawn), in which the domain of validity of this effective evolution equation was overestimated.

We would like to emphasize here that our configuration-space approach to decoherence is effectively equivalent to the usual one based on the momentum-space considerations, provided that the momentum distributions of the individual neutrino states are properly taken into account. However, as we shall see, it allows new insights into the space-time development of decoherence.

The paper is organized as follows. In Section 2 we estimate the spatial length of the WPs of neutrinos produced in core collapse supernovae and discuss the neutrino coherence length with respect to WP separation. In Section 3 we derive and discuss equations of motion (EoMs) satisfied by the neutrino density matrix in dense matter and neutrino backgrounds. We start by reviewing the evolution equation for the neutrino state vector (Section 3.1), the neutrino density matrix in flavour space (Section 3.2) and the density matrix describing neutrino oscillations in vacuum for neutrinos described by Gaussian WPs (Section 3.3). In Section 3.3 we also derive the EoMs satisfied by the coordinate-averaged and un-averaged neutrino density matrices in vacuum and show that these equations essentially coincide, provided that in the latter case the total time derivative is understood as the Liouville operator. Moreover, we demonstrate that the evolution equation of the density matrix of a neutrino as a whole (as contrasted to the equations for fixed-momentum components of the density matrix), as well as that of the coordinate-space density matrix, explicitly contains a damping term which describes decoherence by WP separation. In Section 3.4 we consider the EoM for the neutrino density matrix in normal matter and dense neutrino backgrounds, both in the adiabatic (Section 3.4.1) and in the non-adiabatic (Section 3.4.2) regimes. For the latter, we consider two specific models of adiabaticity violation in Sections 3.4.3 and 3.4.4. We show that, just like for vacuum oscillations, the EoM of the neutrino density matrix in the adiabatic regime contains a damping term which leads to decoherence by WP separation at sufficiently late times. However, no such term in general appears in the non-adiabatic regime. We summarize and discuss our results in Section 4. In the Appendix we briefly review the flavour spin formalism for description of neutrino oscillations in dense neutrino backgrounds and derive within this formalism the expression for the non-adiabaticity parameter in the case of 2-flavour neutrino oscillations. We also give there an approximate analytic solution of the EoM for neutrino flavour spin vectors in the case of synchronized neutrino oscillations.

2 Supernova neutrino wave packets and the coherence

length

An important parameter characterizing decoherence phenomena in neutrino propagation is the coherence length $L_{\rm coh}$. This is the distance traveled by neutrinos beyond which the WPs corresponding to different propagation eigenstates composing the produced neutrino flavour eigenstate separate by more than their spatial length $\sigma_{x}$. Once this has happened, the neutrino WPs no longer overlap, and their coherence is lost.555 Decoherence by WP separation can in principle be undone by neutrino detection if the individual detection processes are characterized by a sufficiently high energy resolution (i.e. are sufficiently delocalized in space and time) [26]. However, at least for terrestrial detection of neutrinos from a galactic supernova this would require unrealistically high energy resolution, and so we disregard this possibility. To find the coherence length of supernova neutrinos, let us first estimate $\sigma_{x}$.

2.1 The length of supernova neutrino WPs

In Refs. [22, 23] the length $\sigma_{x}$ of supernova neutrino WPs was estimated assuming that it is dominated by the length of the WPs of the charged leptons (mostly electrons and positrons) participating in the neutrino production process. However, as we show below, while this assumption is correct in the case of neutrino production in $e^{+}e^{-}$ annihilation, it is actually not valid when neutrinos are produced in processes with participation of nucleons.

It has been demonstrated in [27, 28] that the length of the neutrino WP is determined by the spatial and temporal localization of the neutrino production process and is typically dominated by the latter. The temporal localization of neutrino production (i.e. the duration of the production process) is given by the overlap time $\sigma_{t}$ of the WPs of all particles taking part in the process. In ref. [27] it has been shown that this time can usually be estimated as

[TABLE]

where $\vec{v}_{1}$ and $\vec{v}_{2}$ are the velocities of the particles with the shortest and next-to-shortest WPs in configuration space, and $\sigma_{x2}$ is the spatial length of the next-to-shortest WP. This result has a very simple interpretation: since the production process requires the overlap of the WPs of all participating particles, it is over as soon as at least one of these WPs ceases to overlap with all the others. The time scale of neutrino production is therefore typically given by the interval during which the shortest WP “slides” along the second shortest one, which leads to the estimate in eq. (1). Note, however, that this result may need a modification if there is a longer WP which moves much faster than the shortest and the second shortest ones. We will discuss this point in more detail later on.

For neutrino production processes with participation of nucleons, the shortest and the second shortest WPs are those of the involved nucleons, and their lengths can be estimated as $\sigma_{x1}\sim\sigma_{x2}\lesssim r_{0}$, where $r_{0}$ is the average distance between the nucleons in the medium. The relative velocities of the nucleons are of the order of their mean thermal velocity $\bar{v}$. The parameters $r_{0}$ and $\bar{v}$ can be estimated from the relations

[TABLE]

where $\rho$, $n_{bar}$ and $T$ are the matter density, baryon number density and temperature at the neutrino production point, respectively, and $m$ is the nucleon mass. This gives

[TABLE]

For the length $\sigma_{x}$ of the neutrino WP we then find

[TABLE]

where $v_{g}\simeq 1$ is the mean group velocity of neutrino propagation eigenstates. At density $\rho\sim 10^{12}$ g/cm3 and temperature $T\sim 10$ MeV, typical for the neutrino production region in a supernova, we have

[TABLE]

Note that $\sigma_{x}$ depends only rather weakly on the matter density at neutrino production: for neutrinos produced close to the electron neutrino sphere, where $\rho\sim 5\times 10^{10}$ g/cm3 [29], $\sigma_{x}$ is only about a factor of 2.7 larger than our estimate (5).

Eq. (1), on which the above estimates are based, has been obtained in [27] from a more general (but less transparent) relation under the assumptions that the velocities $\vec{v}_{1}$ and $\vec{v}_{2}$ are not nearly equal, and that the magnitudes of the velocities of other particles taking part in neutrino production are not vastly different from $v_{1}$ and $v_{2}$. The latter condition is actually not satisfied for neutrino production in supernovae, as electrons are relativistic at temperatures $T\sim 10$ MeV, while nucleons are non-relativistic. Let us therefore consider the effect of the electron WP on the duration of the neutrino production process.

It has been demonstrated in [22] that the length of the WPs of electrons participating in neutrino production in supernovae can be estimated as $\sigma_{xe}\simeq(4\pi\alpha^{2}n_{e})^{-1/3}$, where $\alpha$ is the fine structure constant and $n_{e}$ is the electron number density. For $\rho\sim 10^{12}$ g/cm3 and the electron fraction per baryon $Y_{e}\sim 1/2$, this gives $\sigma_{xe}\sim 10^{-11}$ cm, which is much larger than the lengths of the nucleon WPs. However, since electrons move much faster than nucleons, it is conceivable that the overlap time of the electron and nucleon WPs may be shorter than the nucleon–nucleon overlap time. The duration of the neutrino production process, and thus the length of the neutrino WP, would then be determined by the electron–nucleon overlap time. To see if this is indeed the case, let us compare electron–nucleon and nucleon–nucleon overlap times. For the ratio of these quantities we have

[TABLE]

Thus, the overlap time of electron and nucleon WPs is typically larger than the nucleon-nucleon overlap time by about a factor of four or more. This means that the effect of the electron WPs on the temporal localization of the neutrino production processes with participation of nucleons (and therefore on the length of the WPs of produced neutrinos) is small. For Gaussian WPs its relative contribution to the length of neutrino WPs is suppressed by a factor  $\lesssim e^{-8}\sim 3\times 10^{-4}$. Neutrino interactions on the way from their birthplace in the supernova core to the neutrinosphere (from which they essentially stream out freely) do not increase the length of the neutrino WPs.

Interestingly, even though the considerations on which our estimate of $\sigma_{x}$ was based were quite different from those employed by Kersten and Smirnov in [22, 23], numerically our result is rather close to theirs, $\sigma_{x}\sim 10^{-11}$ cm.

2.2 Coherence length of supernova neutrinos

Let us now turn to the coherence length of supernova neutrinos $L_{\rm coh}$. The propagation eigenstates $\nu_{i}$ composing the initially produced flavour eigenstate neutrino in general move with different group velocities $v_{gi}$; as a result, after having traveled a distance

[TABLE]

their WPs separate by a distance exceeding their length $\sigma_{x}$ and therefore cease to overlap. To find the coherence length $L_{\rm coh}$, one therefore needs to know, in addition to $\sigma_{x}$, the difference of their group velocities $\Delta v_{g}$. The group velocities are well defined for the propagation eigenstates and are given by the derivatives of their energies with respect to momentum. Therefore,

[TABLE]

where $\Delta E$ is the difference of the energies of the corresponding neutrino propagation eigenstates, i.e. the difference of the eigenvalues of the effective neutrino Hamiltonian. For relativistic neutrinos in vacuum (free propagation) this gives $\Delta v_{g}\simeq\Delta m^{2}/(2p^{2})$, where $p$ is the neutrino momentum, and eq. (7) yields

[TABLE]

It should be stressed that in a system with more than two neutrino flavours there is a separate coherence length $L_{\rm coh}^{jk}$ for each pair of propagation eigenstates $\nu_{j}$, $\nu_{k}$. It is the distance beyond which the WPs of $\nu_{j}$ and $\nu_{k}$ cease to overlap. Similarly, the quantities $\Delta v_{g}$, $\Delta E$ and $\Delta m^{2}$ in eqs. (7), (8) and (9) should also carry indices $jk$. We have suppressed these indices to simplify the notation. Complete decoherence occurs when the distance $L$ traveled by neutrinos satisfies $L>\max\{L_{\rm coh}^{jk}\}$ for all pairs $\{j,k\}$; however, partial decoherence may also be of interest.

For neutrino oscillations in a medium, one has to take into account that the effective neutrino Hamiltonian depends on the properties of the medium. Consider for definiteness the 2-flavour case, in which the effective Hamiltonian is a $2\times 2$ Hermitian matrix

[TABLE]

Its eigenvalues are

[TABLE]

In this case

[TABLE]

In particular, for neutrino oscillations in ordinary matter, one finds from the standard expression for the neutrino Hamiltonian ${\cal H}$ that [43]

[TABLE]

where $G_{F}$ is the Fermi constant, $s_{20}=\sin 2\theta_{0}$, and $c_{20}=\cos 2\theta_{0}$, with $\theta_{0}$ being the leptonic mixing angle in vacuum. It is interesting to note that, when the density is such that $\sqrt{2}G_{F}n_{e}c_{20}=\Delta m^{2}/(2p)$, the quantity $\Delta v_{g}/v_{g}$ vanishes. For small values of the vacuum mixing angle $\theta_{0}$, the density required to satisfy this condition is very close to the one at the Mikheyev–Smirnov–Wolfenstein (MSW) resonance, where $\sqrt{2}G_{F}n_{e}=c_{20}\Delta m^{2}/(2p)$. Note that for small $\theta_{0}$, the density interval where $\Delta v_{g}/v_{g}$ is suppressed is very narrow; outside this region the absolute value of the last factor on the right hand side of eq. (13) is of order one, i.e. $|\Delta v_{g}|/v_{g}$ is of the same order of magnitude as $\Delta m^{2}/(2p^{2})$.

Consider next the case of 2-flavour oscillations in a dense uniform and isotropic neutrino gas. In this case it is convenient to study neutrino flavour evolution in the framework of the “flavour spin” formalism [30]. In this picture, for each neutrino mode $\omega=\Delta m^{2}/(2p)$, the traceless part of the neutrino density matrix in flavour space $\rho(t,\vec{p})$ is written as a linear combination of the Pauli matrices. The coefficients of the Pauli matrices form a flavour spin vector $\vec{P}_{\omega}$, whose evolution is a precession around the (in general time-dependent) Hamiltonian vector $\vec{H}_{\omega}$. The latter, in turn, is formed by the coefficients of the expansion of the traceless part of the neutrino Hamiltonian in terms of the Pauli matrices (see the Appendix). The vectors $\vec{H}_{\omega}$ depend on the global flavour spin vector $\vec{P}$, which is the sum of all the individual $\vec{P}_{\omega}$. It is easy to show that in the 2-flavour case the difference of the eigenvalues of the effective neutrino Hamiltonian is $\Delta E=|\vec{H}_{\omega}|$, so that

[TABLE]

Substituting here the expression for $|\vec{H}_{\omega}|$ from eq. (A9) yields

[TABLE]

where $\mu=\sqrt{2}G_{F}n_{\nu}$, $P=|\vec{P}|$, and $P_{0}$ is defined as $P(t=0)$. Note that the global flavour spin vector $\vec{P}$ is in general time dependent, and so are $\Delta v_{g}$ and $L_{\rm coh}$. We did not indicate the time dependence in eq. (15) explicitly in order not to overload the notation. Just as in the case of oscillations in ordinary matter, the absolute value of the last factor on the right hand side of eq. (15) is typically of order unity except in a relatively narrow region around the point in the parameter space where $\omega-\mu P_{0}c_{20}=0$.666Indeed, the second term under the square root in eq. (15) vanishes for $\mu\leq\mu_{0}$ (where $\mu_{0}$ is the threshold value of $\mu$ discussed in the Appendix) because in this case $P\to c_{20}P_{0}$; thus the absolute value of the last factor in eq. (15) is essentially equal to one. For $\mu\gg\mu_{0}$ we have $P\to P_{0}$ and the the last factor in eq. (15) is $\simeq-c_{20}$. For intermediate values of $\mu$ the modulus of this factor is also of order unity.

Thus, in most of the parameter space a simple estimate $\Delta v_{g}/v_{g}\simeq\Delta m^{2}/(2p^{2})$ (and therefore the expression for $L_{\rm coh}$ in eq. (9)) should be quite accurate, i.e. the quantities $\Delta v_{g}$ and $L_{\rm coh}$ are to a good approximation time independent.

It should be emphasized that the discussion of the coherence length in this section, strictly speaking, applies only to the oscillations in vacuum and to adiabatic flavour evolution in medium. In the non-adiabatic regime the very notion of the coherence length may lose its sense; we shall discuss this issue in more detail in Sections 3.4.4 and 4.

3 Equation of motion for neutrino density matrix

3.1 Evolution of the neutrino state vector

We start with considering the evolution of an individual neutrino in a background of other neutrinos and ordinary matter. The coordinate-space neutrino state vector $|\nu(t,\vec{x})\rangle$ can be Fourier-expanded as

[TABLE]

where $|\nu(t,\vec{p})\rangle$ is the time-dependent state vector for a neutrino of a given momentum $\vec{p}$. We will assume that $|\nu(t,\vec{p})\rangle$ satisfies the standard Schrödinger-like evolution equation

[TABLE]

In the following we will often be using the matrix representation, where ${\cal H}_{\vec{p}}(t)$ is an $N_{f}\times N_{f}$ Hermitian matrix with $N_{f}$ being the number of neutrino flavours. Then in the flavour basis [30]

[TABLE]

Here $M$ is the diagonal matrix of neutrino masses, $U_{0}$ is the leptonic mixing matrix in vacuum, $\hat{V}(t)$ is the matrix of potentials due to coherent forward scattering of neutrinos on ordinary matter, $G_{F}$ is the Fermi constant, $\rho(t,\vec{p})$ and $\bar{\rho}(t,\vec{p})$ are the density matrices of, respectively, neutrinos and antineutrinos of momentum $\vec{p}$ (to be discussed in more detail below), and $\vec{v}_{\vec{p}}$ is the velocity of neutrinos with momentum $\vec{p}$. The last term on the right hand side of eq. (18) is due to coherent forward scattering of the test neutrino on background neutrinos. Note that the first term in (18) depends only on $p=|\vec{p}|$. In isotropic matter and neutrino backgrounds, the same is true for the whole Hamiltonian; in particular, the last term in (18) reduces to $\sqrt{2}G_{F}[\rho(t)-\bar{\rho}(t)]$, where

[TABLE]

and similarly for $\bar{\rho}(t)$.

Let us introduce the neutrino propagation eigenstates and leptonic mixing matrix in medium. At any time $t$ the Hamiltonian ${\cal H}_{\vec{p}}(t)$ can be diagonalized by a unitary transformation

[TABLE]

where $U(t)$ is the leptonic mixing matrix in medium and ${\cal H}_{\vec{p}}^{d}(t)$ is the diagonal matrix of the eigenvalues $E_{j}(t,\vec{p})$ of ${\cal H}_{\vec{p}}(t)$, i.e.

[TABLE]

Because of the time dependence of ${\cal H}_{\vec{p}}(t)$, its eigenstates and eigenvalues are also $t$-dependent; at a given $t$ they are called instantaneous eigenstates and eigenvalues.777Since the Hamiltonian ${\cal H}_{\vec{p}}(t)$ depends on the density matrices $\rho$ and $\bar{\rho}$, which themselves must be obtained from the evolution equations based on eq. (17), the instantaneous eigenstates and eigenvalues of ${\cal H}_{\vec{p}}(t)$ are not known a priori, unlike in the case of neutrino oscillations in ordinary matter. For our discussion, however, it is sufficient to know that such eigenstates and eigenvalues exist. At any fixed instant of time the neutrino flavour states $|\nu_{\alpha}(\vec{p})\rangle$ can be represented as linear superpositions of the propagation eigenstates $|\nu_{j}(t,\vec{p})\rangle$:

[TABLE]

Since the mixing matrix in medium $U(t)$ diagonalizes the $\vec{p}$-dependent Hamiltonian ${\cal H}_{\vec{p}}(t)$, it is also $\vec{p}$-dependent; we do not indicate this dependence explicitly to simplify the notation.

Let a neutrino be produced at a time $t_{0}$ in the flavour state $\alpha$ ($\alpha=e,\mu$ or $\tau$). According to eq. (22), its state vector can be written as a linear combination of the state vectors of the propagation eigenstates, so that

[TABLE]

We will be assuming that the propagation eigenstates $\nu_{j}$ composing the initially produced flavour state are described by WPs with momentum-space wave functions (momentum distribution amplitudes) $f_{\vec{p}_{j}}(\vec{p})$, where $\vec{p}_{j}$ is the centroid of the momentum distribution. For each $\vec{p}$-component of the WP one can then write

[TABLE]

where $|\nu_{j}^{(0)}(t,\vec{p})\rangle$ are the state vectors of the propagation eigenstates satisfying

[TABLE]

For the amplitudes $f_{\vec{p}_{j}}(\vec{p})$ we adopt the normalization condition

[TABLE]

The momentum distribution is characterized by its effective width $\sigma_{p}$. An example we will often be using below is the Gaussian WP:

[TABLE]

The plane-wave limit is recovered for $\sigma_{p}\to 0$, which yields $f_{\vec{p}_{j}}(\vec{p})=[(2\pi)^{3}/\sqrt{V}]\delta^{3}(\vec{p}-\vec{p}_{j})$, where $V$ is the normalization volume.

Combining eqs. (23) and (24), we find for the initial state

[TABLE]

The neutrino state vector $|\nu(t,\vec{p})\rangle$ at $t>t_{0}$ can then be found by solving the EoM (17) with the initial condition (28). The coordinate space neutrino state vector in the WP picture is then obtained from eq. (16).

Let us now return to the neutrino flavour evolution equation. In non-uniform matter and neutrino backgrounds, the potential $\hat{V}$ and the neutrino density matrix $\rho$ should obviously depend also on the neutrino coordinate; then how can the coordinate-independent eqs. (17) and (18) be obtained? In deriving these equations, neutrinos are implicitly assumed to be essentially pointlike objects. For a pointlike neutrino moving along a classical trajectory $\vec{x}=\vec{x}(t)$ its coordinate is fully determined by its propagation time, and one can consider $\hat{V}$ and $\rho$ as functions of $t$ only, rather than of $(t,\vec{x})$. In other words, what matters is only the values of $\hat{V}$, $\rho$ and $\bar{\rho}$ at the point reached by the neutrino at time $t$.

The question remains, of course, as to how this can be reconciled with the fact that eqs. (17) and (18) are supposed to describe flavour evolution of neutrinos of definite momenta, which correspond to plane waves rather than to pointlike states. The answer is that these ‘plane waves’ should be considered as the limit of wave packets (WPs) when their momentum spread $\sigma_{p}$, though still small compared to the momentum $p$, is nevertheless large compared to the inverse of the distances over which $\hat{V}$, $\rho$ and $\bar{\rho}$ vary significantly.

Let us discuss this point in more detail. Eqs. (17) and (18) (as well as the corresponding equations in the case of the usual MSW effect [31, 32]) are usually obtained under the assumptions that neutrinos are relativistic and that the potentials they feel are small compared to the neutrino energy and vary little over spatial distances of the order of the neutrino de Broglie wavelength $1/p$ (see, e.g., [33, 34]). When neutrinos are described by WPs of finite length, one has to add the conditions that $\hat{V}$ and $\rho$ vary little also over distances of the order of the spatial size of the neutrino WP $\sigma_{x}$ and over time scales of the order of the neutrino passage time $\sim\sigma_{x}/v_{\nu}$. Since usually $\sigma_{x}\sim 1/\sigma_{p}\gg 1/p$, these conditions are more stringent than the standard ones; still, for the potential $\hat{V}$ they are satisfied very well for all cases of practical interest. As for the density matrices $\rho$ and $\bar{\rho}$ entering into eq. (18), in the case when neutrinos are described by WPs, their definition must include averaging over spatial volumes whose linear size is large compared to the length of the neutrino WP [35, 36]. Under this condition they also vary little over distances $\sim\sigma_{x}$ and times $\sim\sigma_{x}/v_{\nu}$. Once the above conditions are satisfied, neutrinos can be described by finite-length WPs (for which $\vec{x}$ and $t$ are independent variables), while their individual momentum components satisfy the EoM in eqs. (17) and (18).

3.2 Neutrino density matrix

A convenient quantity for discussing neutrino flavour evolution in dense media is the neutrino density operator (density matrix) in flavour space. The 1-particle density operator describing flavour transitions of an individual neutrino in a medium consisting in general of both ordinary matter and background neutrinos is defined as

[TABLE]

where the evolution of the state vector $|\nu^{a}(t,\vec{x})\rangle$ is described by eqs. (16)–(18). We have introduced here the superscript $a$ which labels the individual neutrinos ($a=1,...,N$, where $N$ is the total number of neutrinos in the system). The density operator of the whole neutrino system, whose Fourier transform enters into the definition of the effective Hamiltonian ${\cal H}_{\vec{p}}(t)$ in eq. (18), is the sum of the 1-particle ones:888The one-particle density operators $\hat{\rho}^{a}(t,\vec{x})$, as well as the relation in eq. (30), are well defined provided that dynamical multiparticle correlations in the system can be neglected. The latter become sizeable only at extremely high densities, at which the mean distance between the neutrinos in the system is comparable with the range of the weak interactions. Such densities are far above those present in the supernovae.

[TABLE]

The density matrix $\hat{\rho}(t,\vec{x})$ along with the evolution equation (17) and eqs. (16) and (18) describe the flavor content of the neutrino ensemble in the mean field approximation. If the overall neutrino spectrum correctly takes into account the spectra of the individual neutrino WPs, both the standard approach, in which the system is treated as an ensemble of neutrinos of well defined energies, and the wave packet approach result in the same $\hat{\rho}(t,\vec{x})$ and therefore describe the same physics.

The coordinate dependence of the neutrino density operator $\hat{\rho}(t,\vec{x})$ has two sources. First, the density operator of each individual neutrino has a non-trivial $\vec{x}$-dependence because it corresponds to a finite-size WP. This type of coordinate dependence disappears when one averages the density operator over spatial volumes of a linear size $\bar{r}\gg\sigma_{x}$ around each point $\vec{x}$. ($\bar{r}$ should, however, be small compared to the distances over which the macroscopic characteristics of the medium change significantly). We shall refer to this type of averaging as coarse graining. The second source of $\vec{x}$-dependence is related to the fact that the properties of the medium are in general coordinate dependent, and in particular different regions of space may be characterized by different concentration and different flavour composition of neutrinos. This type of coordinate dependence, which is the only one that survives upon coarse graining, is, however, absent in the case of uniform matter and neutrino backgrounds.

Likewise, also the momentum dependence of the neutrino density matrix has two sources. The first one is related to the fact that each individual neutrino is described by a WP characterized by a certain momentum distribution amplitude. The density matrix of the neutrino as a whole is obtained by integrating over its individual momentum modes. The resulting 1-particle density matrix will depend on the mean momentum of the neutrino state, which is given by the weighted average of the centroids of the momentum distributions of the propagation eigenstates composing the produced favour state:

[TABLE]

The individual neutrino state vector (16) and the corresponding expression for the density matrix given in eq. (29) actually depend on this mean momentum. Therefore the summation over the index $a$ should include the summation (integration) over the spectrum of the mean momenta $\langle\vec{p}\rangle_{a}$.

What characterizes the state of the $a$th neutrino? Given the shape of its momentum-space wave function, each neutrino is characterized, in addition to the mean momentum discussed above, by its production time $t_{0}^{a}$, production coordinate $\vec{x}_{0}^{a}$, and initial flavour $\alpha$, as well as by the evolution time $t$. The latter no longer uniquely characterizes the state of the $a$th neutrino, which depends also on the initial conditions, which are different for different neutrinos. In non-uniform matter and neutrino backgrounds, the state vectors of different neutrinos will depend on their individual histories and will no longer be given by the solutions of the same evolution equation (17) – at a given time $t$ the Hamiltonians ${\cal H}_{\vec{p}}(t)$ will be different for different neutrinos.

The situation is much simpler in the case of isotropic and uniform neutrino backgrounds. Isotropy means that the velocity terms in eq. (18) can be omitted and eq. (19) can be utilized. Uniformity means that it is sufficient to consider the coordinate-independent EoM (17), which for a given neutrino production time $t_{0}^{a}$ fully determines the individual neutrino state. If all neutrinos are produced at the same time $t_{0}$ (as it is often assumed in toy models), then the density operators of the individual neutrinos are identical, and the density operator of the complete neutrino system is just given by the total neutrino number $N$ times the 1-particle density operator. In reality, different neutrinos are usually produced at different times $t_{0}^{a}$, and therefore averaging of the 1-particle density operators over the production time has to be performed. We will discuss this issue in Section 3.4.3.

We shall now address the question if and when EoMs describing the evolution of the density operators of the individual neutrinos as wholes (rather than of their specific momentum components) can be found. Our goal, in particular, is to examine if such EoMs would exhibit terms responsible for the wave packet separation, thus facilitating studies of decoherence effects. We start with the simplest case of neutrino oscillations in vacuum.

3.3 Neutrino oscillations in vacuum

For oscillations in vacuum, the Hamiltonian of the neutrino system is given by just the first, time-independent, term on the right hand side of eq. (18). The mixing matrix $U$ in eq. (23) then coincides with the leptonic mixing matrix in vacuum $U_{0}$, and the propagation eigenstates $\nu_{j}$ coincide with the neutrino mass eigenstates. In this subsection we will deal only with 1-particle density matrices which are fully adequate to neutrino oscillations in vacuum.

Let a flavour eigenstate $\nu_{\alpha}$ be produced at $t_{0}=0$ and $\vec{x}_{0}=0$.999 By this we mean that the centroid of the neutrino WP is located at the production time $t_{0}=0$ at the point with the coordinate $\vec{x}_{0}=0$. Here and below we neglect the finite duration of the neutrino production process, which is important for the formation of the neutrino WP (see Section 2), but, once its properties are defined, is irrelevant for its subsequent evolution, which we are mostly interested in. The neutrino state vector at $(t,\vec{x})$ can then be written as

[TABLE]

where $\psi_{j}(t,\vec{x})$ is the coordinate-space wave function (wave packet) of the $j$th neutrino mass eigenstate. From eq. (29) it then follows that for neutrino oscillations in vacuum the 1-particle neutrino density matrix in the mass-eigenstate basis can be written as

[TABLE]

We shall now assume for definiteness that the mass-eigenstate neutrinos constituting the produced flavour state are described in momentum space by Gaussian WPs of width $\sigma_{p}$, see eq. (27). Expanding the neutrino energies $E_{j}(p)=(p^{2}+m_{j}^{2})^{1/2}$ about the peak momenta $p_{j}$ and retaining only the first two terms in the expansion (which amounts to neglecting WP spreading effects),101010 Effects of spreading of WPs are strongly suppressed for ultra-relativistic neutrinos. Although they may have observable consequences for neutrinos propagating over astrophysical distances, they can be safely neglected for neutrino oscillations in supernovae, and in any case they do not affect decoherence issues [23]. one finds that the coordinate space WPs also have Gaussian form. For the density matrix in the mass eigenstate basis this yields [37]

[TABLE]

where

[TABLE]

and similarly for $E_{k}$, $\vec{v}_{k}$. Note that $\vec{v}_{j}$ and $\vec{v}_{k}$ are the group velocities of the $j$th and $k$th neutrino mass eigenstates, respectively.

In the realistic situations we are interested in, the density matrix is averaged over spatial volumes of a linear size $\bar{r}$ that is small compared to the distances over which macroscopic characteristics of the system change significantly, but large compared to the size $\sigma_{x}$ of the individual neutrino WPs. Since the WP amplitude decreases very quickly at distances much greater than $\sigma_{x}$ from their center, one can formally extend the integration over the coordinate to infinity and consider the quantity

[TABLE]

For the Gaussian WPs from eq. (34) we find

[TABLE]

where $\vec{v}_{g}\equiv\frac{\vec{v}_{j}+\vec{v}_{k}}{2}$ is the average group velocity of the states $\nu_{j}$ and $\nu_{k}$. We did not indicate explicitly the $j,k$ dependence of $\vec{v}_{g}$ because for relativistic neutrinos this dependence is small and because $v_{g}$ enters into eq. (37) multiplied by $|\vec{p}_{j}-\vec{p}_{k}|$, which is small. Indeed, coherent production of flavour eigenstates is only possible for $|\vec{p}_{j}-\vec{p}_{k}|\ll\sigma_{p}$. Keeping the $j,k$ dependence of $\vec{v}_{g}$ would therefore introduce terms $\ll\!(\Delta m^{2}/2E^{2})\sigma_{p}t$. At the same time, we keep the last term in the exponent in eq. (37), which is of the order of $(\Delta m^{2}/2E^{2})^{2}\sigma_{p}^{2}\,t^{2}$ because it can become important at sufficiently large times.

Next, we notice that for ultra-relativistic neutrinos

[TABLE]

where $\Delta E\equiv E_{j}-E_{k}$, $\Delta\vec{p}\equiv\vec{p}_{j}-\vec{p}_{k}$, $\Delta m^{2}=m_{j}^{2}-m_{k}^{2}$, and $E$ is the average energy of $\nu_{j}$ and $\nu_{k}$. Introducing the notation

[TABLE]

we can rewrite eq. (37) as

[TABLE]

Differentiating this expression, we find the EoM satisfied by $\rho_{jk}(t)$:

[TABLE]

Here

[TABLE]

is the coherence length, i.e. the distance over which the WPs of different neutrino mass eigenstates separate by a distance of the order of their spatial length $\sigma_{x}$. Eq. (41) resembles the standard evolution equation for the fixed-momentum components of the density matrix, except that the usual energies $E_{j}$ are replaced by the effective ones $\tilde{E}_{j}$ and that an extra term (the last term in the square brackets) appears. It should be stressed, however, that we are now considering the density matrix of the neutrino as a whole, which includes integration over all its momentum modes. It is easy to see that the extra term in (41) leads to damping with time of the off-diagonal elements of $\rho_{jk}(t)$, which is a consequence of the spatial separation of the WPs of different mass eigenstates. Note also that eq. (38) implies

[TABLE]

which is just the quantity relevant for neutrino oscillations.

It is also interesting to study the EoM for the un-averaged density matrix $\rho_{jk}(t,\vec{x})$. From eq. (34) we obtain

[TABLE]

[TABLE]

Combining these two equations, we find

[TABLE]

where we have used eq. (39). Recall now that the total time derivative, when applied to a classical function of time, coordinate and momentum, can be written as the Liouville operator

[TABLE]

If we replace here $\dot{\vec{x}}$ with the average group velocity $\vec{v}_{g}$ and take into account that for free neutrinos (as well as for neutrinos experiencing only forward scattering) the neutrino momentum is conserved, we find that this operator formally coincides with the one acting on $\rho_{jk}(t,\vec{x})$ on the left hand side of (46). Modulo the identification of this operator with the total time derivative, eq. (46) would essentially coincide with eq. (41), i.e. the averaged and un-averaged density matrices would satisfy the same EoM. This fact is rather curious: strictly speaking, the Liouville operator is only applicable for pointlike particles moving along classical trajectories, so that $\vec{x}=\vec{x}(t)$. But from our considerations it follows that it still applies to the case of particles represented by finite-size WPs, provided that their motion is described in terms of the average group velocity.

It is easy to show that our evolution equation (41) actually follows from the standard momentum-space EoM of the fixed-momentum components of the neutrino density matrix. (We shall demonstrate this in Section 3.4.1 in the more general case of neutrino propagation in medium in the adiabatic regime.) The advantage of eq. (41) for studying propagation decoherence is that it allows one to immediately see the effects of the spatial separation of the WPs of the different mass eigenstates, without integration over different neutrino momentum modes.

At this point, some comments are in order.

We have derived the EoM (41) for $\rho_{jk}$ simply by differentiating the already known solution. Of course, if the solution is known, we do not need any EoM. Our point here was just to demonstrate that the fact that the off-diagonal elements of $\rho_{jk}$ decrease with time is reflected in the appearance of an extra term in the EoM, provided that we consider as our object of interest the density matrix of an individual neutrino as a whole, rather than the density matrices of its specific momentum modes. The density matrix for a neutrino ensemble can then be obtained by summing over all neutrinos in the ensemble. 2. 2.

The form of the damping term in eq. (41) depends on the assumed shape of the neutrino WPs. For example, for exponential rather than Gaussian shape of the coordinate-space neutrino WPs,111111Neutrinos with such WPs are produced, e.g., in decays of free or quasi-free parent particles [38, 39, 40]. the extra term would have the form $-1/L_{\rm coh}$ rather than $-2t/L_{\rm coh}^{2}$, as in eq. (41) [25]. This leads to a dependence of the solutions of the neutrino evolution equation on the shape of the WPs at short and intermediate times; the solutions at asymptotically large times, however, exhibit vanishing off-diagonal elements of $\rho_{jk}$ whenever the damping term is present in the EoM, irrespectively of the exact form of this term. 3. 3.

The suppression with time of the off-diagonal elements of the neutrino density matrix in the mass eigenstate basis $\rho_{jk}$ is the direct consequence of the fact that its elements with $j\neq k$ contain products of the wave functions of the neutrino mass eigenstates propagating with different group velocities (see eq. (33)). Therefore, their overlap, which is maximal at neutrino production, decreases with time. This suppression, which can be immediately seen from the expression for the coordinate-integrated density matrix $\rho(t)$ in eq. (37), is actually present in the un-integrated matrix (34) as well. This can be readily seen if we rewrite in the exponent of (34)

[TABLE] 4. 4.

Physically, the reason for the suppression of $\rho_{jk}$ with $j\neq k$ with time is the separation of the WPs corresponding to different propagation-eigenstate components of each neutrino. One might consequently expect this phenomenon to occur also in more general situations – for neutrinos propagating in ordinary matter or in dense neutrino gases. As we shall see in Sections 3.4.1 and 3.4.4, this is indeed the case in the adiabatic neutrino propagation regime, but in general is not true when adiabaticity is violated.

3.4 Oscillations in normal matter

and neutrino backgrounds

We now turn to neutrino oscillations in normal matter and dense neutrino environments in the WP picture. Let us first consider evolution of neutrinos of a given momentum $\vec{p}$ in the propagation eigenstate basis.121212Oscillations in dense neutrino gases have been considered from the viewpoint of the propagation eigenstate basis in [41, 42], but the issue of WPs and their separation has not been addressed there. As was discussed in Section 3.1, the propagation eigenstates $\nu_{j}$ are the instantaneous eigenstates of the Hamiltonian of the neutrino system ${\cal H}_{\vec{p}}(t)$, i.e. the states that diagonalize ${\cal H}_{\vec{p}}(t)$ at the time $t$. We shall be using the matrix representation (18) of the neutrino Hamiltonian, and the column vector $\Psi$, the elements of which are the probability amplitudes of finding the corresponding flavor eigenstate neutrino at time $t$, to describe the neutrino state. The element $\Psi_{\beta}$ of the vector $\Psi$ is related to the state vector in the Hilbert space $|\nu(t,\vec{p})\rangle$ by

[TABLE]

and similarly for $|\nu^{p}\rangle$ and $\Psi^{p}$ describing the propagation eigenstates. Taking into account that the two bases are related by the leptonic mixing matrix $U(t)$, we find from EoM (17) the neutrino evolution equation in the propagation eigenstate basis:

[TABLE]

Here the diagonal matrix ${\cal H}^{d}(t)$ was defined in eqs. (20) and (21), and we suppressed the dependence of the Hamiltonians on the neutrino momentum to simplify the notation.

Note that, unlike ${\cal H}^{d}(t)$, the effective Hamiltonian $\tilde{{\cal H}}(t)$ that governs the evolution of the propagation eigenstates is in general not diagonal because $U^{\dagger}(t)\dot{U}(t)$ is not. This is related to the fact that the Hamiltonian ${\cal H}(t)$ cannot be diagonalized by the same unitary transformation at all times. If, however, the properties of the medium vary slowly enough along the neutrino path, the neutrino system has enough time to ‘adjust’ itself to the changes of these properties. In this case the term $U^{\dagger}(t)\dot{U}(t)$ in (50) can be neglected, i.e. the system evolves adiabatically. The effective Hamiltonian $\tilde{\cal H}(t)$ then essentially coincides with ${\cal H}^{d}(t)$, i.e. is diagonal. This means that transitions between different propagation eigenstates are strongly suppressed, and they evolve independently of each other.

3.4.1 Adiabatic evolution

Let us first consider neutrino flavour evolution in normal matter and neutrino backgrounds in the adiabatic regime. Let a neutrino be produced at a time $t_{0}$ and coordinate $\vec{x}_{0}$ in the flavour state $\alpha$ ($\alpha=e,\mu$ or $\tau$). According to eq. (23), its state vector can be written as a linear combination of the state vectors of the propagation eigenstates. Consider now the time evolution of the produced neutrino state. In the adiabatic regime the propagation eigenstates evolve independently, and their time evolution is very simple – they just acquire phase factors:

[TABLE]

This equation should be understood as follows. In general, the state to which $|\nu_{j}^{(0)}(t_{0},\vec{p})\rangle$ evolves as time advances from $t_{0}$ to $t$ can be expanded in the basis of the eigenstates $|\nu_{k}^{(0)}(t,\vec{p})\rangle$ of the Hamiltonian ${\cal H}(t)$. In the adiabatic regime, this expansion contains only one term (with $k=j$), the expansion coefficient being the exponential factor in (51).

The evolved neutrino state can now be obtained by applying eq. (51) to (28):

[TABLE]

Substituting this into (16), we obtain

[TABLE]

where

[TABLE]

For the 1-particle neutrino density matrix in the propagation eigenstate basis we then obtain

[TABLE]

Notice that the quantities $U_{\alpha j}(t_{0})$ and $U_{\alpha k}(t_{0})$ here implicitly depend on $\vec{p}$ and $\vec{p}\,^{\prime}$, respectively.

It is interesting to note that in the special case when the initially produced neutrino flavour state coincides at $t=t_{0}$ with one of the instantaneous propagation eigenstates, i.e. $U_{\alpha j}(t_{0})=\delta_{\alpha j}$, the off-diagonal elements of the propagation-basis density matrix $\rho_{jk}(t,\vec{x})$ vanish at all times. This, however, is only true in the adiabatic approximation.

Let us expand $E_{j}(t,\vec{p})$ and $E_{k}(t,\vec{p}\,^{\prime})$ in the integrand of eq. (55) near the peaks of the corresponding momentum distribution functions and keep only the first two terms in the expansion.131313 The higher-order terms are responsible for the spreading of the wave packets, which we neglect here. This yields

[TABLE]

where

[TABLE]

and similarly for $E_{k}(t)$ and $\vec{v}_{k}(t)$. Note that $\vec{v}_{j}(t)$ is the group velocity of the $j$th neutrino propagation eigenstate. In general, it is time-dependent in medium because so is the effective Hamiltonian ${\cal H}_{\vec{p}}(t)$. In the following, it will be convenient for us to introduce the mean group velocities $\langle\vec{v}_{j}(t)\rangle$ according to

[TABLE]

Substituting the expansion (56) into eq. (55), we find

[TABLE]

[TABLE]

Here the integration is performed over the shifted momenta $\vec{q}=\vec{p}-\vec{p}_{j}$ and $\vec{q}\,^{\prime}=\vec{p}\,^{\prime}-\vec{p}_{k}$. We have also taken into account that the leptonic mixing matrix in medium $U(t)$ varies little over momentum intervals comparable to the width $\sigma_{p}$ of the momentum distribution function $f_{\vec{p}_{j}}(\vec{p})$; this allowed us to replace $U_{\alpha j}^{*}(t_{0})$ and $U_{\alpha k}(t_{0})$ by their values taken, respectively, at $\vec{p}=\vec{p}_{j}$ and $\vec{p}\,^{\prime}=\vec{p}_{k}$ and pull them out of the momentum integrals.

The integrals over momenta in (59) factorize, and we obtain

[TABLE]

where

[TABLE]

is the envelope function of the coordinate-space WP of $\nu_{j}$. In particular, for Gaussian momentum-space neutrino WPs (27), the envelope functions are also Gaussian:

[TABLE]

Let us now go to the averaged neutrino density matrix $\rho(t)$ by integrating $\rho(t,\vec{x})$ over spatial regions of linear size $\bar{r}\gg\sigma_{x}$ around $\vec{x}$, i.e. by performing coarse graining. As was discussed in Section 3.3, since the coordinate-space WPs quickly decrease at distances from their centers that are large compared with $\sigma_{x}$, one can formally extend the integration over $\vec{x}$ to an infinite spatial volume. Eq. (60) then yields

[TABLE]

where $\rho_{jk}(t)$ was defined in eq. (36) and

[TABLE]

Note that, while $\rho_{jk}(t,\vec{x})$ corresponds to a pure state (i.e. its $j$- and $k$-dependence factorizes), this is not true for the spatially averaged quantity $\rho_{jk}(t)$. It is important to emphasize that, strictly speaking, the ‘coarse-grained’ quantity (63) depends solely on time only in the case of uniform ordinary matter and neutrino backgrounds; otherwise, it will also depend on the coordinate of the point around which the coarse graining with $\bar{r}\gg\sigma_{x}$ is performed. In what follows, we will concentrate on the case of uniform media, though, with minimal modifications, some of our results will apply to the case of non-uniform backgrounds as well.

Let us introduce the notation

[TABLE]

Then we can rewrite eq. (64) as

[TABLE]

[TABLE]

and $\rho_{jk}(t)$ takes the form

[TABLE]

Here, in analogy with eq. (39), we have defined

[TABLE]

The envelope function $g_{j}(\vec{z})$ is peaked at or close to the zero of its argument and decreases rapidly when the argument becomes larger than the spatial width of the WP $\sigma_{x}$. The splitting of the arguments of $g_{j}$ and $g_{k}^{*}$ in the integrand of (67) is $\Delta v(t)(t-t_{0})$; when this splitting exceeds $\sigma_{x}$ (i.e., $t-t_{0}$ exceeds the coherence length), the overlap of $g_{j}$ and $g_{k}^{*}$ gets suppressed, leading to the suppression of the off-diagonal elements of $G_{jk}(t)$.

Let us now find the EoM satisfied by $\rho_{jk}(t)$. From eq. (68) we obtain

[TABLE]

where as usual dots denote time derivatives. The last term in eq. (70) describes the suppression with time of the off-diagonal elements of the neutrino density matrix in the propagation eigenstate basis. In particular, for Gaussian WPs (62) we find

[TABLE]

Eq. (70) then yields

[TABLE]

where the coherence length is given by

[TABLE]

Note that eq. (72) has a form similar to that in the case of vacuum oscillations, eq. (41) (where $t_{0}$ was set equal to zero). Interestingly, the expression for $L_{\rm coh}(t)$ in eq. (73) contains, along with the difference of the mean velocities $\langle\vec{v}_{j}(t)\rangle-\langle\vec{v}_{k}(t)\rangle$, the difference of the instantaneous time-dependent velocities $\vec{v}_{k}(t)-\vec{v}_{k}(t)$. It was pointed out in [43] that in the case of neutrino oscillations in ordinary matter this velocity difference changes its sign close to the MSW resonance, which may in principle lead to the complete or partial restoration of propagation coherence. In Section 2.2 of the present paper it was shown that a similar phenomenon takes place in the case of oscillations in dense neutrino backgrounds. This, however, happens in a very narrow region of parameter space, outside of which the time dependence of $\vec{v}_{j,k}(t)$ is very weak and the difference between $\langle\vec{v}_{j}(t)\rangle-\langle\vec{v}_{k}(t)\rangle$ and $\vec{v}_{k}(t)-\vec{v}_{k}(t)$ essentially disappears (see Section 2.2). In this case the expression for $L_{\rm coh}(t)$ in eq. (73) simplifies to

[TABLE]

which is similar to eq. (42).

For future discussion and to compare our results with those obtained in the standard approach,141414We recall that by the standard approach we mean the formalism in which individual neutrino momentum modes are evolved in time without reference to wave packets. In momentum space, introducing wave packets corresponds to coarse graining in the momentum variable. it will be convenient to express $\rho_{jk}(t)$ and the time derivatives of $\rho_{jk}(t,\vec{x})$ and $\rho_{jk}(t)$ also directly from eq. (55), taking the integral over the coordinates in the expression for $\rho_{jk}(t)$ and $\dot{\rho}_{jk}(t)$ before the momentum integrals. For $\rho_{jk}(t)$ we find

[TABLE]

Differentiating (55) and (75) with respect to time yields, respectively,

[TABLE]

[TABLE]

Note that $\rho_{jk}(t)$ can be obtained either as the integral of $\rho_{jk}(t,\vec{x})$ over the coordinate or, equivalently, as the integral of its Fourier transform $\rho_{jk}(t,\vec{p})$ over momentum; it is then clear that the integrand in eq. (75) is just $\rho_{jk}(t,\vec{p})$. The integrand of eq. (77) therefore represents the right hand side of the standard EoM for the fixed-momentum components of the density matrix

[TABLE]

written in the neutrino propagation-eigenstate basis. Indeed, in the adiabatic regime the Hamiltonian $\tilde{\cal H}_{\vec{p}}(t)$ is diagonal in the propagation eigenstate basis, and the commutator in eq. (78) leads to the factor \big{[}E_{j}(\vec{p})-E_{k}(\vec{p})\big{]} in the integrand of eq. (77), whereas the rest of the integrand is just $\rho_{jk}(t,\vec{p})$. Although at first sight eq. (77) does not have a damping term leading to the suppression of the off-diagonal elements of $\rho_{jk}(t)$, such a damping is actually present there. Indeed, neutrino propagation decoherence can be described either in momentum space or in configuration space. In momentum space it stems from the dephasing of different neutrino modes at late times and is related to the fast oscillations of the integrand in the integral over the neutrino spectrum in (75). In configuration space decoherence is related to the spatial separation of the WPs of different propagation eigenstates after they have traveled long enough distance. The momentum space and configuration space descriptions are equivalent (see, e.g., [28]). In the case we consider here, this follows from the fact that eq. (75) is equivalent to eq. (68), whose time derivative (70) contains the damping term $\propto\dot{F}_{jk}/F_{jk}$.

3.4.2 Beyond the adiabatic approximation

Consider now non-adiabatic neutrino flavour evolution. It is important to note that, unlike for neutrino oscillations in ordinary matter, in dense neutrino environments adiabaticity violation may take place even in the case of constant neutrino density. This comes about because in this case the Hamiltonian of the system depends not only on the overall density of the neutrino gas, but also on its flavour composition, which changes with time.

If neutrino evolution is non-adiabatic, the evolution of different propagation eigenstates is not independent. This means that a propagation eigenstate $\nu_{i}$ produced at a time $t_{0}$ will at a later time $t$ become a linear superposition of different propagation eigenstates:

[TABLE]

The matrix $S_{ji}(t,t_{0};\vec{p})$ is the neutrino evolution matrix in the propagation eigenstate basis, i.e. the evolution matrix that solves eq. (50). It is related to the corresponding evolution matrix in the flavour basis $S^{\rm fl}(t,t_{0};\vec{p})$ by

[TABLE]

Applying eq. (79) to (28), we obtain the evolved neutrino state vector in momentum space:

[TABLE]

Following then the same steps as in Section 3.4.1, for the coordinate-averaged 1-particle neutrino density matrix in the propagation eigenstate basis we find

[TABLE]

In the adiabatic limit, $S_{ji}(t,t_{0};\vec{p})=e^{-i\int_{t_{0}}^{t}E_{j}(t^{\prime},\vec{p})dt^{\prime}}\delta_{ji}$, and we recover eq. (75).

It is important that, unlike in the adiabatic case, the off-diagonal elements of the density matrix (82) do not in general vanish when the initially produced neutrino state coincides with one of the instantaneous propagation eigenstates, i.e. when $U_{\alpha i}(t_{0})=\delta_{\alpha i}$. Indeed, in this case

[TABLE]

and its $j\neq k$ elements are obviously non-zero in the non-adiabatic case.

Let us now consider neutrino evolution in the non-adiabatic regime in more detail. An important point is that in the presence of neutrino self-interactions the Hamiltonian ${\cal H}(t)$ is in general complex (though is still, of course, Hermitian) even in the absence of fundamental CP violation. Indeed, as follows from eq. (18), ${\cal H}(t)$ depends on the neutrino density matrix $\rho$. Since the off-diagonal elements of $\rho$ are in general complex, so is ${\cal H}(t)$.This, in particular, means that ${\cal H}(t)$ is diagonalized by a complex unitary transformation rather than by a real orthogonal one. This applies even to the 2-flavour case, in which in the absence of neutrino self-interactions the mixing matrix $U$ contains no physical complex phases. Indeed, since the phase stemming from neutrino-neutrino interactions is time-dependent, it cannot be rotated away by a simple time-independent redefinition of one of the neutrino fields; therefore, it affects their time-dependent phase difference and hence the probabilities of the flavour transitions.

For simplicity, we will confine our consideration to the 2-flavour case. The unitary matrix $U(t)$ diagonalizing ${\cal H}(t)$ can then be taken in the form

[TABLE]

where $c\equiv\cos\theta(t)$ and $s\equiv\sin\theta(t)$. The angle $\theta(t)$ and the phase $\phi(t)$ can be expressed through the matrix elements of ${\cal H}(t)$:

[TABLE]

Next, we recall that for the purposes of studying neutrino flavour transitions one can always choose ${\cal H}(t)$ to be traceless by using the transformation ${\cal H}(t)\to{\cal H}(t)-\tfrac{1}{2}tr[{\cal H}(t)]\!\cdot\!\mathbbm{1}$. Eq. (20) then yields

[TABLE]

where

[TABLE]

Let us now find the effective Hamiltonian in the propagation eigenstate basis $\tilde{{\cal H}}(t)$. From the definition of $\tilde{{\cal H}}(t)$ in eq. (50) and eq. (84) we obtain

[TABLE]

As usual, the probability of transitions between different propagation eigenstates (i.e. the degree of adiabaticity violation) can be characterized by their mixing, i.e. by the modulus of the ratio of twice the off-diagonal element of $\tilde{{\cal H}}(t)$ to the difference of its diagonal elements. The non-adiabaticity parameter is thus

[TABLE]

For $\lambda\ll 1$, transitions between the propagation eigenstates are strongly suppressed and they evolve practically independently, whereas for $\lambda\gtrsim 1$ strong transitions between them are possible.

The appearance of $\dot{\phi}$ in the expression for $\lambda$ is easy to interpret in the framework of the flavour spin description of neutrino flavour transitions. As discussed in Section 2 and in the Appendix, for each $\omega$-mode neutrino flavour evolution is described by the precession of the flavour spin vector $\vec{P}_{\omega}$ around the (in general time-dependent) Hamiltonian vector $\vec{H}_{\omega}$. The angle between $\vec{H}_{\omega}$ and the 3rd axis in flavour space (the polar angle) is given by 2$\theta(t)$, whereas $\phi(t)$ is the azimuthal angle characterizing the direction of the projection of $\vec{H}_{\omega}$ on the 1–2 plane (see the Appendix). Adiabaticity means that the flavour spin vectors $\vec{P}_{\omega}$ are fast precessing around their corresponding $\vec{H}_{\omega}$, whereas the direction of each $\vec{H}_{\omega}$ changes relatively slowly in flavour space. In that case the flavour spins ‘track’ the motion of the $\vec{H}_{\omega}$. For neutrino oscillations in ordinary matter the effective Hamiltonian ${\cal H}_{\omega}(t)$ is real, so that $\phi(t)\equiv 0$. In that case the vector $\vec{H}_{\omega}(t)$ always lies in the 1–3 plane in flavour space, and its direction is fully characterized by the angle $\theta(t)$. With changing matter density the direction of $\vec{H}_{\omega}$ changes in the 1–3 plane with speed $2\dot{\theta}$, and the adiabaticity is good when this speed is small compared to the angular velocity $\Delta$ of precession of the flavour spin vector $\vec{P}_{\omega}$ around $\vec{H}_{\omega}$. As a result, the (non)adiabaticity parameter is given by the standard expression $\lambda=2|\dot{\theta}|/|\Delta|$.

In the case of oscillations in neutrino backgrounds, each vector $\vec{H}_{\omega}$ performs a 3-dimensional motion, and the speed of the change of its direction depends both on $\dot{\theta}$ and on the angular velocity $\dot{\phi}$ of the rotation of the projection of $\vec{H}_{\omega}$ onto the 1–2 plane. Therefore, knowing just $\dot{\theta}$ would not be sufficient to characterize the degree of adiabaticity violation during neutrino evolution. For example, for $\vec{H}_{\omega}$ precessing around the 3rd axis, the angle $\theta$ remains constant, while $\phi$ changes. Note that such a situation can be approximately realized for synchronized oscillations in dense uniform and isotropic neutrino gas. As shown in the Appendix, in this case $\dot{\theta}$ oscillates around zero or small average, taking both positive and negative values, whereas $\dot{\phi}$ oscillates around a non-zero value, remaining for most of the parameter space sign-definite (see eqs. (A27) and (A28)). Therefore, the situation when $\theta\simeq\text{const.}$ while $\phi$ changes with time obtains for $\dot{\theta}$ and $\dot{\phi}$ averaged over their fast oscillations.

The evolution matrix of propagation eigenstates $S_{ji}(t,t_{0};\vec{p})$ cannot be found in a model independent way. To find out if damping due to WP separation is operative in the non-adiabatic regime, we consider two illustrative models of adiabaticity violation.

3.4.3 A model with short-time non-adiabaticity

In this section we adopt the following simple model of adiabaticity violation: a spatially uniform and isotropic 2-flavour neutrino system in which neutrinos are produced in a flavour $\alpha$ at a time $t_{0}$, propagate adiabatically until $t_{1}-\varepsilon$, experience extremely non-adiabatic evolution during a short time interval $[t_{1}-\varepsilon,\,t_{1}+\varepsilon]$, and then evolve again adiabatically until time $t$.151515Such a situation may be realized e.g. in supernovae when neutrinos cross a propagating shock wave: the adiabaticity is then expected to be good everywhere except inside the shock wave front. We thank Alexei Smirnov for this observation. As follows from eq. (89), extreme non-adiabaticity means that during the interval $[t_{1}-\varepsilon,\,t_{1}+\varepsilon]$ the condition

[TABLE]

is satisfied. Barring possible fine-tunings, this requires, in particular, that during this interval either $|\dot{\phi}|\ll|\dot{\theta}|$, or $c_{2}\ll 1$ (i.e. $\theta(t)\simeq\pi/4=\text{const.}$ and so $|\dot{\theta}|\varepsilon\ll 1$). In the rest of this subsection we will assume that the first of these two possibilities is realized, i.e. $\dot{\phi}$ can be neglected; the second possibility can be considered quite analogously.

The eigenvalues and eigenstates of the effective Hamiltonian ${\cal H}_{\vec{p}}(t)$ before and after the non-adiabatic transition may be drastically different. Let us denote

[TABLE]

[TABLE]

where $\theta(t_{1}-\varepsilon)$ and $\theta(t_{1}+\varepsilon)$ are the values of the in-medium mixing angle immediately before and immediately after the non-adiabatic transition. The evolution matrix for the propagation eigenstates can be written as

[TABLE]

Here the rightmost, central and leftmost matrices describe, respectively, the adiabatic evolution in the interval $[t_{0},t_{1}-\varepsilon]$, extremely non-adiabatic evolution during $[t_{1}-\varepsilon,t_{1}+\varepsilon]$, and then again the adiabatic evolution in the interval  $[t_{1}+\varepsilon,t]$. Direct calculation gives

[TABLE]

Note that due to the smallness of $\varepsilon$ one can actually set $\varepsilon\to 0$ in the integration limits in (91), so that

[TABLE]

Let us now concentrate on the off-diagonal elements of $\rho_{jk}(t)$. We will first consider the question of whether they exist and are damped with time when the initially produced neutrino state coincides with one of the instantaneous propagation eigenstates or when different neutrinos are produced at different times and so the averaging over $t_{0}$ has to be performed when going from the 1-particle neutrino density matrix to the complete one. Recall that in the adiabatic regime in both these cases the density matrix of the neutrino would have no off-diagonal terms in the propagation eigenstate basis.

Consider first the case when the flavour eigenstate $\nu_{\alpha}$ produced at $t_{0}$ coincides with one of the instantaneous propagation eigenstates, which we take to be $\nu_{1}$ for definiteness, i.e. $U_{\alpha i}(t_{0})=\delta_{\alpha 1}$. Then eqs. (82) and  (95) yield, for $t>t_{1}$,

[TABLE]

Next, we assume that the produced flavour eigenstates do not coincide with one of the propagation eigenstates, but the production time $t_{0}$ of different neutrinos is different, and averaging over it results in the vanishing of the phase factors depending on $t_{0}$. Then we obtain from (82)

[TABLE]

where we have assumed $\nu_{\alpha}=\nu_{e}$ and used $|U_{e1}(t_{0})|^{2}=\cos^{2}\theta(t_{0})$, $|U_{e2}(t_{0})|^{2}=\sin^{2}\theta(t_{0})$. The bars over $\sin^{2}(\theta_{0})$ and $\cos^{2}(\theta_{0})$ in (97) denote averaging of these quantities over the production time $t_{0}$.

Let us now discuss eqs. (96) and (97). We can see that in both considered cases ($U_{\alpha i}(t_{0})=\delta_{\alpha i}$ or averaging over $t_{0}$) the off-diagonal elements of $\rho_{jk}(t)$ are non-zero. Moreover, the integration in the exponents in eqs. (96) and (97) is over the time interval $[t_{1},t]$ rather than over the full time interval $[t_{0},t]$. This means that in both cases the off-diagonal elements of $\rho_{jk}(t)$ do not depend on the evolution of the system before the time $t_{1}$. It can also be readily seen from eqs. (96) and (97) that these elements should exhibit the usual damping with time due to WP separation, but the damping will become important only for $t-t_{1}\gtrsim L_{\rm coh}$ rather than for $t-t_{0}\gtrsim L_{\rm coh}$.

Thus, we have demonstrated that in the non-adiabatic regime the off-diagonal elements of the density matrix in the propagation eigenstate basis may be present and be damped with time due to WP separation even when no such off-diagonal elements are initially present (and therefore no WP separation is operative) in the case of pure adiabatic evolution.

3.4.4 A model with extended adiabaticity violation

The matrix $S_{ji}(t,t_{0};\vec{p})$ can also be readily found in a model in which non-adiabatic transitions between propagation eigenstates proceed at a near constant pace. Let us elaborate on this. Let $\tau$ be the characteristic time over which the probability for the propagation eigenstates $\nu_{1}$ and $\nu_{2}$ to convert into each other takes its maximum value allowed by the degree of adiabaticity violation. We will be assuming that the matrix elements of the neutrino Hamiltonian in the propagation eigenstate basis $\tilde{\cal H}$ (88) change very little during this time interval. In this case we can consider them as approximately time-independent. Up to an irrelevant overall phase factor, we then obtain for the evolution matrix $S_{ji}(t,t_{0};\vec{p})$

[TABLE]

where

[TABLE]

In the considered case, transitions between the propagation eigenstates are oscillatory; the characteristic time $\tau$ mentioned above is then just half the oscillation period, $\tau=\pi/(2\Omega)$. The maximal amplitude of the $\nu_{1}\leftrightarrow\nu_{2}$ transitions is $|2\dot{\theta}\pm is_{2}\dot{\phi}|/(2\Omega)=\lambda/\sqrt{1+\lambda^{2}}$, with $\lambda$ being the non-adiabaticity parameter defined in eq. (89). Note that in the adiabatic limit $\dot{\theta}=\dot{\phi}=0$ eqs. (98) and (99) yield $S_{ji}(t,t_{0};\vec{p})=\operatorname{diag}(e^{i\frac{\Delta}{2}t},\,e^{-i\frac{\Delta}{2}t})$, in agreement with the results of Section 3.4.1.

Let us now assume that the initially produced neutrino is a $\nu_{e}$ and consider the general expression for the element $\rho_{12}(t)$ of the density matrix in the 2-flavour case. From eqs. (82) and (84) we have

[TABLE]

Here we have taken into account that, since the neutrinos are produced at $t=t_{0}$ as flavour eigenstates, the phase $\phi(t_{0})$ vanishes. Note that in the adiabatic limit the evolution matrix $S_{ji}(t,t_{0};\vec{p})$ is diagonal, and only the second term in the curly brackets survives.

We shall now evaluate eq. (83) with the model (98) for the evolution matrix $S(t,t_{0};\vec{p})$. It will be convenient for us to introduce the angle $\tilde{\theta}$ which characterizes the mixing of the propagation eigenstates, and thus the adiabaticity violation. The non-adiabaticity parameter $\lambda$ introduced in eq. (89) coincides with $\tan 2\tilde{\theta}$.161616The mixing angle of the propagation eigenstates $\tilde{\theta}$ must not be confused with the angle $\theta$ describing flavour mixing in medium or its in-vacuum analogue $\theta_{0}$. Denoting $\tilde{s}_{2}\equiv\sin 2\tilde{\theta}$, $\tilde{c}_{2}\equiv\cos 2\tilde{\theta}$, $\varphi\equiv\text{arg}[\tilde{\cal H}_{12}(t)]$, from eq. (88) we find

[TABLE]

With this notation, we can rewrite eq. (98) as

[TABLE]

where

[TABLE]

For the off-diagonal element $\rho_{12}(t)$ of the density matrix in the propagation-eigenstate basis from eq. (83) we then find

[TABLE]

where for simplicity we neglected the small difference between the momentum distribution amplitudes of the WPs of the two propagation eigenstates $f_{\vec{p}_{2}}(\vec{p})$ and $f_{\vec{p}_{1}}(\vec{p})$.

Next, we note that our current assumption of near-constant elements of $\tilde{\cal H}(t)$ in eq. (88) actually implies that $c\equiv\cos\theta$ and $s\equiv\sin\theta$ are practically constant, so that $\dot{\theta}\simeq 0$. Eq. (101) then gives $\varphi\simeq 0$, and from eqs. (83) and (102) we get

[TABLE]

Notice that the assumption $\varphi\simeq 0$ adopted here actually implies $|\dot{\theta}|\ll s_{2}|\dot{\phi}|$, which is opposite to the condition used in Section 3.4.3. Thus, the model in this subsection is in a sense complementary to the one in Section 3.4.3.

We are interested in the question of whether the off-diagonal elements of the neutrino density matrix in the propagation eigenstate basis get suppressed at late times. Consider therefore the regime of very large times, for which the elements of $\rho_{jk}$ take time-independent asymptotic values. It can be seen from eq. (105) that this corresponds to the averaging of the oscillating terms in the integrand of the momentum integral.

Indeed, the width of the momentum interval that contributes significantly to the integral in (105) is given by the width $\sigma_{p}$ of the momentum distribution functions $f_{\vec{p}_{i}}(\vec{p})$. If the phase $\Phi=\Omega(t-t_{0})$ of the oscillating factors in (105) changes significantly over this interval, i.e. if

[TABLE]

the $\Phi$-dependent terms in (105) exhibit fast oscillations and can be replaced by their average values. Therefore the characteristic time after which the asymptotic regime sets in is

[TABLE]

where $\sigma_{x}\sim\sigma_{p}^{-1}$ is the spatial length of the neutrino wave packet. In the adiabatic regime we have $2\Omega=\Delta$, and $t_{\rm char}$ essentially coincides with the adiabatic coherence length $(L_{\rm coh})_{\rm adiab}\equiv|\partial\Delta/\partial\vec{p}\,|^{-1}\sigma_{x}$; however, in the non-adiabatic case the characteristic time $t_{\rm char}$ may differ significantly from the naive adiabatic coherence length $(L_{\rm coh})_{\rm adiab}$.

Consider now $\rho_{12}$ at asymptotically large times. Averaging the $\Phi$-dependent oscillating factors, from eq. (105)) we find

[TABLE]

Let us consider this expression in some limiting cases. In the adiabatic limit $\lambda\to 0$ we have $\tilde{s}_{2}\to 0$ (no mixing between the propagation eigenstates) and, as expected, $\rho_{12}$ vanishes at asymptotically large times due to the WP separation. In the opposite, extremely non-adiabatic, limit $\lambda\to\infty$, we have $\tilde{s}_{2}=1$, which corresponds to maximal mixing between the propagation eigenstates. In this case eq. (108) yields $(\rho_{12})_{as}=s(t_{0})c(t_{0})$, i.e. the asymptotic value of $\rho_{12}(t)$ coincides with its initial value at $t=t_{0}$. The overlap of the WPs of the different propagation eigenstates does not reduce with time, i.e. no decoherence by WP separation occurs. In the intermediate case of a moderate adiabaticity violation, partial decoherence takes place.

These results can be readily understood. In the adiabatic limit ($\lambda\to 0$), neutrino propagation eigenstates evolve independently and propagate with well defined group velocities, which are different for different eigenstates. Therefore, at asymptotically large times their overlap essentially vanishes, leading to vanishing off-diagonal elements of the neutrino density matrix and complete decoherence.

In the opposite (extremely non-adiabatic) limit $\lambda\gg 1$, the propagation eigenstates, though mathematically well defined at any given $t$, are strongly mixed and therefore are not a physically meaningful notion. Their group velocities are then not well defined either, and therefore no WP separation can occur. The WP describing the initially produced neutrino flavour eigenstate propagates as a single conglomerate, without separating into its different propagation-eigenstate components. To put it another way, on a time scale $\tau=\pi/(2\Omega)$ which is short compared to the coherence length expected in the adiabatic regime, the propagation eigenstates convert into each other with maximal amplitude. Each swap interchanges the WPs and so flips the sign of their group velocity difference. As a result, the (small) separation of the WPs over one period $\tau$ is completely compensated after the next such period, and no noticeable wave packet separation ever occurs. Absence of WP separation in the extremely non-adiabatic regime means that no decoherence takes place. The off-diagonal elements of the neutrino density matrix at asymptotic times coincide with their initial values at $t=t_{0}$, in agreement with the conclusions drawn above from eq. (108).

In the case of moderate adiabaticity violation ($\lambda\sim 1$), the mixing between the different propagation eigenstates is sizeable but not maximal. The shuffling of the propagation eigenstates still occurs, but with an amplitude smaller than one. Therefore, only a fraction of the faster eigenstate converts into the slower one and vice versa. As a result, at asymptotically large times the separation of WPs still occurs, but the strengths of the separated WPs are uneven: the probability of finding a neutrino in one of them is larger than in the other.

It should be stressed that the above discussion applies only to the case of continuous adiabaticity violation considered here; in the case when adiabaticity is strongly violated only during a short fraction of the overall neutrino evolution time, the situation is different, as illustrated by the sample model discussed in Section 3.4.3.

4 Summary and discussion

We have considered the effects of decoherence by WP separation on collective neutrino oscillations in dense neutrino backgrounds. In particular, we have studied in detail the question of whether decoherence effects can be described by an extra term in the evolution equation of the density matrix of a neutrino as a whole (as contrasted to the EoM of fixed-momentum components of the density matrix) or of the neutrino density matrix in coordinate space.

Propagation decoherence occurs when the WP separation distance exceeds their length. An important quantity characterizing decoherence effects is therefore the length of WPs of neutrinos produced e.g. in supernovae, where collective neutrino oscillations are expected to occur. We considered this issue in Section 2. Our estimates show that, for neutrinos produced in the processes with participation of nucleons (such as $\beta$-processes or neutrino bremsstrahlung in nucleon collisions), the length $\sigma_{x}$ of the neutrino WPs is dominated by the overlap time of the nucleon WPs and is typically $\lesssim 4\times 10^{-12}$ cm. Interestingly, this estimate is rather close to the value $\sigma_{x}\sim 10^{-11}$ cm obtained in [22, 23] from completely different considerations, based on the assumption that for all production processes with participation of electrons $\sigma_{x}$ is dominated by the length of the electron WP.

In Section 3.3 we derived EoMs satisfied by both the coarse grained (coordinate-averaged) and un-averaged neutrino density matrices in vacuum. We have found that these equations essentially coincide, provided that in the un-averaged case the total time derivative is understood as the Liouville operator, with the average group velocity of the WPs of the neutrino mass eigenstates playing the role of the neutrino transport velocity. Thus, modulo this identification, the system is described by a Liouville-type equation despite the fact that neutrinos are described by quantum WPs rather than being considered as classical pointlike particles.

An important feature of the EoMs of the considered neutrino density matrices in vacuum is that these equations contain a damping term, which leads to suppression with time of the off-diagonal elements of the neutrino density matrix in the mass eigenstate basis. This suppression is a consequence of the spatial separation of the WPs of the different neutrino mass eigenstates composing the produced neutrino flavour state. The presence of a term directly describing WP separation in the EoMs of the coordinate-space or coordinate-averaged neutrino density matrices simplifies studies of propagation decoherence, as in this case it is not necessary to explicitly sum over neutrino momentum modes to check if their dephasing occurs at late times.

Equipped with these results, we then studied the EoMs for the density matrices describing neutrino flavour transitions in dense matter and neutrino backgrounds, both in the adiabatic and in the non-adiabatic regimes. We have shown that, just like for vacuum oscillations, in the adiabatic regime the EoM of the neutrino density matrix contains a damping term which describes decoherence by WP separation at late times. However, no such term in general appears in the non-adiabatic regime, and whether or not propagation decoherence occurs depends on the properties of the system. We have considered two specific models of adiabaticity violation – one with strong short-term non-adiabaticity and another with extended non-adiabaticity of an arbitrary strength. For the first model, we found that WP separation does occur, but the degree of decoherence depends on the time elapsed after the short period of adiabaticity violation rather than on the overall time elapsed since neutrino production. For the second model, whether or not WP separation occurs depends on the extent of adiabaticity violation.

The obtained results allow us to give a quantitative interpretation of the well known result [8, 24] that synchronized oscillations that can occur in dense uniform and isotropic neutrino gases are not affected by decoherence effects provided that the neutrino ‘self-interaction’ parameter $\mu=\sqrt{2}G_{F}n_{\nu}$ is sufficiently large (i.e. the neutrino density is high enough). The interpretation depends on the value of $\mu P_{0}$ relative to the mean ‘energy’ of the neutrino WP, $\omega_{0}=[\Delta m^{2}/(2p)]_{0}$, and its mean width in the $\omega$ variable, $\sigma_{\omega}$ (recall that $P_{0}$ is the initial length of the global flavour spin vector). As is shown in the Appendix, for

[TABLE]

we have $\lambda\sim s_{20}\omega_{0}/(\mu P_{0})\ll 1$, i.e. the neutrino system evolves adiabatically. In this case one could in principle expect decoherence by WP separation, which could destroy the synchronized oscillations. This, however, does not happen. Indeed, eq. (A11) means that at the neutrino production time $t_{0}$ one has $s_{2}(t_{0})\simeq s_{20}\omega_{0}/(\mu P_{0})\ll 1$, that is the mixing in medium is suppressed, and the produced neutrino flavour eigenstate coincides with one of the propagation eigenstates. The produced state therefore consists of just one WP, which cannot ‘separate with itself’. As a result, no decoherence occurs.

If $\mu P_{0}$ satisfies

[TABLE]

the non-adiabaticity parameter $\lambda$ is typically large (provided that $\sigma_{\omega}\ll\omega_{0}$, so that from (110) it follows that $\mu P_{0}\sim\omega$ for all $\omega$ in the neutrino spectrum). This means that the propagation eigenstates are strongly mixed, $\tilde{s}_{2}\simeq 1$. It then follows from eq. (108) that the WP overlap at late times is the same as it is at $t_{0}$, i.e. no decoherence by WP separation takes place. The reason for this was discussed in Section 3.4.2. Strong mixing between the propagation eigenstates means that they are not physically meaningful quantities, and neither are their group velocities; therefore there is no WP separation. Another way of saying this is that the propagation eigenstates go into each other (shuffle) on a very short time scale $\tau$, which precludes any observable WP separation. Thus, coherence is maintained in that case. If, however, $\sigma_{\omega}$ is comparable to $\omega_{0}$, eq. (110) does not guarantee that $\lambda$ is large for all $\omega$-modes, and partial or full decoherence by WP separation takes place at late times, leading to partial or complete de-synchronization of neutrino oscillations.

It has been shown in [8, 24] that synchronized neutrino oscillations can occur (and so no decoherence by WP separation takes place) also in the limit

[TABLE]

In this case $\lambda\sim s_{20}\mu P_{0}/\omega_{0}\ll 1$, i.e. adiabaticity is good, and mixing at neutrino production essentially coincides with vacuum mixing. Therefore, the initially produced neutrino flavour state is again a nontrivial mixture of different propagation eigenstates, and naively one could expect decoherence by WP separation and loss of synchronization at late times. It is not difficult to understand why this actually does not happen. Eq. (111) means $\sigma_{\omega}\lll\omega_{0}$, so that the neutrino spectrum is practically monochromatic. Therefore, all neutrinos oscillate with essentially the same frequency, i.e. experience synchronized oscillations.171717One could still expect decoherence to take place in this case, but only at extremely late times $\sim\sigma_{\omega}^{-1}$. The absence of WP separation can also be readily understood at a formal level by considering the off-diagonal elements of the neutrino density matrix. In this case the spectrum shape function $f_{\vec{p}_{i}}(\vec{p})$, when transformed to the $\omega$-variable, is essentially $\delta(\omega-\omega_{0})$, and no time averaging should be performed. Therefore, instead of eq. (108) one should use eq. (105) or (104). Good adiabaticity means $\tilde{s}_{2}\to 0$, $\tilde{c}_{2}\to 1$, and eq. (105) (or (104)) yields

[TABLE]

which up to the phase factor coincides with $\rho_{12}(t_{0})=s(t_{0})c(t_{0})$.

As shown in the Appendix, with $\mu P_{0}$ decreasing and approaching $\mu_{0}P_{0}\sim\sigma_{\omega}$, the non-adiabaticity parameter $\lambda$ decreases, which leads to partial decoherence by WP separation at late times. For $\mu\leq\mu_{0}$ adiabaticity becomes perfect ($\lambda=0$). In this regime the in-medium mixing at neutrino production is of the order of the vacuum mixing and therefore in general is not suppressed. Thus, unlike in the case (109), the initially produced neutrino flavour eigenstate is a nontrivial superposition of different propagation eigenstates. At sufficiently large times their WPs fully separate, destroying synchronized neutrino oscillations.

We have addressed only one aspect of late-time decoherence of oscillations in dense neutrino gases – decoherence by WP separation. To this end, we considered a simplified model – a uniform and isotropic gas of neutrinos described by identical WPs. In realistic situations, one should also take into account that the neutrino ensemble contains WPs with different centroid momenta and perform a statistical averaging over this ensemble, which will also contribute to decoherence effects. Our evolution equations for the density matrix will take this into account automatically, provided that by the momentum distribution functions $|f(\vec{p})|^{2}$ we understand the overall momentum spectrum (that is, the sum of the spectra of all the individual WPs). Then equations like our eq. (70) (in the adiabatic regime) or (76), (77) (in the general case) will describe the space-time evolution of the system and the loss of coherence in it. In particular, some distance or time (coherence time) after which the coherence is lost can be found. However, they in general cannot be interpreted as the distance or time after which the neutrino WPs separate. Thus, interpretation in terms of WP separation will be lost in this case, or will only be a part of the story.

In conclusion, the results of this paper provide novel analytical tools for describing effects of decoherence by wave packet separation in neutrino oscillations. In particular, in Sections 3.3 and 3.4 we have derived effective equations of motion which describe the evolution of a neutrino wave packet as a whole, as opposed to the evolution of its individual momentum modes. These equations could be used, for instance, in simulations of supernova neutrinos to incorporate decoherence effects without the need to simulate a large number of momentum modes.

While this paper was in the stage of the manuscript preparation, the paper [44] appeared, in which, within the wave packet formalism, the Liouville equation for the neutrino density matrix describing neutrino flavour evolution in matter was derived. An extra term in the equation of motion responsible for decoherence due to wave packet separation was obtained and analyzed. However, the nonadiabatic neutrino evolution and the oscillations in neutrino backgrounds were not considered there. Where the results of the present paper and of [44] overlap, they are in agreement with each other.

Acknowledgments. The authors are grateful to Alexei Smirnov and Irene Tamborra for useful discussions and to Antonino Di Piazza for pointing out a typo in eq. (9).

Appendix A:

The flavour spin formalism and the adiabaticity

                              condition

Consider a uniform and isotropic neutrino gas in the 2-flavour framework. The evolution of such a system is conveniently described within the flavour spin formalism [30]. For an isotropic system, one can use the absolute value of the neutrino momentum $p\equiv|\vec{p}|$ rather than the momentum itself in order to label the kinematic characteristics of neutrinos. It proves to be more convenient, however, to use instead the vacuum oscillation frequency $\omega\equiv\Delta m^{2}/(2p)$. The $2\times 2$ density matrix in flavour space can be expanded in terms of the unit matrix and Pauli matrices as

[TABLE]

Here $n_{\nu\omega}$ is the number density of the $\omega$-mode neutrinos, and $\vec{P}_{\omega}$ is the corresponding flavour spin vector. The 2-flavour effective neutrino Hamiltonian can be similarly expanded as

[TABLE]

where the “Hamiltonian vector” $\vec{H}_{\omega}$ is given by

[TABLE]

Here

[TABLE]

where, as discussed in Section 2, $s_{20}\equiv\sin 2\theta_{0}$, $c_{20}\equiv\cos 2\theta_{0}$, with $\theta_{0}$ being the leptonic mixing angle in vacuum, and $n_{\nu}$ is the total density of neutrinos. The quantity $\vec{P}$ is the global flavour spin vector which is the sum (integral) of the individual flavour spin vectors $\vec{P}_{\omega}$ corresponding to fixed $\omega$:

[TABLE]

Here positive and negative $\omega$ correspond to neutrinos and antineutrinos, respectively. The terms in $\rho_{\omega}$ and ${\cal H}$ that are proportional to the unit matrix do not affect neutrino evolution and will hereafter be omitted, i.e. we will consider only the traceless parts of $\rho_{\omega}$ and ${\cal H}$.

The evolution equation for the density matrix (78) yields the following EoM for the individual flavour spin vectors:

[TABLE]

i.e. the vectors $\vec{P}_{\omega}$ precess around their corresponding $\vec{H}_{\omega}$ with frequencies $H_{\omega}\equiv|\vec{H}_{\omega}|$. Integrating eq. (A6) over $\omega$, we obtain the EoM for the global flavour spin vector:

[TABLE]

From EoMs (A6) and (A7) it follows that the quantities $|\vec{P}_{\omega}|$ and $\vec{B}\vec{P}$ are conserved by neutrino evolution.

The flavour content of a given $\omega$-mode is defined by the projection of $\vec{P}_{\omega}$ onto the 3rd axis in flavour space. We will be assuming that all neutrinos are produced at $t=0$ in the electron flavour, which means that the initial conditions for $\vec{P}_{\omega}$ and $\vec{P}$ are

[TABLE]

where $g_{\omega}$ is the neutrino spectrum in the $\omega$ variable normalized according to $\int\,d\omega g_{\omega}=1$, and $\vec{n}_{3}$ is the unit vector in the third direction in flavour space. Then $|\vec{P}_{\omega}|=P_{0}g_{\omega}$, $\vec{B}\vec{P}=-c_{20}P_{0}$, and the length of the vector $\vec{H}_{\omega}$ is given by

[TABLE]

The traceless Hermitian $2\times 2$ Hamiltonian matrix ${\cal H}(t)$ can be written in terms of the components of the real vector $\vec{H}_{\omega}=(H_{\omega 1},\,H_{\omega 2},\,H_{\omega 3})$ as

[TABLE]

Comparing eqs. (A10) and (86), one can find the quantities $c_{2}\equiv\cos 2\theta(t)$, $s_{2}\equiv\sin 2\theta(t)$ and $\phi(t)$, through which the Hamiltonian ${\cal H}$ in eq. (86) is expressed, in terms of the components of the Hamiltonian vector $\vec{H}_{\omega}$ and $\Delta\equiv E_{2}(t)-E_{1}(t)$:

[TABLE]

Here we used the notation $H_{\omega\perp}\equiv\sqrt{H_{\omega 1}^{2}+H_{\omega 2}^{2}}$. Note that

[TABLE]

The angles $\theta(t)$ and $\phi(t)$ have simple interpretation in the flavour spin formalism: $2\theta$ is the angle between the Hamiltonian vector $\vec{H}_{\omega}$ and the 3rd direction in flavour space, i.e. it is the polar angle of $\vec{H}_{\omega}$, whereas $\phi(t)$ is the azimuthal angle characterizing the direction of the projection of $\vec{H}_{\omega}$ onto the 1–2 plane.

The non-adiabaticity parameter $\lambda$ defined in eq. (89) depends on the derivatives $\dot{\theta}$ and $\dot{\phi}$. Let us express them in terms of the components of $\vec{H}_{\omega}$ and their derivatives. From the definitions in eq. (A11) we find

[TABLE]

It proves convenient to introduce the angular velocity of the vector $\vec{H}_{\omega}$, i.e. the derivative of the unit vector in its direction:

[TABLE]

Eq. (A14) then takes the form $2\dot{\theta}=(H_{\omega}/H_{\omega\perp})(\Omega_{H_{\omega}})_{3}$. It is not difficult to show that the expression under the square root in the numerator of the formula (89) for $\lambda$ is just the squared length of $\vec{\Omega}_{H_{\omega}}$:

[TABLE]

The non-adiabaticity parameter $\lambda$ can therefore be written as

[TABLE]

Note that it is similar (but not identical) to the non-adiabaticity parameter

[TABLE]

introduced in [10]. The difference between the two parameters is due to the fact that they have somewhat different meaning. The one in eq. (A18) is the ratio of the speed $\Omega_{H_{\omega}}$ with which the direction of $\vec{H}_{\omega}$ changes to the angular velocity $H_{\omega}$ of precession of $\vec{P}_{\omega}$ around $\vec{H}_{\omega}$: when this ratio is small, the individual flavour spin vectors, while precessing around their $\vec{H}_{\omega}$, track the motion of the latter. On the other hand, our parameter $\lambda$ is defined as the tangent of twice the mixing angle between the propagation eigenstates; it describes to what extent the independence of the evolution of these eigenstates is violated. Both parameters characterize the degree of adiabaticity violation in the neutrino system.

The evolution of the system under consideration exhibits a threshold behaviour in the neutrino self-interaction parameter $\mu$: there exists a critical value $\mu_{0}$ which depends on the neutrino spectrum $g_{\omega}$ and is of the order of the width of this spectrum $\sigma_{\omega}$ [8, 24]. For $\mu\leq\mu_{0}$ complete decoherence is achieved at asymptotic (late) times, whereas for $\mu>\mu_{0}$ the system maintains partial or complete coherence at all times. For $\mu>\mu_{0}$, the late time evolution of the global flavour spin vector $\vec{P}$ is to a good approximation a simple precession around the vector $\vec{B}$ with a frequency $\omega_{s}$:

[TABLE]

In this regime the length of the vector $\vec{P}$ is (approximately) conserved; however, it may change significantly at short and intermediate times, before the asymptotic regime sets in. For $\mu\gg\mu_{0}$ the EoM in (A19) is practically exact, and the asymptotic regime described by this equation sets in immediately at $t=0$; the length of the vector $P$ is therefore preserved by neutrino evolution, i.e. $P=P_{0}$. The precession frequency $\omega_{s}$ in this case coincides with the spectral average $\omega_{0}\equiv\int g_{\omega}\omega\,d\omega$. The flavour spin vectors $\vec{P}_{\omega}$ of all the individual modes precess around $\vec{B}$ with the same frequency $\omega_{0}$ and in phase with each other, i.e. synchronized oscillations take place.

For $\mu\gtrsim\mu_{0}$, i.e. $\mu$ exceeding $\mu_{0}$ but not by much, the late-time EoM (A19) is only approximate; the length of the component of $\vec{P}$ orthogonal to $\vec{B}$ is actually not exactly conserved, but exhibits small oscillations. Still, averaging over these small oscillations (which leads to eq. (A19)) is a good first approximation. In this case the length of the vector $\vec{P}$ shrinks during the time period before the asymptotic regime sets in, i.e. the asymptotic value of $P$ satisfies $P_{\text{min}}<P<P_{0}$. The existence of a non-zero minimum value $P_{\text{min}}=c_{20}P_{0}$ follows from the conservation of $\vec{B}\vec{P}$. The shrinkage of $P$ is a consequence of partial dephasing of the individual $\omega$-modes. In configuration-space description, it is related to decoherence by WP separation, which, however, is not complete in this case (see Section 3.4.4). The precession frequency $\omega_{s}$ is in general different from $\omega_{0}$.

For $\mu\leq\mu_{0}$, at asymptotically large times the flavour spin vector $\vec{P}$ aligns with $\vec{B}$ and its evolution stops. Its length $P$ takes its minimum possible value $P=c_{20}P_{0}$, which means complete dephasing of the individual $\omega$-modes, or equivalently complete decoherence by WP separation in configuration space.

Let us now consider late-time evolution of $\vec{P}$ for $\mu\geq\mu_{0}$ assuming the EoM in eq. (A19) to be exact. It is not difficult to solve this equation. To fix the initial conditions, we have to specify the value of $P=|\vec{P}|$, as well as the direction of the vector $\vec{P}$ at a given time $t_{0}$. Let us introduce the parameter $R$ describing the length of this vector in the asymptotic regime:

[TABLE]

The length of the component of $\vec{P}$ that is orthogonal to $\vec{B}$ is then given by181818The notation $\vec{P}_{\perp\vec{B}}$ was introduced here to stress the difference between the vector components orthogonal to $\vec{B}$ from those orthogonal to the 3rd axis in flavour space, cf. eq. (A11).

[TABLE]

Note that the quantity $R/s_{20}$ can be considered as the order parameter characterizing decoherence in the system. Indeed, for $\mu\gg\mu_{0}$ we have $P\simeq P_{0}$ and $R\simeq s_{20}$, whereas at the threshold $\mu=\mu_{0}$ we obtain $P=P_{\text{min}}=c_{20}P_{0}$, $R=0$.

We shall adopt the initial condition such that at a time $t_{0}$ the vector $\vec{P}$ lies in the 1–3 plane in flavour space. Direct integration of eq. (A19) then yields

[TABLE]

where

[TABLE]

The components of the vector $\vec{H}_{\omega}$ can then be found by substituting (A22) into (A3):

[TABLE]

It should be remembered that in the above expressions for the components of $\vec{H}_{\omega}$ the frequency $\omega$ in general does not coincide with $\omega_{s}$, though for narrow neutrino spectra ($\sigma_{\omega}\ll\omega_{s}$) they must be close to each other.

For the lengths of $\vec{H}_{\omega}$ and of its component $\vec{H}_{\omega\perp}$ orthogonal to the 3rd axis in flavour space we then find

[TABLE]

[TABLE]

Note that eq. (A21) means that (A25) is consistent with eq. (A9).

The mixing angle in neutrino background $\theta(t)$ and the angle $\phi(t)$ characterizing the projection of $\vec{H}_{\omega}$ onto the 1–2 plane in flavour space are defined through the components of the vector $\vec{H}_{\omega}$ by the relations in eq. (A11). Substituting there the expressions for these components from eq. (A24) we obtain, in the case when the asymptotic evolution of $\vec{P}$ is described by eq. (A19),

[TABLE]

[TABLE]

In calculating $2\dot{\theta}$ we have taken into account that eqs. (A3) and (A19) imply that $\vec{H}_{\omega}$ satisfies $\dot{\vec{H}}_{\omega}=\omega_{s}\vec{B}\times\vec{H}_{\omega}$, i.e. in the regime described by EoM (A19), $\vec{H}_{\omega}^{2}=\text{const.}$ and $\dot{\vec{H}}_{\omega}\vec{H}_{\omega}=0$. Note that, while $\dot{\theta}$ oscillates around zero or small average, taking both positive and negative values, $\dot{\phi}$ oscillates around a non-zero value, remaining for most of the parameter space sign-definite.

Let us calculate the numerator of expression (89) for $\lambda$ in the regime described by EoM (A19). Making use of eqs. (A28) and (A27), we find

[TABLE]

which is time-independent. For the non-adiabaticity parameter $\lambda$ we find

[TABLE]

Consider first the limit $\mu P_{0}\gg\mu_{0}P_{0}\sim\sigma_{\omega}$. In this case $P\to P_{0}$, $R\to s_{20}$, $\omega_{s}\to\omega_{0}$. From eqs. (A28), (A11) and (A9) it then follows that one can neglect $c_{2}\dot{\phi}$ compared to $H_{\omega}$ in eq. (A30) provided that $\mu P_{0}\gg\omega$ or $\mu P_{0}\ll\omega$. This gives $\lambda\sim s_{20}\omega_{0}\mu P_{0}/H_{\omega}^{2}$. For $\mu P_{0}\gg\omega_{0}$ we have $H_{\omega}\approx\mu P_{0}$ and $\lambda\sim s_{20}\omega_{0}/(\mu P_{0})\ll 1$, whereas for $\mu P_{0}\ll\omega_{0}$ we obtain $H_{\omega}\approx\omega_{0}$ and $\lambda\sim s_{20}\mu P_{0}/\omega_{0}\ll 1$. Thus, in the limits when $\mu P_{0}$ is vastly different from $\omega_{0}$, neutrino evolution is adiabatic. If $\mu P_{0}\sim\omega_{0}$, one cannot in general neglect $c_{2}\dot{\phi}$ compared to $H_{\omega}$ in the denominator in eq. (A30). In this case $\lambda$ is either $\sim 1$ or $\gg 1$, i.e. adiabaticity is either moderately or strongly violated.

For $\mu$ approaching $\mu_{0}$ from above, we have $R\to 0$, so that $\lambda\to 0$, and neutrino evolution at late times is fully adiabatic. The same is also true for $\mu<\mu_{0}$.