Fast Neutrino Flavor Conversion: Collective Motion vs. Decoherence

TL;DR

This paper investigates the dynamics of neutrino flavor conversion, revealing the coexistence of collective and non-collective modes, and identifying conditions that lead to tachyonic instabilities in dense neutrino gases.

Contribution

It provides a detailed analysis of the dispersion relation for fast neutrino flavor modes, including explicit eigenfunctions and conditions for instabilities based on angular mode crossings.

Findings

Existence of a continuum of non-collective modes alongside discrete collective modes.

Single angular crossings lead to tachyonic instabilities, while multiple crossings may not.

Explicit eigenfunctions for both collective and non-collective modes are derived.

Abstract

In an interacting neutrino gas, flavor coherence becomes dynamical and can propagate as a collective mode. In particular, tachyonic instabilities can appear, leading to "fast flavor conversion" that is independent of neutrino masses and mixing angles. On the other hand, without neutrino-neutrino interaction, a prepared wave packet of flavor coherence simply dissipates by kinematical decoherence of infinitely many non-collective modes. We reexamine the dispersion relation for fast flavor modes and show that for any wavenumber,there exists a continuum of non-collective modes besides a few discrete collective ones. So for any initial wave packet, both decoherence and collective motion occurs, although the latter typically dominates for a sufficiently dense gas. We derive explicit eigenfunctions for both collective and non-collective modes. If the angular mode distribution of electron-lepton number crosses between positive and negative values, two non-collective modes can merge to become a tachyonic collective mode. We explicitly calculate the interaction strength for this critical point. As a corollary we find that a single crossing always leads to a tachyonic instability. For an even number of crossings, no instability needs to occur.

Figures

Equations

$$({\partial }_t+\mathbf{v}\cdot {\partial }_{\mathbf{r}}){\varrho }_{\mathbf{p}}=-i[{\mathsf{H}}_{\mathbf{p}},{\varrho }_{\mathbf{p}}]+\mathcal{C}({\varrho }_{\mathbf{p}}),$$

$${\mathsf{H}}_{\mathbf{p}}=\sqrt{2}{G}_{\mathrm{F}}\int \frac{d^3{\mathbf{p}}^{\prime }}{(2\pi {)}^3}(1-\mathbf{v}\cdot {\mathbf{v}}^{\prime })({\varrho }_{{\mathbf{p}}^{\prime }}-{\bar{\varrho }}_{{\mathbf{p}}^{\prime }}),$$

$${\varrho }_{\mathbf{p}}=\frac{f_{{\nu }_e,\mathbf{p}}+f_{{\nu }_x,\mathbf{p}}}{2}\text{mathbbm}1+\frac{f_{{\nu }_e,\mathbf{p}}-f_{{\nu }_x,\mathbf{p}}}{2}\left(\begin{array}{cc}{s}_{\mathbf{p}}& S_{\mathbf{p}}\\ S_{\mathbf{p}}^{\ast }& -s_{\mathbf{p}}\end{array}\right),$$

$$i({\partial }_t+\mathbf{v}\cdot {\partial }_{\mathbf{r}})S_{\mathbf{p}}=-\sqrt{2}{G}_{\mathrm{F}}\int \frac{d^3\mathbf{p}}{(2\pi {)}^3}(1-\mathbf{v}\cdot {\mathbf{v}}^{\prime })\left(S_{{\mathbf{p}}^{\prime }}{f}_{{\mathbf{p}}^{\prime }}-{\bar{S}}_{{\mathbf{p}}^{\prime }}{\bar{f}}_{{\mathbf{p}}^{\prime }}\right).$$

$$i({\partial }_t+\mathbf{v}\cdot {\partial }_{\mathbf{r}})S_{\mathbf{v}}=-\sqrt{2}{G}_{\mathrm{F}}\int \frac{d{\mathbf{v}}^{\prime }}{4\pi }(1-\mathbf{v}\cdot {\mathbf{v}}^{\prime })g_{{\mathbf{v}}^{\prime }}{S}_{{\mathbf{v}}^{\prime }}.$$

$$g_{\bf v}=\int_{0}^{\infty}\frac{E^{2}dE}{2\pi^{2}}\,\bigl{(}f_{\nu_{e},{\bf p}}-f_{\bar{\nu}_{e},{\bf p}}-f_{\nu_{x},{\bf p}}+f_{\bar{\nu}_{x},{\bf p}}\bigr{)}$$

$$i({\partial }_t+u{\partial }_r)S_u=-\mu {\int }_{-1}^{+1}d{u}^{\prime }(1-u{u}^{\prime })G_{u^{\prime }}{S}_{u^{\prime }},$$

$$G_u$$

$$\mu$$

$$S_u(t,r)=Q_u(\Omega ,K)e^{-i(\Omega t-Kr)},$$

$$(\Omega -uK)Q_u(\Omega ,K)=-\mu {\int }_{-1}^{+1}d{u}^{\prime }(1-u{u}^{\prime })G_{u^{\prime }}{Q}_{u^{\prime }}(\Omega ,K).$$

$$Q_u=\frac{a+bu}{w-u}$$

$$a+bu=-\frac{\mu }{K}{\int }_{-1}^{+1}d{u}^{\prime }\frac{(1-u{u}^{\prime })(a+b{u}^{\prime })}{w-u^{\prime }}{G}_{u^{\prime }}.$$

$$\left(\begin{array}{cc}+K& \\ & -K\end{array}\right)\left(\begin{array}{c}a\\ b\end{array}\right)=-\mu \left[{\int }_{-1}^{+1}du\frac{G_u}{w-u}\left(\begin{array}{cc}1& u\\ u& u^2\end{array}\right)\right]\left(\begin{array}{c}a\\ b\end{array}\right).$$

$$\Vert \left(\begin{array}{cc}K& \\ & -K\end{array}\right)+\mu \left(\begin{array}{cc}⟨1{⟩}_w& ⟨u{⟩}_w\\ ⟨u{⟩}_w& ⟨u^2{⟩}_w\end{array}\right)\Vert =0,$$

$$⟨u^n{⟩}_w={\int }_{-1}^{+1}du{G}_u\frac{u^n}{w-u}.$$

$$K_w=-\mu \frac{{\kappa }_w\pm \sqrt{{\Delta }_w}}{2},$$

$${\mathcal{G}}_n={\int }_{-1}^{+1}du{G}_u{u}^n.$$

$${\Omega }_{w\to \pm \infty }=-\mu \frac{{\mathcal{G}}_0-{\mathcal{G}}_2\pm \sqrt{{\left({\mathcal{G}}_0+{\mathcal{G}}_2\right)}^2-4{\mathcal{G}}_1^2}}{2}.$$

$$(w-u)Q_{w,u}=-\tilde{\mu }{\int }_{-1}^{+1}d{u}^{\prime }(1-u{u}^{\prime })G_{u^{\prime }}{Q}_{w,u^{\prime }},$$

$$i{\partial }_{t}S=\mathcal{H}(S),$$

$${\cal H}(S)\big{|}_{u}=u\,S_{u}(t)-\tilde{\mu}\int_{-1}^{+1}du^{\prime}\,(1-uu^{\prime})\,G_{u^{\prime}}\,S_{u^{\prime}}(t)\,.$$

$$\tilde{\mu }=-\frac{2}{{\kappa }_w\pm \sqrt{{\Delta }_w}},$$

$$Q_{w,u}=A_w\left[\frac{1}{w-u}+{\alpha }_w\delta (w-u)\right]+B_w$$

$$A_w+B_w(w-u)=-\tilde{\mu }{A}_w\left[\text{vbox}-{\int }_{-1}^{+1}d{u}^{\prime }\frac{1-u{u}^{\prime }}{w-u^{\prime }}{G}_{u^{\prime }}+{\alpha }_w{G}_w(1-uw)\right]-\tilde{\mu }{B}_w{\int }_{-1}^{+1}d{u}^{\prime }(1-u{u}^{\prime })G_{u^{\prime }}.$$

$$\text{vbox}-{\int }_{-1}^{+1}d{u}^{\prime }\frac{1-u{u}^{\prime }}{w-u^{\prime }}{G}_{u^{\prime }}=(1-uw)\underset{F_w}{\underbrace{\text{vbox}-{\int }_{-1}^{+1}d{u}^{\prime }\frac{G_{u^{\prime }}}{w-u^{\prime }}}}+u\underset{{\mathcal{G}}_0}{\underbrace{{\int }_{-1}^{+1}d{u}^{\prime }{G}_{u^{\prime }}}}.$$

$$A_w+B_w(w-u)+\tilde{\mu }{A}_w\left[(1-uw)(F_w+{\alpha }_w{G}_w)+u{\mathcal{G}}_0\right]+\tilde{\mu }{B}_w\left({\mathcal{G}}_0-u{\mathcal{G}}_1\right)=0.$$

$$\left(\begin{array}{cc}\tilde{\mu }(F_w+{\alpha }_w{G}_w)+1& \tilde{\mu }{\mathcal{G}}_0+w\\ \tilde{\mu }(F_w+{\alpha }_w{G}_w)w-\tilde{\mu }{\mathcal{G}}_0& \tilde{\mu }{\mathcal{G}}_1+1\end{array}\right)\left(\begin{array}{c}{A}_w\\ B_w\end{array}\right)=0.$$

$${\alpha }_w{G}_w=-F_w-\frac{1+\tilde{\mu }{\mathcal{G}}_0(\tilde{\mu }{\mathcal{G}}_0+w)+\tilde{\mu }{\mathcal{G}}_1}{\tilde{\mu }(1-w^2-\tilde{\mu }w{\mathcal{G}}_0+\tilde{\mu }{\mathcal{G}}_1)}.$$

$${\beta }_w\equiv \frac{B_w}{A_w}=\frac{w+\tilde{\mu }{\mathcal{G}}_0}{1-w^2-\tilde{\mu }(w{\mathcal{G}}_0-{\mathcal{G}}_1)},$$

$${\alpha }_w$$

Full text

1 Introduction

Beginning fifty years ago with the “solar neutrino problem” of an apparent $\nu_{e}$ flux deficit from the Sun and twenty years ago with the up-down asymmetry of atmospheric neutrinos, the issue of flavor conversion has developed to a main-stream field of experimental and theoretical investigation. One intruiging ingredient is the refractive effect on neutrino propagation by a background medium [1] that can lead to large flavor conversion even for small mixing angles (MSW effect [2]). It remains forever fascinating that for many cases of interest, the vacuum oscillation frequency $\Delta m^{2}/2E$ and the matter potential $\sqrt{2}G_{\rm F}n_{e}$ are similar, one depending on the small neutrino mass differences and the other on their weak interaction strength.

One theoretical frontier that remains somewhat unsettled is neutrino propagation in a medium consisting of other neutrinos [3, 4, 5, 6, 7, 8, 9]. In particular, this question concerns collapsing stars, often leading to core-collapse supernovae, and neutron-star mergers. Neutrinos play a dominant dynamical role in such events and their flavor composition likely shapes nucleosynthesis in the concomitant neutrino-driven matter outflows. Moreover, interpreting the flavor composition of the neutrino signal from the next nearby supernova will depend on a better theoretical understanding of neutrino-neutrino refraction in the source region.

In the mean-field approximation, neutrino radiation is represented by density matrices $\varrho_{{\bf p}}(t,{\bf x})$ in flavor space for each momentum mode ${\bf p}$, where the diagonal entries are the usual occupation numbers, while the off-diagonal elements encode flavor coherence. In the ultrarelativistic limit, the space-time evolution follows from the kinetic equation [10, 11, 12, 13, 14, 15, 16, 17, 18, 19]

[TABLE]

where the neutrino velocity vector ${\bf v}={\bf p}/|{\bf p}|$ is taken to be a unit vector describing the direction of motion. The left-hand side (lhs) is the advection term provided by the Liouville operator and describes free streaming, assuming we ignore momentum changes by coherent forces and notably gravitational bending of trajectories. The collision term on the right-hand side (rhs) accounts for scattering and emission or absorption. Coherent flavor evolution is governed the Hamiltonian matrix ${\sf H}_{\bf p}$ that depends on neutrino masses and the flux densities of background particles including other neutrinos.

Ignoring coherent flavor conversion means replacing ${\sf H}_{\bf p}$ by the unit matrix, in which case flavor evolution can arise only by collisions such as pair annihilation and charged-current processes that are encoded in the collision term. In numerical supernova (SN) simulations, this assumption is always used and justified by the large matter effect which effectively demixes neutrinos, i.e., propagation eigenstates are very nearly identical to flavor eigenstates. Coherent flavor evolution, e.g. by MSW conversion, has usually been studied by means of post-processing. It remains to be understood if it is necessary to incorporate coherent flavor evolution in numerical SN simulations to obtain reliable results regarding the neutrino-driven explosion mechanism and nucleosynthesis.

One recent focus of attention concerns “fast flavor conversion,” i.e., nontrivial dynamical solutions of equation (1.1) that do not depend on the vacuum oscillation frequency and would exist even for vanishing neutrino masses. Some initial off-diagonal seed in the $\varrho_{{\bf p}}(t,{\bf x})$ distribution could lead to self-induced flavor conversion [20, 21, 22, 23, 24, 25, 26]. The relevant length scale would be governed by a typical neutrino-neutrino refractive energy $\mu=\sqrt{2}G_{\rm F}n_{\nu\bar{\nu}}$ that far exceeds the vacuum oscillation frequency. So these effects could be important on short length scales and modify neutrino flavor evolution near the neutrino decoupling region in a SN core or in neutron-star binary mergers.

Several questions need answering. What are the required conditions for the neutrino gas to support fast collective modes? If these conditions are satisfied, what exactly triggers these modes? And if they are started, will unstable modes indeed reach the nonlinear phase and cause tangible effects relevant for the physics of compact objects?

The first of these questions can be nicely addressed in terms of a normal-mode analysis of the linearised version of the kinetic equation [27, 28, 29, 30, 31, 32, 33, 34, 35]. The crucial ingredient for supporting possible fast tachyonic instabilities appears to be the angle distribution of the electron lepton number (ELN) carried by neutrinos that seems to require a “crossing” between positive and negative values, although the concept of the ELN distribution would need to be formulated more generally if the other flavors, notably the muon neutrinos, also carry a flavor lepton number flux. Several groups have studied the ELN angle distribution in various types of SN and binary merger simulations to look for astrophysical environments that would show such crossings and thus would support fast-flavor conversion [36, 37, 38, 39, 40, 41].

We here return to more basic questions of the normal-mode analysis for the linearised version of equation (1.1) to clarify the emergence of collective modes from the multitude of non-collective ones that exist in the absence of neutrino refraction and otherwise co-exist with collective modes. In the “fast flavor limit” of vanishing neutrino masses, and ignoring neutrino-neutrino refraction, any prepared initial condition of flavor coherence would be dissipated by the multi-directional motion of the various neutrino modes given by the advection term. In other words, the dispersion relation corresponding to the kinetic equation should always include all modes, collective and non-collective ones, even if it is only the former that ultimately are of interest for fast flavor conversion.

To this end we begin in section 2 with the linearised equation of motion and reduce its dimensionality to an axially symmetric system that offers the simplest sandbox for our demonstration. In section 3 we explicitly consider the non-collective modes and derive the explicit form of their singular eigenfunctions, notably in the presence of nonvanishing neutrino-neutrino refraction and we show how collective modes emerge from the continuum of non-collective ones. As a somewhat unexpected bonus, the explicit identification of the non-collective modes allows us in section 4 to identify the critical points where a crossed ELN angle spectrum produces a branch point for a tachyonic solution to appear from the coalescence of two non-collective modes. In section 5 we return to our original goal and show, in a simple example, how collective motion takes over from the non-collective modes with increasing interaction strength. Section 6 is finally given over to conclusions.

2 Linearised equation of motion

2.1 Kinetic equation for fast flavor modes

The starting point for our study is the kinetic equation 1.1 where we neglect the collision term and work in the “fast flavor limit” where neutrino masses and mixing are ignored. Moreover, the background medium is taken to be homogeneous, isotropic, and stationary, in which case the matter effect can be “rotated away.” Therefore, the Hamiltonian matrix on the right-hand side (rhs) includes only neutrino-neutrino interactions and thus has the form

[TABLE]

where $\bar{\varrho}_{{\bf p}}$ is the density matrix for antineutrinos.

We linearise the kinetic equation, so we can limit our discussion to a two-flavor system consisting of $\nu_{e}$ and some other flavor $\nu_{x}$. We write the density matrices as

[TABLE]

where $s_{\bf p}$ is real, $S_{\bf p}$ is complex, and $s_{\bf p}^{2}+|S_{\bf p}|^{2}=1$. To linear order $s_{\bf p}=1$, so we ask only for the space-time evolution of $S_{\bf p}$ which holds the information about flavor coherence. Rotating away also the diagonal neutrino-neutrino matter effect, the linearised EOM is

[TABLE]

An analogous equation applies to the antineutrino flavor coherence $\bar{S}_{\bf p}$. Here we use the difference spectra $f_{\bf p}=(f_{\nu_{e},{\bf p}}-f_{\nu_{\mu},{\bf p}})$ and $\bar{f}_{\bf p}=(f_{\bar{\nu}_{e},{\bf p}}-f_{\bar{\nu}_{\mu},{\bf p}})$.

The EOM does not depend on neutrino energy and is the same for neutrinos and antineutrinos. Therefore, assuming identical initial conditions, we may integrate over energies and sum over neutrinos and antineutrinos to obtain an EOM that depends only on the direction ${\bf v}$ of a given mode,

[TABLE]

The integration is over all directions ${\bf v}$, i.e., over the unit sphere in the space of velocities. The angle distribution of the effective density of electron-lepton number (ELN) is

[TABLE]

with ${\bf p}=E\,{\bf v}$.

2.2 Axial symmetry

We study a restricted class of backgrounds and solutions, where the neutrino angle distribution is axially symmetric relative to some direction ${\bf r}$ that could be the radial direction in the supernova context. Moreover, we consider only those flavor modes that have wave vectors along that same direction, and we assume that the solution itself is axially symmetric, ignoring those modes that spontaneously break axial symmetry. So finally we consider the 1+1 dimensional problem

[TABLE]

where $u=\cos\theta$ is the velocity component along the symmetry direction ${\bf r}$, and the ELN angular distribution and effective interaction strength are

[TABLE]

In this way, $G_{u}$ is a dimensionless function with values of the order of unity, whereas $\mu$ (units of energy) is an effective interaction strength between neutrinos. We have arbitrarily normalised these quantities to the sum of the $\nu_{e}$ and $\bar{\nu}_{e}$ densities, but of course other definitions are possible. There is no entirely natural way to normalise the dimensionless angle distribution.

2.3 Dispersion relation

As equation (2.6) is linear, we can solve it in Fourier space and consider solutions of the form

[TABLE]

where the eigenfunction obeys the equation

[TABLE]

For a given angle distribution $G_{u}$ and interaction strength $\mu$ we can solve this equation and find, for given $K$, the corresponding $\Omega$ as well as the eigenfunction $Q_{u}$ for $(\Omega,K)$.

The solution is trivial in the absence of interactions ($\mu=0$). For given $K$, any $\Omega=wK$ with $-1\leq w\leq+1$ is a solution with the eigenfunctions $Q_{u}(\Omega,K)=\delta(u-w)$. Here $w=\Omega/K$ has the interpretation of the phase velocity for the given wave, which here is less than the speed of light, commensurate with the picture that these modes describe the flavor coherence of a single neutrino mode (or rather a cone of modes) $u=\cos\theta$ relative to the radial direction in a supernova and thus travels with $-1\leq u\leq+1$ along that direction. If we were to set up a wave, or a wave packet, over many modes, it would quickly dissipate by decoherence because of the continuum of frequencies with which it would oscillate. Notice that these modes are not dynamical—the oscillation with frequency $\Omega$ seen by an observer in the laboratory frame arises because a wave with wave number $K$ simply drifts by. This is a purely kinematical effect. One focus of our paper is to understand what happens to these modes once the interaction is turned on.

First, however, we recapitulate the dispersion relation for collective modes that appear in the presence of a nonvanishing $\mu$ and that involve all directions (or cones) $u$ simultaneously, in contrast to the kinematical modes just described. To this end we observe that the rhs of equation (2.9) is of the form $a+bu$ with unknown coefficients $a$ and $b$. Therefore, the eigenfunctions must have the form

[TABLE]

so long as $w=\Omega/K$ is outside of the interval $[-1,+1]$ and thus the denominator does not become singular for $u$ within this interval. In other words, for given $K$, this ansatz applies to real $\Omega$ “outside of the light cone” with phase velocity $|w|>1$ or to unstable modes where $\Omega$ has an imaginary part. Inserting this ansatz on both sides in equation (2.9) yields

[TABLE]

This equality must apply for any $u$ and thus represents two equations, one consisting of the terms that are independent of $u$ and the other linear in $u$. Therefore, after dropping the prime in the integration variable $u^{\prime}$, one finds

[TABLE]

This equation has nontrivial solutions for $a$ and $b$ if the determinant vanishes,

[TABLE]

where we use the notation

[TABLE]

This is a quadratic equation with two solutions

[TABLE]

where \kappa_{w}=\bigl{\langle}1-u^{2}\bigr{\rangle}_{w} and \Delta_{w}=\bigl{\langle}(1-u)^{2}\bigr{\rangle}_{w}\,\bigl{\langle}(1+u)^{2}\bigr{\rangle}_{w}. For those values of $w$ where $\Delta_{w}\geq 0$, we find the parametric solutions $(\Omega,K)=(wK_{w},K_{w})$ of propagating waves.

In the limit $w\to\pm\infty$, where we can neglect $u$ in the denominator of equation (2.14), we find $\langle u^{n}\rangle_{w}\to{\cal G}_{n}/w$, where we use the moments of the $G$-distribution

[TABLE]

In the $w\to\infty$ limit, $\kappa_{w}$ and $\Delta_{w}$ scale with $1/w$ so that $K_{w}\to 0$, i.e., this is the case of vanishing wave number. The corresponding $\Omega_{w}$ involve an extra power of $w$ and become

[TABLE]

We find the same result if we consider the homogeneous case $K=0$ from the start and solve the eigenvalue equation for $\Omega$, as it was done previously in reference [42].

If $G_{u}$ is either positive or negative on the interval $-1\leq u\leq+1$ (no crossing between positive and negative values) and for $|w|>1$, the averages in the definition of $\Delta_{w}$ have equal signs so that $\Delta_{w}\geq 0$. Therefore, in the absence of $G_{u}$-crossings we find two real branches of the dispersion relation and no tachyonic instabilities. For the simplest case of an isotropic medium with $G_{u}=1$ we show the two branches of the collective dispersion relation in the left panel of figure 1 as thick blue lines. Under the light cone (gray shaded region) we show in addition a grid of non-collective modes, assuming they are not modified by the interaction. They provide dense coverage under the light cone ($-1<w<+1$).

We will find it useful to display the same information in the form of the right panel of figure 1, where on the vertical axis we show the phase velocity $w=\Omega/K$ for a given mode as a function of $\tilde{\mu}=\mu/K$ on the horizontal axis. The light cone is now simply the shown gray horizontal band and the non-collective modes are horizontal lines, assuming they are not modified by the interaction. The case of vanishing interaction $\tilde{\mu}=0$ is at the center of the plot. It corresponds to the asymptotic solution $K\to\pm\infty$ in the left panel. While $\mu$ is defined as a positive parameter, $\tilde{\mu}$ can be both positive and negative.

2.4 Schrödinger-like equation

The representation of the right panel of figure 1 corresponds to a somewhat different perspective on the dispersion relation. The only dimensionful parameter of our problem is the interaction strength $\mu$, so all frequencies and wave numbers can be expressed in units of $\mu$. However, we can also imagine that we pick one fixed wave number $K$ and ask for the spectrum of eigenfrequencies $\Omega$, depending on the chosen $\mu$, and express the eigenfrequencies in units of $K$. In this way, what was the phase velocity $w=\Omega/K$ can be seen as the eigenfrequency in units of $K$ and we ask how it behaves as a function of the interaction strength $\tilde{\mu}$. In this picture, equation (2.9) is

[TABLE]

where the chosen wave number is hidden in $\tilde{\mu}$. Returning to an equation for the time evolution, this corresponds to

[TABLE]

where dimensionless time is measured in units of $1/|K|$, but otherwise $K$ does not appear. Here $S$ is a time-dependent complex function of $u$ on the interval $-1\leq u\leq+1$ and the linear operator on the rhs is

[TABLE]

So the right panel of figure 1 shows the eigenvalues of the operator ${\cal H}$ as function of $\tilde{\mu}$. For vanishing interaction, the spectrum of eigenvalues is the continuum $-1\leq w\leq+1$, whereas for nonvanishing $\tilde{\mu}$, in addition two discrete collective modes emerge.

We can think of equation (2.19) as a Schrödinger equation with Hamiltonian ${\cal H}$ and eigenvalues $w$. However, unless $G_{u}={\rm const.}$, the real operator ${\cal H}$ is not symmetric and thus not Hermitean and indeed can have complex eigenvalues. The latter are, of course, the main interest concerning collective modes. Our linearised EOMs are not complete—we have ignored the evolution of $s_{{\bf p}}$, so the evolution of our sub-system need not be unitary.

We should solve the eigenvalue equation and find the energies $w$ as a function of the interaction strength $\tilde{\mu}$. For real $w$ we can turn this question around and find

[TABLE]

where we have used the results of section 2.3. It is then straightforward to find the coefficients $a$ and $b$ for the normalised eigenfunctions of the form (2.10), although the expressions get cumbersome.

3 Non-collective modes

3.1 Eigenfunctions

To determine the non-collective modes (the propagating ones within the light cone), we expect that for any $\tilde{\mu}$ the spectrum of eigenfrequencies densely fills the interval $-1<w<+1$. So for a given $\tilde{\mu}$ we know the spectrum of non-collective eigenvalues and all that is missing is the form of the eigenfunctions. In the noninteracting case they are proportional to $\delta(w-u)$, so switching on a nonvanishing $\tilde{\mu}$ should also involve a $\delta$ function. On the other hand, it should also have the form of equation (2.10) by the same logic that the rhs of equation (2.18) is of the form $a+bu$. If we add $-bw+bw$ in the numerator of equation (2.10) we see that it is $(a+bw)/(w-u)+b$ or $A/(w-u)+B$ in terms of other parameters $A$ and $B$. So overall we propose that the eigenfunctions should be of the form

[TABLE]

with the stipulation that in equation (2.18) the principal value of the integral should be taken. Or turning this point around, the integral by itself is not defined and can be made meaningful by the principal-part prescription up to an unknown constant that must be fixed by the unknown coefficient $\alpha_{w}$.

Inserting this ansatz on both sides of equation (2.18) and using $(w-u)\delta(w-u)=0$, analogous to the noninteracting case, we find the condition

[TABLE]

We use $(1-uu^{\prime})/(w-u^{\prime})=(1-uw)/(w-u^{\prime})+u$ under the integral to write

[TABLE]

The first integral on the rhs is a certain transform of $G(u)$ on the interval $-1\leq u\leq+1$ which we call $F(u)$ and is here to be taken at the location $w$, whereas the second integral is the 0th moment of $G_{u}$ as defined in equation (2.16). So equation (3.2) is

[TABLE]

By the same argument that was used earlier, this equation represents two independent equations for the coefficients $A_{w}$ and $B_{w}$, one from the terms that do not depend on $u$ and one from those that are linear in $u$,

[TABLE]

For any nontrivial solution, the determinant of the matrix must vanish, implying

[TABLE]

Inserting this solution in the matrix, one easily finds

[TABLE]

leaving open the overall normalisation.

The simplest example is an isotropic system with $G_{u}=1$, implying ${\cal G}_{0}=2$, ${\cal G}_{1}=0$, and $F_{u}=\log[(1+u)/(1-u)]$. In this case the explicit coefficients are

[TABLE]

In general, the coefficient $\alpha_{w}$ is fixed by equation (3.6), except for the special case when $G_{w}=0$, i.e., at an angular crossing. In this case, $\delta(w-u)$ is the eigenfunction, so our ansatz would not apply. We will see later that in this case, for certain values of $\tilde{\mu}$ a tachyonic solution branches off.

3.2 Normalisation

These eigenfunctions cannot be normalised, whereas the collective modes with $|w|>1$ have eigenfunctions that can be normalised. This is analogous to free vs. atomic bound electrons. One solution is to confine the electron to some large box instead of infinite space. In our case, we may use discrete angles $u_{i}=\cos\theta_{i}$ with equal spacing $\Delta u$ for the neutrino directions, which is equivalent to using some large box for the neutrinos in the direction of our symmetry axis, leading to quantized neutrino momenta in that direction. For arbitrarily small $\Delta u$, this means $\int_{-1}^{+1}du\,f(u)=\sum_{i=1}^{N}f_{i}\,\Delta u$ for some function $f(u)$ with $f_{i}=f(u_{i})$. For the function $\delta(u)$ with $\int du\,\delta(u)=1$, this implies $\int du\,|\delta(u)|^{2}=1/\Delta u$ because the $\delta$ function becomes a Kronecker $\delta$ at one grid point where it must be represented by the function value $1/\Delta u$ to provide unit sum.

Next we consider the function $1/(w-u)$. Its principal part integrates to a finite value, whereas its normalisation $\int_{-1}^{+1}du\,(w-u)^{-2}$ diverges at $u=w$. To realise the principal-part prescription in the discrete case, we need to fix the grid points $u_{i}$ symmetrically around $w$. For illustration we use $w=0$ and $u_{i}=({\textstyle\frac{1}{2}}\pm i)\Delta u$ with $i$ the integer grid-point indices. Noting that $\sum_{i=-\infty}^{+\infty}({\textstyle\frac{1}{2}}-i)^{-2}=\pi^{2}$ we find that $\int_{-1}^{+1}du\,(w-u)^{-2}=\pi^{2}/\Delta u$.

Finally we consider the normalisation of $Q_{w,u}$ given in equation (3.1). We notice that in $\int du\,|Q_{w,u}|^{2}$ the mixed terms either vanish or are finite, leaving us with

[TABLE]

So $B_{\omega}$ does not contribute to the normalisation and if all eigenfunctions are to be normalised in the same way, we may fix $A_{w}=a_{w}\,\sqrt{\Delta u}$ with

[TABLE]

Of course, the overall sign (or phase) of the eigenfunctions is not fixed and we introduced $s_{w}$ for later convenience. Notice that $B_{w}$ does not vanish even though it does not contribute to the overall normalisation of the singular eigenfunctions. Rather $B_{w}$ is derived from $A_{w}$ through the ratio $\beta_{w}$ given in equation (3.7). We may write it in the form $B_{w}=b_{w}\sqrt{\Delta u}$ and with equation (3.10) one finds

[TABLE]

to achieve a common normalisation for all eigenfunctions.

We finally notice that $\alpha_{w}$ measures the relative contributions of the functions $\delta(w-u)$ and $1/(w-u)$ to the eigenfunction, which thus may be written in the form

[TABLE]

where $A$ is now a chosen global normalisation factor, for example $A=\sqrt{\Delta u}$, that does not depend on the eigenvalue $w$. The angle $\varphi_{w}$ is given by the relations

[TABLE]

For vanishing interaction strength, one finds $\alpha_{w}\to-\infty$, $\varphi_{w}=0$, and $B_{w}=0$ so that we have only the $\delta$ function.

For an isotropic system we show $\varphi_{w}$ and $b_{w}$ for several values of $\tilde{\mu}$ in figure 2. The mixing angle varies in the range $-\pi<\varphi_{w}<+\pi$ on the interval $-1<w<+1$. In the isotropic case, the Hamiltonian in equation (2.20) is symmetric and thus Hermitean, so the eigenfunctions must be orthogonal. It is a simple exercise to show that indeed \mathchoice{{\vbox{\hbox{\textstyle-}}\kern-4.86108pt}}{{\vbox{\hbox{\scriptstyle-}}\kern-3.25pt}}{{\vbox{\hbox{\scriptscriptstyle-}}\kern-2.29166pt}}{{\vbox{\hbox{\scriptscriptstyle-}}\kern-1.875pt}}\!\int_{-1}^{+1}du\,Q_{w_{1},u}Q_{w_{2},u}=0 for $-1<w_{1}<w_{2}<+1$ if we use the expressions $\alpha_{w}$ and $\beta_{w}$ for the isotropic case given in equation (3.8). In addition, the non-collective eigenfunctions are also orthogonal to the two collective ones. For a non-isotropic system, the eigenfunctions are linearly independent, but not orthogonal, so an explicit cross-check is less obvious.

3.3 Discrete angles

It is instructive to solve the Schrödiger-like equation (2.19) on a discrete grid of angles. For numerical solutions, this would have to be done anyway and we already had to take conceptual recourse to discrete angles to make sense of the continuum eigenfunctions and their normalisation. Therefore, we may determine the eigenvalues and eigenvectors of the “Hamilton operator” ${\cal H}$. In the discretised case it is an $N{\times}N$ matrix

[TABLE]

without summation of repeated indices on the rhs, where the angular spectrum is $G_{i}=G(u_{i})$. The grid points are

[TABLE]

so each represents the center of a bin of width $\Delta u$. In the noninteracting case the eigenvalues are $w_{i}=u_{i}$ and the normalised eigenvectors are $Q_{i,j}=\delta_{ij}$, which is the $j$th component of the $i$th eigenvector. In figure 3 we show the eigenvalues as a function of $\tilde{\mu}$ for the isotropic system, in analogy to the right panel of figure 1 where the horizontal lines were examples for the continuous distribution of eigenvalues, whereas here we show solutions of the discrete system with $N=6$ and 20 grid points. So here the eigenvalues $w_{i}$ are shifted relative to the undisturbed ones $u_{i}$. Even for a small number of bins, the two collective modes are very close to those found from the continuous solution. The two collective modes develop from the two modes at the edge of the light cone and peel off for increasing interaction strength.

While the shift of the eigenvalues may seem small for large $N$, this is not the case if we write them in the form

[TABLE]

In units of the splitting $\Delta u$ between eigenvalues, the shift is not small and approaches an asymptotic value for large $N$ for a fixed value of $i/N$.

By the same logic that allowed us to determine the form of the eigenfunctions in the continuum case, we find that the discrete eigenvector $k$ has components

[TABLE]

They are normalised as $\sum_{i=1}^{N}Q_{k,i}^{2}=1$. This sum is convergent because $\epsilon_{k}\not=0$, so the eigenvalues usually do not fall on exact grid points. They do fall on grid points in the noninteracting case and can do so also for special values of $\tilde{\mu}$, in which case the eigenvector is simply $Q_{k,i}=\delta_{ki}$. Starting from equation (3.17) we can make the transition to the continuous case with $\Delta u\to 0$, where the main trick is to split the term proportional to $A_{k}$ in the form

[TABLE]

In the continuous limit, the first term in a summation turns to the principal-part integration of $1/(w-u)$, whereas the second term becomes $\delta(w-u)$. An integral over $u$ turns to a summation over $i$, so in the second term we need to evaluate

[TABLE]

where we have extended the sum to $\pm\infty$ in the continuum limit $\Delta u\to 0$.

So more precisely, the second term in equation (3.18) turns to $\pi\,{\rm cotan}(\pi\epsilon_{w})\,\delta(w-u)$. Multiplying both terms with $\sin(\pi\epsilon_{w})/\pi$ reveals $\sin(\pi\epsilon_{w})/[\pi(w-u)]+\cos(\pi\epsilon_{w})\delta(w-u)$ for the relative contributions as in equation (3.12). We conclude that $\varphi_{w}=\pi\epsilon_{w}$, so the mixing angle between the two singular functions is equivalent, in the discrete case, to the offset $\epsilon_{i}$ between the grid points $u_{i}$ and the eigenvalues $w_{i}$ or equivalently to the shift of the eigenvalues relative to the ones of the noninteracting system.

4 Crossings of the angle distribution

4.1 Single crossing

For a general non-isotropic function $G_{u}$, the solutions are qualitatively similar to what we have shown so far if $G_{u}$ does not change sign on the interval $-1<u<+1$. However, in the presence of a crossing, tachyonic solutions can appear, which here mean complex eigenvalues of the Hamiltonian ${\cal H}$. We first consider a simple example of a schematic singly-crossed spectrum of the form $G_{u}={\textstyle\frac{1}{2}}-u$ that varies from $\frac{3}{2}$ to $-{\textstyle\frac{1}{2}}$ on the interval $-1\leq u\leq+1$. We show the discrete eigenvalues in figure 4 for $N=8$ and 32 for both real and imaginary part. We find both a real and a complex collective solution. As observed earlier, the collective solutions depend only mildly on the number of grid points. Fast flavor oscillations do not produce spurious instabilities that plague the numerical analysis of slow flavor oscillations.

To identify the points where the complex solutions branch off within the light cone, we note that for $\tilde{\mu}=0$ the only solutions are non-collective ones. The branch point appears for increasing $\tilde{\mu}$ and requires the merging of two real eigenvalues into a pair of complex conjugate ones. In the continuous case, every real $-1<w<+1$ is an eigenvalue with a unique eigenfunction determined by the value of $\alpha_{w}$ that follows from the determinant equation (3.6). The only exception is the location $w=u_{\rm c}$ where the function $G(u)$ has a crossing and thus $G_{w}=G_{u_{\rm c}}=0$ so that $\alpha_{u_{\rm c}}$ is not determined. Therefore, the branch point must occur for $w=u_{\rm c}$. In our example, this is for $w=u_{\rm c}=1/2$, shown as a horizontal dashed line in the upper panels of figure 4.

At the branch point we need two eigenfunctions for the two degenerate eigenvalues which merge at $w=u_{\rm c}$ for some value of $\tilde{\mu}$. This is only possible if the rhs of equation (3.6) vanishes, in which case this equation is fulfilled because both sides are separately zero. So the interaction strength for the branch points follows from

[TABLE]

This is a quadratic equation for $\tilde{\mu}$ with the solutions

[TABLE]

In our example we find $\tilde{\mu}_{\rm c}=-1/2$ and $3/2$, shown as vertical dashed lines.

This condition becomes more transparent if we write the angle spectrum, assumed to have a crossing at $u=u_{\rm c}$, in the form

[TABLE]

where $P(u)$ is some other function that could have crossings at other values of $u$. In analogy to equation (2.16) we introduce the moments of the $P(u)$ distribution

[TABLE]

This definition implies $F_{u_{\rm c}}=-{\cal P}_{0}$ and ${\cal G}_{n}={\cal P}_{n+1}-u_{\rm c}{\cal P}_{n}$. The interaction strength at the branching points is then given entirely by the $P$-moments as

[TABLE]

It is remarkable that $\tilde{\mu}_{\rm c}$ does not depend on the eigenvalue $w_{\rm c}$, i.e., $\tilde{\mu}_{\rm c}$ depends only on the first moments of the function $P(u)$. An alternative derivation will be presented in appendix A.

This result becomes yet clearer if we introduce averages \hbox{{\big{\langle}}\kern-3.5pt{\big{\langle}}}f(u)\hbox{{\big{\rangle}}\kern-3.6pt{\big{\rangle}}}=\int_{-1}^{+1}du\,P(u)\,f(u) that imply

[TABLE]

If the function $G(u)$ has exactly one crossing, then the function $P(u)$ has no crossing at all, implying that \hbox{{\big{\langle}}\kern-3.5pt{\big{\langle}}}(1+u)^{2}\hbox{{\big{\rangle}}\kern-3.6pt{\big{\rangle}}} and \hbox{{\big{\langle}}\kern-3.5pt{\big{\langle}}}(1-u)^{2}\hbox{{\big{\rangle}}\kern-3.6pt{\big{\rangle}}} has the same sign because the integrand $(1\pm u)^{2}>0$. So the discriminant is positive and $\tilde{\mu}_{\rm c}$ is real. We conclude that a single crossing indeed leads to branching points and thus to eigenmodes with complex $\Omega$.

The crossing at $u=u_{\rm c}$ could be of the form $G(u)=(u-u_{\rm c})^{n}Q(u)$ with positive integer $n$. If $n$ is odd, the previous argument does not change except that $P(u)$ would be non-negative instead of positive (or non-positive instead of negative) and thus would have no crossing. For even $n$ there would be no crossing of $G(u)$ because it would not change sign at $u_{\rm c}$. In other words, it is enough that $G(u)$ changes sign at $u_{\rm c}$, the change need not be a simple crossing of the form $u-u_{\rm c}$.

4.2 Multiple crossings

In principle, the function $G_{u}$ could change sign in several places. For $n$ crossings it can then be written in terms of some other function $P(u)$, assumed to have no crossings, as

[TABLE]

assuming also $u_{1}<u_{2}<\ldots<u_{n}$.

A simple example with two crossings is $G(u)=(u-u_{1})(u-u_{2})$ with $-1<u_{1}<u_{2}<+1$. Possible critical points are at $w_{\rm c}=u_{1,2}$. Considering $w_{\rm c}=u_{1}$, the discriminant in equation (4.6) is $\frac{16}{9}(4u_{2}^{2}-1)$ and thus positive for $|u_{2}|>\frac{1}{2}$. Everything is symmetric under $1\leftrightarrow 2$, so we have two instabilities if $|u_{1,2}|>\frac{1}{2}$, one instability if $|u_{j}|>\frac{1}{2}$ for $j=1$ or 2, but not both, and no instability at all if $|u_{1,2}|<\frac{1}{2}$. This simple example proves that a multi-crossed spectrum need not have any tachyonic solutions.

For a similar example with three crossings of the form $G(u)=(u-u_{1})(u-u_{2})(u-u_{3})$ we find cases with 1, 2 or 3 instabilities, depending on the values of $u_{1,2,3}$, but there is always at least one instability.

We suspect that an odd number of crossings would guarantee at least one unstable solution so that it would be enough for $G_{u}$ to have opposite signs at $u=\pm 1$. This would then be a more general definition of what we mean with a crossing of the spectrum. However, we have not pursued the question of how to prove this conjecture.

5 Collective motion vs. dissipation

5.1 General expressions

For the question of flavor conversion in an astrophysical environment with dense neutrinos, the existence of unstable solutions is the main question that can be addressed with a linearised stability analysis. If the initial conditions provide a seed for an exponentially growing solution, it will eventually dominate and push the system into the nonlinear regime. However, here we are interested in the role of non-collective modes for conditions without an exponentially growing solution and ask what will happen to a nonvanishing initial condition of flavor coherence. This might be a wave packet, but we here consider an even more elementary case of a plane wave with wave vector ${\bf K}$.

So effectively we consider solutions of the Hamiltonian equation (1.1) with a chosen interaction strength $\tilde{\mu}$. One way of defining the overall flavor coherence of the ensemble is by simply averaging over all modes

[TABLE]

The initial condition $S_{u}(0)$ can be expanded in the normalised eigenfunctions $Q_{w,u}$, each contributing an amplitude $T_{w}$ such that

[TABLE]

where $j=1,2$ refers to the two collective eigenmodes with $|w_{j}|>1$, while the non-collective modes are indexed with their eigenfrequency $-1\leq w\leq+1$. Therefore, the overall flavor coherence evolves as

[TABLE]

where

[TABLE]

for both collective and non-collective modes.

In the noninteracting case, there are no collective modes and the eigenfunctions of the non-collective ones are $Q_{w,u}=\delta(w-u)$, where we use $A=1$ for the global normalisation factor. Assuming $S_{u}(0)=1$ implies $T_{w}=1$ and $\overline{Q}_{w}=G_{w}$ so that

[TABLE]

As expected, the flavor coherence dissipates as $t^{-1}$ by kinematical decoherence (de-phasing) of the continuum of energy eigenstates $w$.

5.2 Isotropic system

In the interacting case, the Hamiltonian of equation (2.20) is usually not Hermitean and the eigenfunctions not orthogonal, so determining the amplitudes $T_{w}$ requires the dual basis. The one exception is the isotropic system ($G_{u}=1$), where the Hamiltonian is symmetric and thus the eigenfunctions are orthogonal as discussed in section 3.2, so for simplicity we use this case to illustrate the transition from decohering to collective behavior. In addition we assume $S_{u}(0)=1$, implying for the non-collective modes

[TABLE]

and similar for the collective ones, so

[TABLE]

For the non-collective modes, we use the eigenfunctions of the first line of equation (3.12) with the coefficients $\alpha_{w}$ and $\beta_{w}$ for the isotropic case given in equation (3.8), providing

[TABLE]

For vanishing $\tilde{\mu}$ this is $T_{w}^{2}=1$ and we recover equation (5.5). For large $\tilde{\mu}$, one finds

[TABLE]

so in this limit the contribution of the non-collective modes decreases with $\tilde{\mu}^{-4}$ and the evolution becomes entirely collective.

The contributions $T_{1,2}^{2}$ of the collective modes are straightforward to calculate from the collective eigenfunctions and eigenvalues, but the expressions become cumbersome. In figure 5 we show ${\textstyle\frac{1}{2}}T_{1,2}^{2}$ as a function of $\tilde{\mu}$ as well as the total contribution from the non-collective modes. The three contributions add up to unity. In our example of an isotropic medium, the collective mode No. 1 has the eigenfunction $Q_{1,u}=1/\sqrt{2}$ in the large-$\tilde{\mu}$ limit, identical to the initial condition of $S_{u}$, leaving no projection on the other modes, so that for large $\tilde{\mu}$ only one eigenfrequency contributes. (The other collective mode has the large-$\tilde{\mu}$ eigenfunction $Q_{2,u}=\sqrt{3/2}\,u$ and thus has no overlap with a constant function on the interval $-1<u<+1$.)

We conclude that for $\tilde{\mu}\lesssim 0.3$, the initial flavor coherence dissipates by decoherence, whereas for somewhat larger interaction strength one of the collective modes takes over. For an initial wavepacket, the components with $K$ exceeding a few times $\mu$ completely dissipate, whereas the smaller-$K$ components propagate as collective flavor excitations. In an intermediate range, the initial flavor coherence partially dissipates.

We mention in passing the similarity to the question of synchronized oscillations vs. decoherence at intermediate densities for ordinary collective oscillatoins (driven by the mass term) that was addressed in reference [43]. In that case, the dephasing effect came from different energies in $\Delta m^{2}/2E$, whereas the angle degree of freedom was intergrated out.

6 Conclusions

We have re-examined the dispersion relation for the mean field of flavor coherence in a dense neutrino gas, considering in particular “fast modes,” where neutrino masses do not enter and the question of which modes are supported depends on the neutrino angle distribution, or rather, on the angle distribution of the electron lepton number (ELN) carried by neutrinos. The purpose of our discussion was primarily to illuminate the connection between non-collective modes that exist even in the absence of neutrino-neutrino interactions and collective modes that appear as a result of neutrino refraction in a background of dense neutrinos.

In the absence of collective effects, the angle-dependent neutrino flux simply carries along any putative flavor coherence that may have been imprinted initially. Unless we consider very special angle distributions, such as two-beam toy models that have sometimes been used in the literature, the angular range of neutrino directions will quickly decohere any initially prepared flavor coherence.

The presence of neutrino-neutrino interactions causes the appearance of collective modes that may or may not show tachyonic instabilities, whereas the non-collective modes get modified. We have derived explicit expressions for the non-collective eigenfunctions and we have also shown that, as the interaction strength is increased, the amplitude of flavor coherence carried by these modes decreases. For a large interaction strength, the emergent collective modes dominate the evolution.

For fast flavor modes, the neutrino-neutrino interaction energy that we have defined as $\mu=\sqrt{2}G_{\rm F}(n_{\nu_{e}}+n_{\bar{\nu}_{e}})$ is the only energy scale of the problem, so all wave vectors and frequencies are naturally measured in units of $\mu$. So for flavor modes with $K\gg\mu$ collective effects are practically irrelevant, whereas for $K\lesssim\mu$ they are important. Of course, for unstable modes, collective effects are always important in the sense that any unstable mode eventually dominates if it has enough time to grow.

We have also clarified that the appearance of tachyonic collective modes can be understood as the coalescence of two non-collective modes, i.e., two real eigenfrequencies merge at a critical point to become two complex-conjugate solutions. This happens at crossings of the neutrino angle distribution and allows us to calculate the required interaction strength. As a corollary we show that a single crossing of the ELN angle distribution guarantees a tachyonic instability, whereas the case of several crossings is less clear at present.

In this way the role of the non-collective and collective modes of the linearised equations of motion of neutrino flavor coherence has been further illuminated. Of course, these dispersion relations only allow us to judge which kind of solutions are supported by a given dense neutrino gas, but do not show which ones will actually be excited and by what. Flavor coherence of the neutrino mean field probably can be sourced ultimately only by the mass term, perhaps in conjunction with density fluctuations of the medium. In this sense, understanding the dispersion relation is probably better suited to understand if tachyonic solutions are or are not supported by the neutrino gas, but if they are actually excited, and with which strength, and how they evolve in the nonlinear regime are questions that go beyond the topics that can be addressed by the dispersion relation alone.

Acknowledgments

We acknowledge partial support by the Deutsche Forschungsgemeinschaft through Grants No. SFB 1258 (Collaborative Research Center Neutrinos, Dark Matter, Messengers) and EXC 2094 (Excellence Cluster Origins), as well as the European Union through Grant No. H2020-MSCA-ITN-2015/674896 (Innovative Training Network Elusives).

Appendix A Critical points: Alternative derivation

We can derive the critical points $(w_{\rm c},\tilde{\mu}_{\rm c})$ derived in section 4.1 where a tachynoic solution appears without having to worry about the non-collective eigenfunctions. To this end we recall that any collective mode fulfills equation (2.13) which reads, with $\tilde{\mu}=\mu/K$,

[TABLE]

Considering a complex eigenvalue, we now interpret $w$ as its real part and $\kappa$ its imaginary part so that the definition of equation (2.14) now reads explicitly

[TABLE]

At the critical point, the imaginary part is infinitesimally small, so we can use the Sokhotski-Plemelj Theorem in the form

[TABLE]

where the integral on the rhs is the Cauchy principal value. Therefore, in equation (A.1) we should insert

[TABLE]

Equation (A.1) then must be fulfilled simultaneously for the real and imaginary part of the determinant, the real part implying

[TABLE]

The imaginary part vanishes if

[TABLE]

One solution is $G_{w_{\rm c}}=0$ which happens exactly at the crossings of $G_{u}$, i.e., when $w_{\rm c}=u_{\rm c}$. If we express again $G(u)=(u-u_{\rm c})P(u)$, then the quantities $p_{n}$, for $w_{\rm c}=u_{\rm c}$, correspond to ${\cal P}_{n}$ defined in equation (4.4), so our new equation (A.5) is identical with the previous result of equation (4.5).

Another solution of equation (A.6) is when the expression in square brackets vanishes, providing another expression for $\tilde{\mu}_{\rm c}$ that must agree with equation (A.5), in principle allowing us to calculate an expression for $w_{\rm c}$. Simple examples do not provide consistent solutions of these equations. The discussion in section 4.1 suggests that critical points only appear at crossings of $G(u)$ where $G(w_{\rm c})=0$, but no proof is available that there could not exist unforeseen special cases of other situations.