# Do Neutrinos Become Flavor Unstable Due to Collisions with Matter in the
Supernova Decoupling Region?

Shashank Shalgar ID Niels Bohr International Academy & DARK, Niels Bohr
Institute,
University of Copenhagen, Blegdamsvej 17, 2100 Copenhagen, Denmark  Irene
Tamborra ID Niels Bohr International Academy & DARK, Niels Bohr Institute,
University of Copenhagen, Blegdamsvej 17, 2100 Copenhagen, Denmark

###### Abstract

In core-collapse supernovae, the neutrino density is so large that neutrino
flavor instabilities, leading to flavor conversion, can be triggered by the
forward scattering of neutrinos among each other, if a crossing between the
angular distributions of electron neutrinos and antineutrinos exists (fast
instability in the limit of vanishing vacuum frequency) or in the presence of
perturbations induced by the neutrino vacuum frequency (slow instability).
Recently, it has been advanced the conjecture that neutrino collisions with
the medium could be another mean to kickstart flavor change (collisional
instability). Inspired by a spherically symmetric core-collapse supernova
model with mass $18.6\ M_{\odot}$, we compute the neutrino angular
distributions solving the kinetic equations for an average energy mode and
investigate the occurrence of flavor instabilities at different post-bounce
times, ranging from the accretion phase to the early cooling phase. We find
that fast and slow flavor instabilities largely dominate over the collisional
ones in the decoupling region for all post-bounce times. While more work is
needed to assess the relevance of collisional instabilities in neutrino-dense
environments, our findings suggest that neutrino collisions with matter affect
the flavor evolution in the decoupling region, but are not responsible for
triggering flavor conversion, if crossings in the neutrino lepton number
angular distribution exist.

##  I Introduction

In dense neutrino astrophysical environments, such as the interior of core-
collapse supernovae (SNe) and neutron star mergers, neutrino flavor evolution
has a rich phenomenology owing to the nonlinear nature of its evolution
Mirizzi  _et al._ (2016); Duan  _et al._ (2010); Chakraborty  _et al._
(2016a); Tamborra and Shalgar (2021); Richers and Sen (2022). The nonlinearity
is a consequence of the coherent forward scattering of neutrinos among each
other. Flavor conversion is classified as fast, if purely driven by the
neutrino-neutrino interaction term for vanishing neutrino vacuum frequency
Sawyer (2005, 2009, 2016); Chakraborty  _et al._ (2016b); Izaguirre  _et al._
(2017); Tamborra and Shalgar (2021); in this case, a necessary but not
sufficient condition for the development of fast flavor instabilities is the
existence of electron lepton number (ELN) crossings in the angular
distribution of electron neutrinos and antineutrinos (note that, under the
assumption of periodic boundary conditions, the existence of an ELN crossing
becomes a sufficient condition for the development of a flavor instability),
albeit a large growth rate for such flavor instabilities does not necessarily
imply large flavor conversion Padilla-Gay  _et al._ (2022a); Morinaga (2022).
In addition, neutrino self-interactions could be triggered at high densities
by the vacuum frequency in the Hamiltonian Shalgar and Tamborra (2023a, b); in
this case, we consider the flavor instability to be slow since it develops
thanks to the vacuum frequency, which is smaller than the neutrino-neutrino
interaction strength; note that, in this case, a system with ELN crossings
which is stable for vanishing vacuum frequency may become unstable.

Owing to the nonlinearity of the system, the momentum changing collisional
processes, which are lower order processes, can have a significant impact on
the flavor evolution, despite the widely different characteristic time scales
characterizing flavor evolution and collisions Shalgar and Tamborra (2021);
Johns (2023); Sasaki and Takiwaki (2022); Hansen  _et al._ (2022); Kato and
Nagakura (2022); Johns and Xiong (2022); Padilla-Gay  _et al._ (2022b); Lin
and Duan (2023); Liu  _et al._ (2023); Kato  _et al._ (2023). The collisional
processes of interest include momentum changing scattering as well as
reactions responsible for the creation and absorption of neutrinos Burrows
_et al._ (2006); Janka (2012); Mezzacappa  _et al._ (2020); Shalgar and
Tamborra (2019). As the neutrino emission properties are affected by the
collisional processes, before and during flavor evolution, it is imperative to
understand their interplay. First attempts coupling collisions and neutrino
flavor evolution self-consistently have been reported in the literature, but
there is still a long journey ahead Shalgar and Tamborra (2023a, b); Nagakura
(2023); Xiong  _et al._ (2023a).

Existing work employing a given initial ELN angular distribution with a
crossing finds that the effect of collisions on flavor evolution is to
isotropize the angular distribution, with a resulting enhancement or
suppression of flavor conversion according to the initial ELN distributions
Shalgar and Tamborra (2021); Sasaki and Takiwaki (2022); Hansen  _et al._
(2022); Martin  _et al._ (2021). On the other hand, a system with no flavor
instability in the absence of the collision term can also become unstable in
the flavor space because of the so-called collisional instabilities Johns
(2023); Johns and Xiong (2022); Martin  _et al._ (2021); Lin and Duan (2023);
Padilla-Gay  _et al._ (2022b); Xiong  _et al._ (2023a, b).

In this paper, we aim to investigate slow 111In this work, slow flavor
instability indicates a flavor instability that vanishes for zero vacuum
frequency, independently on the existence of crossings in the lepton number
angular distribution. For fast instabilities, instead, we mean instabilities
that arise in the limit of vanishing vacuum frequency., fast, and collisional
flavor instabilities, not relying on arbitrary neutrino angular distributions,
rather solving the neutrino kinetic equations with the collision terms for an
average energy mode to obtain the neutrino angular distributions for each
flavor. To this purpose, we rely on static hydrodynamic backgrounds and
thermodynamic properties extrapolated at different post-bounce snapshots of a
spherically symmetric core-collapse SN model with mass $18.6\ M_{\odot}$
Archive . To compute the neutrino angular distributions, in the collisional
kernel, we include emission and absorption of neutrinos through
Bremsstrahlung, pair production, and beta reactions. In addition, we consider
the direction changing effects due to neutral current scattering of neutrinos
on nucleons. We then investigate through the linear stability analysis whether
collisional, slow, or fast flavor instabilities exist.

The manuscript is organized as follows. Section II introduces the neutrino
quantum kinetic equations, then we discuss the procedure adopted to obtain the
neutrino angular distributions solving the kinetic equations in the absence of
flavor conversion in Sec. III. In Sec. IV, we perform the linear stability
analysis for the homogeneous mode, with as well as without the collision term
in order to understand the effect of collisional instabilities, while we
discuss the impact of spatial inhomogeneities in Sec. V. We summarize our
findings and offer an outlook in Sec. VI. In addition, Appendix A provides an
outline of the collision terms responsible for the production, absorption, and
momentum changing interactions; while a comparison of the neutrino number
densities computed in this work with the ones obtained from the hydrodynamical
SN simulation is shown in Appendix B.

##  II Neutrino quantum kinetic equations

Figure 1: Schematic diagram of the shell within which we solve the neutrino
kinetic equations taking into account emission, absorption, and direction
changing processes. The polar angle $\theta$ is defined with respect to the
radial direction. It should be noted that $\theta$ changes with the radius for
a given trajectory.

For the sake of simplicity, we consider two monoenergetic neutrino flavors,
$\nu_{e}$ and $\nu_{x}$ with $x=\mu,\tau$ (or a linear combination of the
two). The equations of motion of neutrinos and antineutrinos are

$\displaystyle i\left(\frac{\partial}{\partial
t}+\vec{v}\cdot\vec{\nabla}\right)\rho$ $\displaystyle=$
$\displaystyle[H,\rho]+i\mathcal{C}[\rho,\bar{\rho}]\ ,$ (1) $\displaystyle
i\left(\frac{\partial}{\partial t}+\vec{v}\cdot\vec{\nabla}\right)\bar{\rho}$
$\displaystyle=$
$\displaystyle[\bar{H},\bar{\rho}]+i\bar{\mathcal{C}}[\rho,\bar{\rho}]\ ,$ (2)

with $\rho$ and $\bar{\rho}$ being $2\times 2$ density matrices representing
the flavor content of neutrinos and antineutrinos, respectively. The diagonal
components of the density matrices ($\rho_{ii}$, with $i=e,x$) represent the
neutrino occupation numbers [in units of $1/\textrm{cm}^{3}$]; the number
density of a given species, $n_{i}$, at a given location is obtained by
$\int_{-1}^{1}\rho_{ii}d\cos\theta$. In Eqs. 1 and 2, the dependence of the
density matrices on the polar angle $\theta$ (measured with respect to the
radial direction, see Fig. 1) and radius $r$ is suppressed for the sake of
brevity.

The velocity vector is represented by $\vec{v}$ while $\vec{\nabla}$ is the
gradient operator. The inner product of the two is given by:

$\displaystyle\vec{v}\cdot\vec{\nabla}=\cos\theta\frac{\partial}{\partial
r}+\left(\frac{\sin^{2}\theta}{r}\right)\frac{\partial}{\partial\cos\theta}\
,$ (3)

where the dependence on the azimuthal angle $\phi$ is ignored due to the
assumption of spherical symmetry.

The Hamiltonians for neutrinos and antineutrinos
($H=H_{\textrm{vac}}+H_{\textrm{mat}}+H_{\nu\nu}$ and
$\bar{H}=-H_{\textrm{vac}}+H_{\textrm{mat}}+H_{\nu\nu}$, respectively) are
made of the vacuum, matter, and the self-interaction terms:

$\displaystyle H_{\textrm{vac}}$ $\displaystyle=$
$\displaystyle\frac{\omega_{\rm{vac}}}{2}\begin{pmatrix}-\cos
2\theta_{\textrm{V}}&\sin 2\theta_{\textrm{V}}\\\ \sin
2\theta_{\textrm{V}}&\cos 2\theta_{\textrm{V}}\end{pmatrix}\ ,$ (4)
$\displaystyle H_{\textrm{mat}}$ $\displaystyle=$
$\displaystyle\begin{pmatrix}\lambda&0\\\ 0&0\end{pmatrix}\ ,$ (5)
$\displaystyle H_{\nu\nu}$ $\displaystyle=$
$\displaystyle\sqrt{2}G_{\textrm{F}}\int
d\cos\theta^{\prime}\left[\rho(\cos\theta^{\prime})-\bar{\rho}(\cos\theta^{\prime})\right]$
(6) $\displaystyle\times$ $\displaystyle(1-\cos\theta\cos\theta^{\prime})\ .$

In the vacuum Hamiltonian, the vacuum frequency is $\omega_{\rm{vac}}=\Delta
m^{2}/2\epsilon_{\nu}$, with $\Delta m^{2}=2.5\times 10^{-3}$, and
$\epsilon_{\nu}$ being the average neutrino energy calculated over all
species, while the mixing angle is $\theta_{\textrm{V}}=10^{-3}$. In the
matter term, the strength of neutrino-matter interactions is given by
$\lambda=\sqrt{2}G_{\textrm{F}}n_{e}$, with the electron number density being
$n_{e}=Y_{e}\rho_{B}/m_{N}$ where $Y_{e}$ is the electron fraction, $\rho_{B}$
is the baryon mass density, and $m_{N}$ the nucleon mass (in the following, we
restrict our analysis to the linear regime in which the matter term can be
ignored because it does not affect the instability of the system). Note that,
when we carry out the linear stability analysis in the following sections, we
neglect the matter term and effectively account for it though a smaller value
for the mixing angle (i.e., $\theta_{\textrm{V}}=10^{-3}$). The strength of
the self-interaction potential, generally denoted by
$\mu=\sqrt{2}G_{\textrm{F}}n_{\nu_{e}}$, characterizes the $\nu$–$\nu$ term of
the Hamiltonian.

The collision term, schematically represented by the functionals $\mathcal{C}$
and $\bar{\mathcal{C}}$, includes the absorption, emission, and direction
changing processes defined as outlined in Appendix A:

$\displaystyle\mathcal{C}[\rho,\bar{\rho}]$ $\displaystyle=$
$\displaystyle\mathcal{C}_{\textrm{emit}}-\mathcal{C}_{\textrm{absorb}}\odot\rho(\cos\theta)$
(7) $\displaystyle+$ $\displaystyle\frac{\mathcal{C}_{\textrm{dir-ch}}}{2}\int
d\cos\theta^{\prime}\left[-\rho(\cos\theta)+\rho(\cos\theta^{\prime})\right]$
$\displaystyle+$ $\displaystyle\cos\theta\ \mathcal{C}_{\textrm{ani}}\int
d\cos\theta^{\prime}\cos\theta^{\prime}\rho(\cos\theta^{\prime})\ ,$
$\displaystyle\bar{\mathcal{C}}[\rho,\bar{\rho}]$ $\displaystyle=$
$\displaystyle\bar{\mathcal{C}}_{\textrm{emit}}-\bar{\mathcal{C}}_{\textrm{absorb}}\odot\bar{\rho}(\cos\theta)$
(8) $\displaystyle+$ $\displaystyle\frac{\bar{\mathcal{C}}_{\textrm{dir-
ch}}}{2}\int
d\cos\theta^{\prime}\left[-\bar{\rho}(\cos\theta)+\bar{\rho}(\cos\theta^{\prime})\right]$
$\displaystyle+$ $\displaystyle\cos\theta\
\bar{\mathcal{C}}_{\textrm{ani}}\int
d\cos\theta^{\prime}\cos\theta^{\prime}\bar{\rho}(\cos\theta^{\prime})\ .$

Here,
$\mathcal{C}_{\textrm{emit}}=\textrm{diag}(\mathcal{C}_{\textrm{emit}}^{ee},\mathcal{C}_{\textrm{emit}}^{xx})$
and
$\bar{\mathcal{C}}_{\textrm{emit}}=\textrm{diag}(\bar{\mathcal{C}}_{\textrm{emit}}^{ee},\bar{\mathcal{C}}_{\textrm{emit}}^{xx})$
govern the emission rate for each neutrino flavor.
$\mathcal{C}_{\textrm{absorb}}$ and $\bar{\mathcal{C}}_{\textrm{absorb}}$ are
matrices and $\odot$ represents element-wise multiplication:

$\displaystyle\mathcal{C}_{\textrm{absorb}}\odot\rho(\cos\theta)$
$\displaystyle=$
$\displaystyle\mathcal{C}_{\textrm{absorb}}^{ee}\rho_{ee}+\mathcal{C}_{\textrm{absorb}}^{ex}\rho_{ex}$
$\displaystyle+$
$\displaystyle\mathcal{C}_{\textrm{absorb}}^{xe}\rho_{xe}+\mathcal{C}_{\textrm{absorb}}^{xx}\rho_{xx}$
$\displaystyle\bar{\mathcal{C}}_{\textrm{absorb}}\odot\bar{\rho}(\cos\theta)$
$\displaystyle=$
$\displaystyle\bar{\mathcal{C}}_{\textrm{absorb}}^{ee}\bar{\rho}_{ee}+\bar{\mathcal{C}}_{\textrm{absorb}}^{ex}\bar{\rho}_{ex}$
$\displaystyle+$
$\displaystyle\bar{\mathcal{C}}_{\textrm{absorb}}^{xe}\bar{\rho}_{xe}+\bar{\mathcal{C}}_{\textrm{absorb}}^{xx}\bar{\rho}_{xx}\
,$

which is linearly proportional to the number of neutrinos at a given location
representing the rate at which neutrinos are absorbed by the surrounding
matter Richers  _et al._ (2019). The off-diagonal components of
$\mathcal{C}_{\textrm{absorb}}$ and $\bar{\mathcal{C}}_{\textrm{absorb}}$ are
the average of the diagonal elements of the diagonal terms. The terms
$\mathcal{C}_{\textrm{dir-ch}}$ and $\bar{\mathcal{C}}_{\textrm{dir-ch}}$
describe the direction changing neutral current scatterings, which do not
change the number of neutrinos and are flavor independent. The
$\mathcal{C}_{\rm{ani}}$ and $\bar{\mathcal{C}}_{\rm{ani}}$ terms depend on
the proton and neutron fractions. Moreover, we take into account the
anisotropy in the direction changing scatter term.

In Eqs. 7 and 8, the emission and absorption terms seem to be independent of
each other, however, that is not the case and they are related by the
Kirchoff’s law as explained in Appendix A. Hence, given either the absorption
or the emission rate, the other can be calculated (for this reason, in
Appendix A, we only provide the absorption rates for all the processes).

##  III Steady state neutrino angular distributions

Figure 2: Radial profiles of the characteristic quantities extracted from our
benchmark hydrodynamical SN simulation (see main text for details) at the
post-bounce times $t_{\rm{pb}}=0.05$, $0.12$, $0.25$, $0.5$, $0.75$, and $1$ s
and then adopted to compute the neutrino angular distributions. The different
panels represent the baryon number density, the electron temperature, the
neutrino and electron chemical potentials, and the neutron and electron
fractions, from top left to bottom right, respectively.

We adopt the outputs of a one-dimensional hydrodynamical simulation of a SN
with mass $18.6\ M_{\odot}$, SFHo nuclear equation of state without muons, and
gravitational mass $1.4\ M_{\odot}$ Archive . We consider static
hydrodynamical backgrounds and thermodynamical quantities, as well as the
local neutrino number density and first energy moment, for selected time
snapshots extracted at post-bounce times $t_{\rm{pb}}=0.05$, $0.12$, $0.25$,
$0.5$, $0.75$, and $1$ s, which are representative of the accretion as well as
the early Kelvin-Helmoltz cooling phases. Then, we use these quantities as
inputs to calculate the collision terms (Eqs. 7 and 8). The relevant
properties of our benchmark hydrodynamical SN simulations adopted to compute
the collision terms are illustrated in Fig. 2 (see also Appendix A).

The neutrino angular distributions are not provided by the hydrodynamic
simulation, hence, we compute them solving Eqs. 1 and 2 and imposing
$H=\bar{H}=0$. Since we are only interested in the angular distributions in
the proximity of the decoupling region, we solve Eqs. 1 and 2 in a radial
shell (see Fig. 1) whose width $[r_{\rm{min}},r_{\rm{max}}]$ differs for each
post-bounce time as summarized in Table 1 and has been chosen such that the
minimum radius $r_{\rm min}$ has a number density of neutrinos that coincides
with the energy integrated Fermi-Dirac distribution extracted from the
hydrodynamical simulation (i.e., at large optical depth, such that neutrinos
are fully coupled with matter), while neutrinos are completely decoupled at
the maximum radius $r_{\rm max}$ (i.e., the backward component of the neutrino
and antineutrino number densities is negligible: $\rho_{ii}(\cos\theta)\approx
0$ for $\cos\theta<0$). We solve the neutrino kinetic equations starting with
no neutrinos and antineutrinos in the shell. As the system evolves in time,
the simulation shell is populated with neutrinos and antineutrinos with the
boundary conditions being guided by the hydrodynamical simulation as explained
above.

Table 1: Values of $r_{\textrm{min}}$ and $r_{\textrm{min}}$ used for each post-bounce time snapshot to take into account the varying location and width of the decoupling region (see also Fig. 1). $t_{\rm{pb}}$ (s) |  $r_{\textrm{min}}$ (km) |  $r_{\textrm{min}}$ (km)
---|---|---
$0.05$ | 25 | 150
$0.12$ | 25 | 100
$0.25$ | 22 | 57
$0.5$ | 20 | 35
$0.75$ | 17 | 32
$1$ | 16 | 31

For each time snapshot, we use a grid with $75$ radial bins and $75$ angular
bins uniformly spread over $\cos\theta$ and solve the neutrino kinetic
equations for the time interval required to achieve a classical steady state
configuration.

Figure 3 shows the radial profiles of the angle-integrated neutrino number
densities computed for all three flavors through the procedure described
above. While involving some approximations, our method allows us to reproduce
a local neutrino density in agreement with the output of our benchmark
hydrodynamical simulation, as shown in Appendix B. We note that, due to the
single energy approximation used in this paper for simplicity, we slightly
overestimate the number density of $\bar{\nu}_{e}$ with respect to the one of
$\nu_{e}$ at large radii when comparing with our benchmark hydrodynamical
simulation (see also Fig. 11). As shown later, the fast flavor instability is
triggered by the presence of ELN crossings that would be still present,
despite the single-energy approximation; hence, the relative importance of
fast flavor instabilities with respect to the collisional ones is not affected
by the single-energy treatment. We stress that the goal of the procedure
illustrated in this section is not to reproduce the neutrino properties
obtained through the hydrodynamical simulation, rather to obtain neutrino
angular distributions that are not arbitrary but inspired by conditions found
in the SN interiors.

Figure 3: Radial evolution of the number densities of $\nu_{e}$ (in red),
$\bar{\nu}_{e}$ (in green), and $\nu_{x}$ (in blue) for selected post-bounce
time snapshots and computed solving the neutrino kinetic equations. We mark
the radius of neutrino decoupling through a vertical line; this is the radius
at which the number density deviates from a Fermi-Dirac distribution (see
Appendix B). One can see that the decoupling surface shrinks as a function of
time. Note that the radial range is not the same for all time snapshots.

We define the decoupling radius as the one at which the neutrino number
density differs from the Fermi-Dirac one up to $1\%$ (see Appendix B) and mark
such radius in Fig. 3 to guide the eye. One can see from Fig. 3 that the
decoupling radius becomes smaller as the post-bounce time increases. Because
of the different nature of their interactions, while the electron flavors
decouple over a broader radial range, the non-electron flavor decouples almost
instantaneously Tamborra  _et al._ (2017).

Figure 4 displays an example of the radial evolution of the neutrino angular
distributions computed through the procedure illustrated above for $t_{\rm
pb}=0.05$ s. As also found in Ref. Tamborra  _et al._ (2017), the angular
distributions become progressively forward peaked as a function of the radial
distance and have an angular spread that is the largest for $\nu_{e}$,
followed by $\bar{\nu}_{e}$ and $\nu_{x}$.

Figure 4: Polar diagram of the radial variation of the spectral intensity of
$\nu_{e}$, $\bar{\nu}_{e}$, and $\nu_{x}$, from top to bottom respectively,
computed solving the neutrino kinetic equations. The selected radii and the
time snapshot $t_{\rm pb}=0.05$ s have been chosen for illustrative purposes.
One can see that the neutrino distributions become progressively forward
peaked as $r$ increases, with an angular spread that is largest for $\nu_{e}$,
followed by $\bar{\nu}_{e}$, and then $\nu_{x}$.

Figure 5 focuses on the radial evolution of the ELN angular distribution for a
representative post-bounce time snapshot ($t_{\rm pb}=0.05$ s), interestingly,
we find ELN crossings for most radii. This finding holds for all post-bounce
times considered in this work for our benchmark SN model and is in contrast
with the results presented in Ref. Tamborra  _et al._ (2017); the latter
concluded, through a systematic analysis of different one-dimensional SN
models, that ELN crossings are not present in one-dimensional hydrodynamical
simulations, albeit ELN crossings were found in multi-dimensional SN
simulations Morinaga  _et al._ (2020); Abbar  _et al._ (2020, 2019); Delfan
Azari  _et al._ (2020); Nagakura  _et al._ (2019); Abbar  _et al._ (2021). The
results presented in Fig. 5 seem to fulfill the criterion of Ref. Shalgar and
Tamborra (2019) that linked the occurrence of ELN crossings to comparable
number densities of neutrinos and antineutrino (see also Refs. Lin and Duan
(2023); Kato  _et al._ (2023)). Yet, the reason why ELN crossings are found in
this work, but not in Ref. Tamborra  _et al._ (2017) requires further
investigation and will be the focus of an upcoming dedicated paper et al. —we
stress that the purpose of this work is to rely on the angular distributions
derived from first principles to investigate the interplay among fast, slow,
and collisional instabilities. The different findings could be connected to
the fact that the proto-neutron star cooling model of our benchmark SN
simulation includes effects from proto-neutron star convection (treated by a
mixing-length approximation Mirizzi  _et al._ (2016); Weiss  _et al._ (2004))
in contrast to the SN models investigated in Ref. Tamborra  _et al._ (2017).
Other reasons might be the different nuclear equation of state (SFHo for our
benchmark SN model vs. Lattimer & Swesty with compressibility modulus $K=220$
MeV in Ref. Tamborra  _et al._ (2017)) or secondary subtle differences in our
neutrino transport results compared to those obtained with the VERTEX
transport used in Ref. Tamborra  _et al._ (2017).

Figure 5: Angular distribution of the ELN for $t_{\rm pb}=0.05$ s extracted at
different radii as indicated. We find ELN crossings for all considered post-
bounce time snapshots.

##  IV Linear stability analysis for the homogeneous mode

The radial range where significant flavor evolution might develop can be
determined through the linear stability analysis Banerjee  _et al._ (2011);
Airen  _et al._ (2018); Izaguirre  _et al._ (2017); Morinaga and Yamada
(2018). In this section, we linearize Eqs. 1 and 2 assuming homogeneity (i.e.,
ignoring the advective term in the neutrino kinetic equations) and perform the
linear stability analysis with and without the collision term to investigate
the impact of collisional instabilities. Since we consider monoenergetic
neutrinos for simplicity, we here assume that the neutrino energy
$\epsilon_{\nu}$ entering $\omega_{\rm{vac}}$ changes as a function of the
radius and it coincides with the average of the first energy moments extracted
from the hydrodynamical SN simulation over all four flavors (see bottom right
panel of Fig. 9).

The linearized version of Eqs. 1 and 2 for the off-diagonal terms of the
density matrices at a given $r$ and $\cos\theta$ is

$\displaystyle i\frac{\partial\rho_{ex}}{\partial t}$ $\displaystyle=$
$\displaystyle\left[(H_{ee}-H_{xx})\rho_{ex}-(\rho_{ee}-\rho_{xx})H_{ex}\right]$
(10) $\displaystyle-$
$\displaystyle\frac{i}{2}(\mathcal{C}^{e}_{\textrm{absorb}}+\mathcal{C}^{x}_{\textrm{absorb}})\rho_{ex}$
$\displaystyle-$
$\displaystyle\frac{i}{4}(\mathcal{C}^{e}_{\textrm{dirch}}+\mathcal{C}^{x}_{\textrm{dir-
ch}})$ $\displaystyle\times$
$\displaystyle\int[\rho_{ex}(\cos\theta)-\rho_{ex}(\cos\theta^{\prime})]d\cos\theta^{\prime}\
,$ $\displaystyle i\frac{\partial\bar{\rho}_{ex}}{\partial t}$
$\displaystyle=$
$\displaystyle\left[(\bar{H}_{ee}-\bar{H}_{xx})\bar{\rho}_{ex}-(\bar{\rho}_{ee}-\bar{\rho}_{xx})\bar{H}_{ex}\right]$
(11) $\displaystyle-$
$\displaystyle\frac{i}{2}(\bar{\mathcal{C}}^{e}_{\textrm{absorb}}+\bar{\mathcal{C}}^{x}_{\textrm{absorb}})\bar{\rho}_{ex}$
$\displaystyle-$
$\displaystyle\frac{i}{4}(\bar{\mathcal{C}}^{e}_{\textrm{dirch}}+\bar{\mathcal{C}}^{x}_{\textrm{dir-
ch}})$ $\displaystyle\times$
$\displaystyle\int[\bar{\rho}_{ex}(\cos\theta)-\bar{\rho}_{ex}(\cos\theta^{\prime})]d\cos\theta^{\prime}\
,$

It should be noted that, due to the direction changing contribution in the
collision term, the linearized equations for $\rho_{ex}(r,\cos\theta)$ and
$\bar{\rho}_{ex}(r,\cos\theta)$ depend on the density matrix at all other
angles, due to the integral in the last term of Eqs. 10 and 11.

Equations 10 and 11 form an eigensystem for which we calculate the eigenvalues
following the procedure described in Ref. Banerjee  _et al._ (2011). In the
absence of the collision term, the eigenvalues of the system are complex
conjugate pairs; however, due to the presence of the collision term this is no
longer true for Eqs. 10 and 11. The eigenmode that has the eigenvalue with the
largest positive imaginary term, which we denote by $\kappa$, dominates and is
responsible for an instability for which we can compute its growth rate.

The growth rate of fast, slow, and collisional flavor instabilities is linked
to three time scales of our system, $\mu^{-1}$, $\omega_{\textrm{vac}}^{-1}$,
and the time between two successive collisions of neutrinos–see also Refs.
Shalgar and Tamborra (2023b, a). The time scales associated with the growth
rate is determined by the fastest time scale that contributes to the flavor
instability.

Figure 6: Radial profile of the growth rate of the flavor instability,
obtained under the assumption of homogeneity, at various time snapshots (from
$t_{\rm pb}=0.05$ s in the top left panel to $t_{\rm pb}=1$ s in the bottom
right panel). Each panel displays results obtained for four different
scenarios. The blue line uses the vacuum and the self-interaction Hamiltonian
(slow instability), while the green line only considers the self-interaction
Hamiltonian (fast instability); for both blue and green lines, collisional
damping is ignored. The same calculation is repeated to obtain the dashed red
and dot-dashed orange lines, but with collisional damping. It should be noted
that, except for $t_{\rm pb}=0.12$ s, the four lines are almost on top of each
other and difficult to distinguish. The location of the onset of ELN crossings
is marked through the black vertical line to guide the eye.

Figure 6 shows the growth rates for the fast instability (obtained by
calculating the eigenvalues of Eqs. 10 and 11 when $\omega_{\textrm{vac}}=0$,
green curve), along with the growth rates obtained when the vacuum and the
collision terms are included. We find regions of flavor instability for all
time snapshots considered in this work. Interestingly, in all cases except for
$t_{\textrm{pb}}=0.12$ s, the largest growth rate is triggered by the ELN
crossing (green curve), with regions (especially at larger radii) where
instabilities due to $\omega_{\textrm{vac}}\neq 0$ dominate (blue and dashed
red curves), but with a growth rate which is about three orders of magnitude
smaller than the first peak. The behavior of the growth rate for $t_{\rm
pb}=0.12$ s represents a special case (top right panel of Fig. 6); in fact,
the fast instability is absent (green line) in spite of the presence of an ELN
crossing. It is the interplay between the vacuum frequency with the ELN
crossing that makes the system unstable and leads to the growth rate of the
blue curve. This is confirmed by the smaller (slow) growth rate that we
observe in this case with respect to all other snapshots. The fact that the
blue and red curves are exactly on top of each other implies a subleading role
of the collisional instability at $r\gtrsim 56$ km. The effect of the
collisional instability is visible for $\omega_{\textrm{vac}}=0$ (dot-dashed
orange curve).

Figure 7: Zoomed in version of Fig. 6 for $t_{\rm pb}=0.05$, $0.12$, and $1$ s
to highlight the radial region where the collisional instability starts to be
present. In the top and the bottom panels, the blue line is not visible as it
is behind the green line.

We conclude that the growth rate does not seem to be dominated by collisional
effects (dot-dashed orange curve) across the investigated radial range and for
all times, except for a few minor exceptions. In very small regions just prior
to the growth of $\kappa$ of the fast flavor instability, the collisional
instability is found, albeit with a very small growth rate compared to the
fast flavor instability, the collisional instability having the effect of
pushing the occurrence of the flavor instability to a smaller radius by a few
kilometers. Hence, collisions seem to modify the effective vacuum frequency,
in agreement with the findings of Ref. Xiong  _et al._ (2023b). In order to
better inspect the growth of the collisional instability, Fig. 7 shows a
zoomed in version of Fig. 6 for $t_{\rm{pb}}=0.05$, $0.12$, and $1$ s. We can
clearly see a small collisional growth rate that appears just prior to the
fast or slow instability; however, the impact of this small instability is
likely to be negligible. The reason for this is that the flavor instability
develops before neutrinos are completely decoupled from matter, and hence
forward peaked. The flavor evolution at larger radii can thus be transported
to smaller radii because of neutrino advection; hence, the fast or slow flavor
instability might still have a larger impact than collisional instability on
the flavor evolution at smaller radii. These findings suggest that, if the
growth rate of slow and fast instabilities should always be larger than the
collisional one, then the only scenario where the collisional instability
might be important should be when both fast and slow flavor instability are
absent.

It should be emphasized that the presence of a flavor instability in a limited
radial range does not imply that significant flavor conversion develops in or
is restricted to that region; in fact, the advective term can cause the flavor
content at other radii to be modified because of flavor conversion, even if
the flavor instability is not initially present in a certain spatial range
Shalgar  _et al._ (2020); Wu  _et al._ (2021); Nagakura and Zaizen (2023);
Shalgar and Tamborra (2023b, a)

##  V Impact of inhomogeneity on the development of the flavor instability

In this section, we address the effect of inhomogeneities on the results
presented in Sec. IV. To this purpose, we generalize Eqs. 10 and 11 to include
the spatial dependence as follows:

$\displaystyle i\frac{\partial\rho_{ex}}{\partial t}$ $\displaystyle=$
$\displaystyle-i\vec{v}\cdot\vec{\nabla}\rho_{ex}+\left[(H_{ee}-H_{xx})\rho_{ex}-(\rho_{ee}-\rho_{xx})H_{ex}\right]$
(12) $\displaystyle-$
$\displaystyle\frac{i}{2}(\mathcal{C}^{e}_{\textrm{absorb}}+\mathcal{C}^{x}_{\textrm{absorb}})\rho_{ex}$
$\displaystyle-$
$\displaystyle\frac{i}{4}(\mathcal{C}^{e}_{\textrm{dirch}}+\mathcal{C}^{x}_{\textrm{dir-
ch}})$ $\displaystyle\times$
$\displaystyle\int(\rho_{ex}(\cos\theta)-\rho_{ex}(\cos\theta^{\prime}))d\cos\theta^{\prime}\
,$ $\displaystyle i\frac{\partial\bar{\rho}_{ex}}{\partial t}$
$\displaystyle=$
$\displaystyle-i\vec{v}\cdot\vec{\nabla}\bar{\rho}_{ex}+\left[(\bar{H}_{ee}-\bar{H}_{xx})\bar{\rho}_{ex}-(\bar{\rho}_{ee}-\bar{\rho}_{xx})\bar{H}_{ex}\right]$
(13) $\displaystyle-$
$\displaystyle\frac{i}{2}(\bar{\mathcal{C}}^{e}_{\textrm{absorb}}+\bar{\mathcal{C}}^{x}_{\textrm{absorb}})\bar{\rho}_{ex}$
$\displaystyle-$
$\displaystyle\frac{i}{4}(\bar{\mathcal{C}}^{e}_{\textrm{dirch}}+\bar{\mathcal{C}}^{x}_{\textrm{dir-
ch}})$ $\displaystyle\times$
$\displaystyle\int(\bar{\rho}_{ex}(\cos\theta)-\bar{\rho}_{ex}(\cos\theta^{\prime}))d\cos\theta^{\prime}\
.$

It can be seen that the growth rate of $\rho_{ex}$ and $\bar{\rho}_{ex}$ as a
function of time at a given spatial location can be modified by the spatial
derivative only if $\vec{\nabla}\rho_{ex}$ and $\vec{\nabla}\bar{\rho}_{ex}$
are significant in comparison to the other terms entering the kinetic
equations.

If periodic boundaries are assumed, in order to perform the linear stability
analysis for the inhomogeneous modes, the ansatz that $\rho_{ex}$ and
$\bar{\rho}_{ex}$ can be represented as a combination of Fourier modes, i.e.
$\rho_{ex}(\vec{x})=\sum_{\vec{k}}e^{i\vec{k}\cdot\vec{x}}\tilde{\rho}_{ex}^{\vec{k}}$
and
$\bar{\rho}_{ex}(\vec{x})=\sum_{\vec{k}}e^{i\vec{k}\cdot\vec{x}}\bar{\tilde{\rho}}_{ex}^{\vec{k}}$,
has been usually considered in the literature Duan and Shalgar (2015); Abbar
_et al._ (2015). Here, $\tilde{\rho}_{ex}^{\vec{k}}$ and
$\bar{\tilde{\rho}}_{ex}^{\vec{k}}$ are the Fourier components of the off-
diagonal components of the density matrices. Plugging the expressions for
$\rho_{ex}(\vec{x})$ and $\bar{\rho}_{ex}(\vec{x})$ in Eqs. 12 and 13, one can
see that the linearized equations decouple for each Fourier mode, if and only
if the Hamiltonian is independent of the spatial location, i.e. $\rho_{ee}$
and $\bar{\rho}_{ee}$ are not functions of $\vec{x}$. If the Hamiltonian has
non-trivial spatial dependence, as in our case or the one investigated in Ref.
Xiong  _et al._ (2023a), this standard approach cannot be applied to
investigate the growth of the inhomogeneous modes 222Note that the method
adopted in the original papers proposing the linear stability analysis for the
inhomogeneous modes was justified since their goal was to explore the
spontaneous breaking of spherical symmetry, without inhomogeneity in radial
direction Duan and Shalgar (2015); Abbar  _et al._ (2015).

To perform the linear stability analysis in the presence of inhomogeneous
modes, we should seek the eigenvalues of Eqs. 12 and 13 for all spatial points
simultaneously. Since this task is challenging, we look for the eigenvalues
through an iterative procedure. We begin with the ansatz of collective
evolution:

$\displaystyle\rho_{ex}(\vec{x})$ $\displaystyle\sim$ $\displaystyle
Q_{\theta}(r)e^{-i\Omega(r)t}\ ,$ (14) $\displaystyle\bar{\rho}_{ex}(\vec{x})$
$\displaystyle\sim$ $\displaystyle\bar{Q}_{\theta}(r)e^{-i\Omega(r)t}\ .$ (15)

Our goal is to find the functions $\Omega(r)$, $Q_{\theta}(r)$, and
$\bar{Q}_{\theta}(r)$ that satisfy Eqs. 12 and 13. We find the zeroth order
approximation for these functions, ignoring the advective term; note that
$\Omega^{0}(r)$, $Q_{\theta}^{0}(r)$, and $\bar{Q}_{\theta}^{0}(r)$ are
independent of the spatial derivative in the zeroth order approximation. For
each given location, the zeroth order estimate of $\Omega(r)$ denoted by
$\Omega^{0}(r)$ is obtained requiring that

$\displaystyle\begin{vmatrix}I[1]-1&-I[\cos\theta]\\\
I[\cos\theta]&-I[\cos^{2}\theta]-1\end{vmatrix}=0\ ,$ (16)

where

$\displaystyle I[f]=\int
d\cos\theta\frac{\rho_{ee}-\bar{\rho}_{ee}}{\Omega^{0}-H_{ee}}\ .$ (17)

Once the eigenvalue is known, the zeroth order estimate of the eigenvectors
$Q^{0}_{\theta}(r)$ and $\bar{Q}^{0}_{\theta}(r)$ can be computed. With the
knowledge of the approximate dependence of $Q_{\theta}(r)$ and
$\bar{Q}_{\theta}(r)$ on $\cos\theta$ and $r$, it is possible to estimate the
magnitude of the advective term in Eqs. 12 and 13. This helps us to obtain an
iterative estimate of the growth rate solving Eq. 16, but with the denominator
that is modified to include the advective term from the previous iteration.

The advective term applied to Eqs. 14 and 15 gives

$\displaystyle\vec{v}\cdot\nabla\rho_{ex}=\left(\frac{1}{Q_{\theta}(r)}\vec{v}\cdot\nabla
Q_{\theta}(r)-i\frac{d\Omega(r)}{dr}\right)\rho_{ex}\ ;$ (18)

and an analogous expression holds for antineutrinos. The first term in the
parenthesis is generally small compared to the second term, but depends on
$\cos\theta$. The effect of the advective term is to replace $\Omega^{0}$ in
Eq. 17 with $(\Omega^{0}-i{d\Omega}/{dr})$. The real part of ${d\Omega}/{dr}$
is equivalent to adding a matter term that does not affect the growth rate of
the instability. The imaginary part of ${d\Omega}/{dr}$ shifts the imaginary
part of the eigenvalue. This correction is not significant if the derivative
is small compared to $\textrm{Im}(\Omega^{0})$. For all the cases, we find
that the results of this iterative approach lead to results that are
comparable to the ones shown in Fig. 6 and hence are not shown.

Figure 8: Growth of $\left|\int\rho_{ex}~{}d\cos\theta\right|$ as a function
of time at $r=21.3$ km for $t_{\textrm{pb}}=1$ s. The growth rate is given by
the slope of $\left|\int\rho_{ex}~{}d\cos\theta\right|$. The green line has
been obtained by solving the quantum kinetic equations, focusing on the linear
regime. The orange line shows the growth predicted by the linear stability
analysis for the homogeneous mode. The excellent agreement between these two
independent approaches in the linear regime suggest that the inhomogeneous
modes might negligibly affect the growth rate.

The impact of inhomogeneities on the growth rate can also be explored solving
the neutrino quantum kinetic equations numerically in the linear regime and
comparing the latter with the growth rate resulting from the linear stability
analysis carried out considering the homogeneous mode. Figure 8 shows such
comparison at $r=21.3$ km for $t_{\rm{pb}}=1$ s. We can see good agreement. We
have carried out this comparison for all snapshots and a range of radii for
each time snapshot. To this purpose, however, we have used a self-interaction
strength three orders of magnitude smaller than the one computed relying on
the neutrino properties. The attenuation of the self-interaction strength has
been employed because of the technical challenges due to the fact that the
numerical simulation has to be performed for very long time to see the
exponential growth for the very small growth rates typical of certain radii.
In all cases, we find good agreement between the numerical solution of the
quantum kinetic equations and the linear stability analysis carried out for
the homogeneous mode (results not shown here as comparable to the ones
reported in Fig. 8). We note that we also find good agreement between the
results of the linear stability analysis carried out considering the
homogeneous mode (cf. the red dashed line of Fig. 6) and the solution of the
neutrino quantum kinetic equations for $t_{\rm{pb}}=0.12$ s . As shown in Fig.
6, for this snapshot fast instabilities are irrelevant for the homogeneous
mode despite the existence of ELN crossings (cf. green line); because of the
interplay between $\omega_{\rm vac}\neq 0$ and the ELN crossing, the system
becomes unstable, with a growth rate represented by the blue curve. We note
that the employment of the attenuated self-interaction strength may
overestimate the relevance of the vacuum term; nevertheless our main
conclusions would not be qualitatively affected.

Our findings seem to suggest that, in the absence of periodic boundary
conditions, there may be a unique growth rate that is of relevance and it is
given by the solution of Eqs. 12 and 13. This is different from what is
considered in Ref. Xiong  _et al._ (2023a), when the growth rate is calculated
for several Fourier modes; in this case, each Fourier mode has its own growth
rate, the largest of them is interpreted as the physically relevant one.
However, such an analysis assumes self-sustaining perturbations at all
possible length scales, which is possible only if periodic boundary conditions
are imposed. While our findings should be taken with caution due to the
approximations involved, they call for dedicated work to assess the interplay
of inhomogeneities with fast and collisional instabilities. To this purpose, a
method to carry out the linear stability analysis in the presence of
inhomogeneities and in the absence of periodic boundaries should be also
worked out.

##  VI Discussion and outlook

The coherent forward scattering of neutrinos onto each other is responsible
for triggering flavor instabilities that could be fast (if driven by an ELN
crossing for $\omega_{\rm vac}=0$) or slow (if the neutrino vacuum frequency
also plays a role). In addition, more recently, it has been advanced the
possibility that also the collision of neutrinos with the background medium
could be responsible for triggering flavor conversion Johns (2023); Johns and
Xiong (2022); Martin  _et al._ (2021); Lin and Duan (2023); Padilla-Gay  _et
al._ (2022b); Xiong  _et al._ (2023a, b); Kato  _et al._ (2023).

References Shalgar and Tamborra (2023b, a) pointed out the possible co-
existence of slow and fast instabilities in the neutrino decoupling region.
Because of the non-negligible collision rate, in this paper, we aim to
investigate the interplay among slow, fast, and collisional instabilities in
the SN core. To this purpose, we rely on static hydrodynamic backgrounds and
thermodynamic properties extracted at post-bounce times between $0.05$ and $1$
s from a spherically symmetric hydrodynamical simulation of a SN with mass
$18.6\ M_{\odot}$ Archive . We compute the neutrino angular distributions in
the absence of flavor conversion solving the neutrino kinetic equations in the
decoupling region for an average energy mode.

We then carry out the linear stability analysis for the homogeneous mode and
investigate the occurrence of slow, fast, or collisional flavor instabilities.
We find that, due to the ubiquitous existence of ELN crossings, the fast
flavor instability dominates virtually in all cases, except for one time
snapshot ($t_{\rm pb}=0.12$ s) where the slow flavor instability dominates.
Under no circumstances did we find the collisional instability to be the
primary contributor to the flavor instability.

Our findings are in contrast with the analysis presented in Ref. Xiong  _et
al._ (2023a). This could be due to several reasons, such as the application in
Ref. Xiong  _et al._ (2023a) of the attenuation method Nagakura and Zaizen
(2022) to the self-interaction potential as well as an arbitrary rescaling of
the collision term for energies above $90$ MeV, and an implementation of the
collision term that is different from ours. This approach may have led to a
fictitious growth of the collisional instabilities, which should not have
occurred otherwise. Note that the growth rates in Ref. Xiong  _et al._ (2023a)
(cf. their Fig. 3) are strongly time dependent due to rapidly changing number
densities, whereas we assume a classical steady state in which the number
densities of neutrinos do not change in the linear regime.

We have ignored the effects of energy dependent collision terms that can be
important in several ways as pointed out in Ref. Kato  _et al._ (2023).
However, the conclusions reached in Ref. Kato  _et al._ (2023) cannot be
directly extended to our model since the analysis in Ref. Kato  _et al._
(2023) does not assume an initial classical steady state as we do. A
systematic study of collisional instabilities with energy dependent collision
terms is left to future work. The role of inhomogeneities in the development
of flavor instabilities, for which we only present preliminary findings in
this paper, should also be subject of future dedicated work.

In conclusion, collisional instabilities remain an intriguing possibility to
trigger flavor conversion; however, it is not clear whether in neutrino-dense
astrophysical environments, flavor instabilities due to collisions could drive
and sustain flavor conversion. Further work will be needed to analyze the
development of flavor instabilities in a large set of SN models, as well as
neutron star merger simulations, exhibiting a wide range of thermodynamic and
neutrino properties.

###### Acknowledgements.

We thank Thomas Janka and Evan O’Connor for insightful discussions. We are
also grateful to Thomas Janka, Lucas Johns, Zewei Xiong, and Meng-Ru Wu for
helpful comments on our manuscript. This project has received support from the
Villum Foundation (Project No. 13164), the Danmarks Frie Forskningsfonds
(Project No. 8049-00038B), and the Deutsche Forschungsgemeinschaft through
Sonderforschungbereich SFB 1258 “Neutrinos and Dark Matter in Astro- and
Particle Physics” (NDM).

##  Appendix A Collisional terms

Figure 9: Radial profile of the first energy moment for $\nu_{e}$ (top left),
$\bar{\nu}_{e}$ (top right) and $\nu_{x}$ (bottom left) for all selected post-
bounce times $t_{pb}=0.05$, $0.12$, $0.25$, $0.5$, $0.75$, and $1$ s extracted
for our benchmark hydrodynamical SN simulation. The bottom right panel
represents the flavor-averaged first energy moment that we use as
characteristic neutrino energy in the linear stability analysis.

In the interior of a core-collapse SN, several processes contribute to the
production, absorption, and interaction of neutrinos with matter. The most
relevant ones are Bremsstrahlung, pair-production, beta reactions, and
direction changing processes. In this appendix, we present the implementation
of these terms in the neutrino kinetic equations (see Eqs. 7 and 8). Our
modeling of the collisional kernel was heavily inspired by the open-source
neutrino transport code NuLib O’Connor (2015).

It should be noted that there are other thermal processes occurring in the SN
interior which we ignore, namely, plasmon decay, recombination, and photo-
production as they are subdominant Itoh  _et al._ (1996). In addition, we
neglect the energy changing non-elastic neutrino-electron scattering because
it is subdominant, especially for the electron type-neutrinos which determine
the development of the instabilities of interest to this work. We do not
consider charge current interactions on muons since our benchmark SN
simulation does not include them for simplicity, despite the fact that muons
could be important Bollig  _et al._ (2017).

For calculating beta processes and direction changing rates, we use the
flavor-dependent average energy of neutrinos extracted from our benchmark
hydrodynamical SN simulation and displayed in Fig. 9, whereas for thermal
processes (i.e. Bremmstrahlung and pair processes) we calculate the production
by integrating over the energy dependent emission rates, and calculate the
total absorption rate by applying Kirchoff’s law to the energy integrated
production rate.

###  A.1 Link between the absorption and emission rates

The absorption and emission rates are not independent in the steady state
configuration, but are related by the Kirchoff’s law of radiation. As a
consequence of the Kirchoff’s law, the emission rate is such that the number
density reaches a thermal distribution asymptotically in the absence of
advection. The Kirchoff’s law is not just valid in thermal equilibrium, but
also in a steady state. Even in the low density region where the distribution
of particles is far from thermal, the Kirchoff’s law can be used to relate the
absorption term to the emission one Thompson (2002):

$\displaystyle\mathcal{L}^{\textrm{int}}[\rho_{ii}]=\eta_{\nu}(1-\rho_{ii})-\chi_{\nu}\rho_{ii}\
,$ (19)

where $\chi_{\nu}$ is the absorption rate and $\eta_{\nu}$ is the emission
rate without Pauli blocking. The right hand side of Eq. 19 is zero in
equilibrium or equivalently:

$\displaystyle\eta_{\nu}=\chi_{\nu}\frac{\rho_{ii}^{\textrm{eq}}}{1-\rho_{ii}^{\textrm{eq}}}\
,$ (20)

with $\rho_{ii}^{\textrm{eq}}$ being the Fermi-Dirac distribution of neutrinos
Substituting, Eq. 20 in Eq. 19, we obtain

$\displaystyle\mathcal{L}^{\textrm{int}}[\rho_{ii}]=\mathcal{C}^{ii}_{\textrm{absorb}}[\rho_{ii}^{\textrm{eq}}-\rho_{ii}]\
,$ (21)

where
$\mathcal{C}^{ii}_{\textrm{absorb}}=\chi_{\nu}/(1-\rho_{ii}^{\textrm{eq}})$ is
the absorption rate with the denominator included to account for Pauli
blocking.

In the following, we discuss the emission and absorption terms. As for
$\nu_{e}$ and $\bar{\nu}_{e}$, charged current interactions (beta reactions)
dominate over the neutral current processes, but neutral current processes are
important for heavy lepton neutrinos.

###  A.2 Emission and absorption processes: Beta reactions

Beta reactions are the dominant channels determining the number density of
electron type neutrinos. The absorption and emission of $\nu_{e}$ and
$\bar{\nu}_{e}$ occurs through the following interactions with free nucleons:

$\displaystyle\nu_{e}+n$ $\displaystyle\leftrightarrow$ $\displaystyle
p+e^{-}$ (22) $\displaystyle\bar{\nu}_{e}+p$ $\displaystyle\leftrightarrow$
$\displaystyle n+e^{+}\ .$ (23)

The correction to the beta reaction rates due to the presence of heavier
nuclei is subdominant and can be neglected for our purposes.

The absorption cross section of electron neutrinos on neutrons is given by
Bruenn (1985):

$\displaystyle\sigma^{a}_{\nu_{e}n}$ $\displaystyle=$
$\displaystyle\sigma_{0}\left(\frac{1+3g_{A}^{2}}{4}\right)\left(\frac{\langle\epsilon_{\nu_{e}}\rangle+\Delta_{np}}{m_{e}c^{2}}\right)^{2}$
(24) $\displaystyle\times$
$\displaystyle\left[1-\left(\frac{m_{e}c^{2}}{\langle\epsilon_{\nu_{e}}\rangle+\Delta_{np}}\right)^{2}\right]^{\frac{1}{2}}(1-f_{e^{-}})W_{M}\
.$

Here, $\sigma_{0}$ is the characteristic neutrino cross section,

$\displaystyle\sigma_{0}=\frac{4G_{\textrm{F}}\left(m_{e}c^{2}\right)^{2}}{\pi(\hbar
c)^{4}}=1.7612\times 10^{-44}\textrm{cm}^{2}\ ,$ (25)

with $g_{A}=-1.254$, $\Delta_{np}$ being the difference in the mass of the
neutron and proton ($m_{n}-m_{p}=1.2933$ MeV), $m_{e}$ the mass of the
electron, $\langle\epsilon_{\nu_{e}}\rangle$ the $\nu_{e}$ average energy
displayed in Fig. 9, $f_{e^{-}}$ the Fermi-Dirac distribution of electrons,
and $W_{M}=1-1.01\langle\epsilon_{\nu_{e}}\rangle/m_{n}$ the weak magnetic
correction Vogel (1984); Thompson (2002).

The related inverse transport mean free path of neutrinos is

$\displaystyle\mathcal{C}_{\textrm{absorb}}=\sigma^{a}_{\nu_{e}n}X_{n}\rho_{B}N_{A}(1-\rho_{ee}^{\textrm{eq}})^{-1}\
,$ (26)

where $\sigma^{a}_{\nu_{e}n}$ is given by Eq. 24, and $X_{n}$ and $\rho_{B}$
are the neutron fraction and the baryon density shown in Fig. 2. In the
expression above, we ignore the Pauli blocking of nucleons as it is
subdominant in the decoupling and free streaming region Bruenn (1985); Raffelt
.

A similar expression holds for $\bar{\nu}_{e}$ absorption with the exception
that we ignore the Pauli blocking for the outgoing positron and with other
appropriate changes following Ref. Bruenn (1985):

$\displaystyle\sigma^{a}_{\bar{\nu}_{e}p}$ $\displaystyle=$
$\displaystyle\sigma_{0}\left(\frac{1+3g_{A}^{2}}{4}\right)\left(\frac{\langle\epsilon_{\bar{\nu}_{e}}\rangle-\Delta_{np}}{m_{e}c^{2}}\right)^{2}$
(27) $\displaystyle\times$
$\displaystyle\left[1-\left(\frac{m_{e}c^{2}}{\langle\epsilon_{\bar{\nu}_{e}}\rangle-\Delta_{np}}\right)^{2}\right]^{\frac{1}{2}}W_{\bar{M}}\
,$

with $W_{\bar{M}}=(1-7.1\langle\epsilon_{\bar{\nu}_{e}}\rangle/m_{p})$. The
inverse of the transport mean free path is defined as

$\displaystyle\bar{\mathcal{C}}_{\textrm{absorb}}=\sigma^{a}_{\bar{\nu}_{e}p}X_{p}\rho_{B}N_{A}(1-\bar{\rho}_{ee}^{\textrm{eq}})^{-1}\
.$ (28)

The absorption mean free path can be used to calculate the emission rate
relying on the Kirchoff’s law:

$\displaystyle\mathcal{C}_{\textrm{emit}}$ $\displaystyle=$
$\displaystyle\mathcal{C}_{\textrm{absorb}}\rho_{ii}^{\textrm{eq}}$ (29)
$\displaystyle\bar{\mathcal{C}}_{\textrm{emit}}$ $\displaystyle=$
$\displaystyle\bar{\mathcal{C}}_{\textrm{absorb}}\bar{\rho}_{ii}^{\textrm{eq}}$
(30)

It should be noted that the equality above holds for a steady state
configuration and it is not necessary to have thermal equilibrium. Hence, we
can use the Kirchoff’s law in the trapped regime, decoupling regime, and the
free streaming regime. The top panels of Fig. 10 display the radial evolution
of the absorption rates due to beta processes for electron neutrinos and
antineutrinos for our selected post-bounce snapshots.

Figure 10: Radial evolution of the absorption rates due to beta reactions for
$\nu_{e}$ and $\bar{\nu}_{e}$ (top row; the emission rate is linked to the
absorption one through the Kirchoff’s law), pair production (second row),
Bremsstrahlung (third row), and direction changing neutral current scattering
processes (bottom row) that we compute at various post-bounce time snapshots
and for $\nu_{e}$ (left panels), $\bar{\nu}_{e}$ (middle panels), and
$\nu_{x}$ (right panels). The difference in the interaction rate for the
different flavors is due to the different average energy of each flavor. The
emission rates are related to the absorption rates by the Kirchoff’s law. At
small $r$, beta processes are dominant. At larger radii, direction-changing
reactions start to become competitive with beta processes and are more
efficient than pair and bremsstrahlung reactions, with the latter rapidly
decreasing as a function of time.

###  A.3 Emission and absorption processes: Pair production

Electrons and positrons annihilate to produce pairs of neutrinos and
antineutrinos:

$\displaystyle e^{-}+e^{+}\leftrightarrow\nu_{i}+\bar{\nu}_{i}\quad i=e,x\ .$
(31)

The $s$-channel Feynman diagrams can create pairs of neutrinos of all flavors,
along with the $t$-channel diagrams that can only produce pairs of electron
type neutrinos. Note that the electron-positron pairs responsible for neutrino
production are relativistic in the SN interior and, hence, the production
cross section depends on the momenta of the initial pair.

To compute the production rate of the neutrinos, the production kernel can be
decomposed in Legendre polynomials as described in Refs. Bruenn (1985).
Assuming that the initial electron-positron pairs are isotropically
distributed, we use Eq. C63 of Ref. Bruenn (1985) to calculate the production
rate for all neutrino species. Note that the produced neutrinos do not have an
isotropic distribution.

The emission rate can be calculated from the absorption rate. However, there
are some caveats that need to be addressed. Unlike the beta processes,
absorption occurs in pairs. This implies that the $\nu_{e}$ absorption rate
due to pair processes is overestimated in the region where there are far more
$\nu_{e}$’s than $\bar{\nu}_{e}$’s. This effect is not very important for
electron type neutrinos as the beta reaction rates dominate over the pair
processes. Hence we ignore the implications of such an approximation, but it
is possible to add a correction to this overestimation of the absorption rate
due to pair process Betranhandy (2022). The second row of Fig. 10 shows the
absorption rates integrated over the neutrino energy distribution for all
three neutrino flavors for different snapshots.

###  A.4 Emission and absorption processes: Bremsstrahlung

The reaction $n+n\leftrightarrow n+n+\nu+\bar{\nu}$ is a subdominant thermal
production process responsible for the creation of all flavors of neutrinos;
as a consequence, it has virtually no impact on the neutrino distributions of
all flavors, nevertheless, we include it in our calculations. Even in the case
of heavy lepton neutrinos, the bremsstrahlung channel can be ignored in most
cases, as visible from Fig. 10.

For estimating the production and absorption rates of neutrinos through
Bremsstrahlung we use the following approximate formula Burrows  _et al._
(2006):

$\displaystyle\frac{dQ^{\prime}_{\nu}}{d\epsilon_{\nu}}$ $\displaystyle=$
$\displaystyle\frac{0.234Q_{nb}}{T}\left(\frac{\epsilon_{\nu}}{T}\right)^{2.4}e^{-1.1\epsilon_{\nu}/T}\
,$ (32) $\displaystyle Q_{nb}$ $\displaystyle=$ $\displaystyle 2.0778\times
10^{30}\zeta
X_{n}^{2}\rho_{B,14}^{2}\left(\frac{T}{\textrm{MeV}}\right)^{5.5}\ ,$ (33)

with $T$ being the temperature, $\rho_{B,14}$ the baryon density in units of
$10^{14}$ g/cm3, $X_{n}$ the neutron fraction (see Fig. 2), and $\zeta=0.5$
the correction factor 333It should be noted we follow Ref. O’Connor (2015) and
Eq. 33 differs from Eq. 142 of Ref. Burrows  _et al._ (2006) because of an
error in that paper O’Connor .  Equation 33 is accurate up to $3\%$ Thompson
(2002). The emission term associated with this process is obtained by
integrating the emissivity over the energy range,

$\displaystyle\mathcal{C}_{\textrm{emit}}=\int
d\epsilon_{\nu}\frac{dQ^{\prime}_{\nu}}{d\epsilon_{\nu}}\ .$ (34)

The absorption term is related to the emission term through the Kirchoff’s
law,

$\displaystyle\mathcal{C}_{\textrm{absorb}}=\frac{\mathcal{C}_{\textrm{emit}}}{n_{i}^{\textrm{eq}}}$
(35)

which is shown in the third row of Fig. 10.

###  A.5 Direction changing processes: Neutral current interactions

The momentum changing interaction of neutrinos with matter can be due to
neutral current interactions as well as charged current interactions. As for
the latter, only the electron flavors ($\nu_{e}$ and $\bar{\nu}_{e}$) can
participate in charged current interactions with electrons; yet, this momentum
changing process is subdominant with respect to the neutral current scattering
of neutrinos with nucleons. Therefore, in the following, we focus on neutral
current processes.

We calculate the momentum changing cross section in the single energy
approximation (relying on the flavor-dependent neutrino average energy shown
in Fig. 9), so the momentum changing cross section can be viewed as a
direction changing one. The direction changing cross section also depends on
the matter composition. However, for our purposes, the contribution of nuclei
(and related nuclear matrix elements) can be ignored. Hence we focus on the
interaction of neutrinos with protons and neutrons. We also ignore the neutral
current interactions of neutrinos with electrons since they are subdominant
with respect to the interaction with nucleons.

As the interaction of neutrinos with nucleons is independent on the neutrino
flavor, the transport cross section with neutrons and protons is given by
Thompson (2002):

$\displaystyle\sigma_{n}^{\textrm{tr}}$ $\displaystyle=$
$\displaystyle\frac{\sigma_{0}}{4}\left(\frac{\langle\epsilon_{\nu}\rangle}{m_{e}c^{2}}\right)^{2}\left(\frac{1+5g_{A}^{2}}{6}\right)\
,$ (36) $\displaystyle\sigma_{p}^{\textrm{tr}}$ $\displaystyle=$
$\displaystyle\frac{\sigma_{0}}{6}\left(\frac{\langle\epsilon_{\nu}\rangle}{m_{e}c^{2}}\right)^{2}\left[(C_{V}^{\prime}-1)^{2}+5g_{A}(C_{A}^{\prime}-1)^{2}\right]\
,$

where $\sigma_{0}$ is defined as in Eq. 25, $g_{A}=1.254$,
$C_{V}^{\prime}=1/2+2\sin^{2}\theta_{W}$ (with $\theta_{W}$ being the Weinberg
angle), and $C_{A}^{\prime}=1/2$. The direction changing collision term
$C_{\textrm{dir-ch}}$ as follows:

$\displaystyle\mathcal{C}_{\textrm{dir-
ch}}=\rho_{B}N_{A}(X_{n}\sigma_{n}^{\textrm{tr}}+X_{p}\sigma_{p}^{\textrm{tr}})\
.$ (38)

It should be noted that the direction changing neutral current interactions
are not isotropic, but backward peaked. The neutrino interactions with protons
have an anisotropy of about 10%
($\mathcal{C}_{\textrm{ani}}=0.1~{}\mathcal{C}_{\textrm{dir-ch}}$) whereas the
interactions with neutrons have an anisotropy of about 33%
($\mathcal{C}_{\textrm{ani}}=0.33~{}\mathcal{C}_{\textrm{dir-ch}}$). This is
taken into account through $\mathcal{C}_{\rm{ani}}$, which is the weighted
average of the two anisotropies in Eqs. 7 and 8. The bottom panels of Fig. 10
illustrate the dependence of the direction-chancing processes on the radius;
note that the flavor dependence of the neutral current interaction rate in
Fig. 10 is due to the fact that the average energy of neutrinos depends on the
neutrino flavor.

##  Appendix B Comparison of the neutrino number densities with the output of
our benchmark supernova model

Relying on the hydrodynamical and thermodynamic quantities, local neutrino
number density, and first energy moment extracted at selected post-bounce
times from our benchmark SN simulation, we compute the neutrino angular
distributions, following the modeling of the collision term outlined in
Appendix A. In this appendix, we compare the flavor-dependent and angle-
integrated neutrino number densities that we obtain from post-processing the
output of our $18.6\ M_{\odot}$ SN model Archive  with the neutrino number
density extracted from the hydrodynamical simulation to assess the degree of
accuracy of our collision term. The comparison is shown in Fig. 11. Despite
the good agreement between the neutrino number densities that we compute in
post-processing and the outcome of hydrodynamical simulations, we stress that
the purpose of this work is not to reproduce the neutrino properties from
hydrodynamical simulations, but rather rely on angular distributions derived
from first principles to investigate the relative importance of collisional,
fast, and slow flavor instabilities.

Figure 11: Radial evolution of the neutrino number densities computed in this
work (solid lines) of $\nu_{e}$ (in red), $\bar{\nu}_{e}$ (in green), and
$\nu_{x}$ (in blue) along with the ones extracted from our benchmark SN
simulation (denoted by “SN model” in the legend) for selected post-bounce time
snapshots. The radius at which the number density of each neutrino flavor
deviates from the Fermi-Dirac one corresponds to the radius of neutrino
decoupling marked by a vertical line. While involving some approximations, our
method allows us to reproduce the local neutrino density which is in general
agreement the output of the hydrodynamical simulation. Note that the radial
range is not the same for all time snapshots.

In the SN core, neutrinos are continuously produced and absorbed at the same
time (see Appendix A). These two processes eventually lead to thermal
equilibration. However, the baryon density falls rapidly as a function of
radius (see Fig. 2), leading neutrinos to go out of thermal equilibrium and
stream freely. The radius of decoupling at which this happens depends on the
magnitude of the production and absorption rates compared to the direction
changing term. As shown in Fig. 11, this implies that the neutrino number
density coincides with the number density obtained by integrating Fermi-Dirac
distribution over energy in the dense region, and it deviates from the energy-
integrated Fermi-Dirac number density at larger radii where the matter density
is lower, and the angular distribution becomes forward peaked. In order to
assess the location of the neutrino decoupling surfaces in Fig. 3, we compare
the local neutrino number density with the one obtained from the energy-
integrated Fermi-Dirac neutrino density and define the neutrino decoupling
radius as the one where the neutrino number density deviates from the Fermi-
Dirac one more than $1\%$ (see vertical lines in Fig. 11).

Figure 12: Top: Comparison between the neutrino number densities obtained
employing multi-energy dependence in the collisional rate (solid lines) and
the ones coming from our benchmark SN simulation (dotted lines) for
$t_{\textrm{pb}}=1$ s. The already small descrepancy between the number
densities seen in Fig. 11 is further alleviated by the multi-energy treatment.
Bottom: Growth rate of the flavor instability obtained employing multi-energy
dependence in the collisional rate. The qualitative feature of collisional
instability not being important remains unchanged by the single energy
approximation.

It should be noted that in Fig. 11 the $\bar{\nu}_{e}$ number density is
slightly larger than the one from the “SN model” because of the single energy
approximation employed in this work. We demonstrate this by performing an
energy-dependent computation and report the results in Fig. 12 for the
snapshot $t_{\textrm{pb}}$=1 s. The top panel of Fig. 12 shows that the
discrepancy between the number densities obtained in this work and the ones
from the Garching simulation is alleviated almost completely. The remaining
negligible residual difference is due to the fact that our implementation of
the collisional term is not identical to the one employed in the Garching
simulation.

The bottom panel of Fig. 12 shows that the growth rate obtained performing a
multi-energy computation remains qualitatively unchanged (cf. bottom right
panel of Fig. 6). It should be noted, however, that there are some
quantitative differences between the single energy and multi-energy
calculations: in the multi-energy simulation the maximum growth rate is
slightly smaller and the region of instability is shifted to slightly larger
radii. However, this does not affect our main conclusions. Figure 12, hence,
demonstrates that the single-energy approximation is justified in our case.

We note that the green dotted line (i.e., $\bar{\nu}_{e}$) is below the red
one (i.e. $\nu_{e}$) in Fig. 11; on the other hand, the solid and red lines
cross at $\sim 22.6$ km and the $\bar{\nu}_{e}$ number density is larger than
the $\nu_{e}$’s one at larger radii. Such crossing between the $\nu_{e}$ and
$\bar{\nu}_{e}$ number densities is also present in our multi-energy
simulations in Fig. 12, although much less prominent, and it might give rise
to differences in the development of ELN crossings. We conclude that the
single energy approximation slightly overestimates the number density
$\bar{\nu}_{e}$ facilitating the presence of a sharper ELN crossing. However,
although shallower, an ELN crossing is also expected in the multi-energy
calculation (cf. Fig. 12). Thus implies that our single-energy treatment does
not affects the relative importance of fast vs. collisional instabilities–the
latter being dominant only in the absence of ELN crossings. In fact, a
comparison between the bottom panel of Fig. 12 and the bottom right panel of
Fig. 6 shows that the growth rate for fast instabilities already develops at
smaller radii (the ELN crossing developing at $20.5$ km for the single-energy
simulation vs. $21.1$ km for the multi-energy simulations) and is
qualitatively unchanged by this feature.

## References

  * Mirizzi  _et al._ (2016) Alessandro Mirizzi, Irene Tamborra, Hans-Thomas Janka, Ninetta Saviano, Kate Scholberg, Robert Bollig, Lorenz Hüdepohl, and Sovan Chakraborty, “Supernova Neutrinos: Production, Oscillations and Detection,” Riv. Nuovo Cim. 39, 1–112 (2016), arXiv:1508.00785 [astro-ph.HE] .
  * Duan  _et al._ (2010) Huaiyu Duan, George M. Fuller, and Yong-Zhong Qian, “Collective Neutrino Oscillations,” Ann. Rev. Nucl. Part. Sci. 60, 569–594 (2010), arXiv:1001.2799 [hep-ph] .
  * Chakraborty  _et al._ (2016a) Sovan Chakraborty, Rasmus Hansen, Ignacio Izaguirre, and Georg G. Raffelt, “Collective neutrino flavor conversion: Recent developments,” Nucl. Phys. B 908, 366–381 (2016a), arXiv:1602.02766 [hep-ph] .
  * Tamborra and Shalgar (2021) Irene Tamborra and Shashank Shalgar, “New Developments in Flavor Evolution of a Dense Neutrino Gas,” Ann. Rev. Nucl. Part. Sci. 71, 165–188 (2021), arXiv:2011.01948 [astro-ph.HE] .
  * Richers and Sen (2022) Sherwood Richers and Manibrata Sen, “Fast Flavor Transformations,” (2022), arXiv:2207.03561 [astro-ph.HE] .
  * Sawyer (2005) Raymond F. Sawyer, “Speed-up of neutrino transformations in a supernova environment,” Phys. Rev. D 72, 045003 (2005), arXiv:hep-ph/0503013 .
  * Sawyer (2009) Raymond F. Sawyer, “The multi-angle instability in dense neutrino systems,” Phys. Rev. D 79, 105003 (2009), arXiv:0803.4319 [astro-ph] .
  * Sawyer (2016) Raymond F. Sawyer, “Neutrino cloud instabilities just above the neutrino sphere of a supernova,” Phys. Rev. Lett. 116, 081101 (2016), arXiv:1509.03323 [astro-ph.HE] .
  * Chakraborty  _et al._ (2016b) Sovan Chakraborty, Rasmus Sloth Hansen, Ignacio Izaguirre, and Georg G. Raffelt, “Self-induced neutrino flavor conversion without flavor mixing,” JCAP 03, 042 (2016b), arXiv:1602.00698 [hep-ph] .
  * Izaguirre  _et al._ (2017) Ignacio Izaguirre, Georg G. Raffelt, and Irene Tamborra, “Fast Pairwise Conversion of Supernova Neutrinos: A Dispersion-Relation Approach,” Phys. Rev. Lett. 118, 021101 (2017), arXiv:1610.01612 [hep-ph] .
  * Padilla-Gay  _et al._ (2022a) Ian Padilla-Gay, Irene Tamborra, and Georg G. Raffelt, “Neutrino Flavor Pendulum Reloaded: The Case of Fast Pairwise Conversion,” Phys. Rev. Lett. 128, 121102 (2022a), arXiv:2109.14627 [astro-ph.HE] .
  * Morinaga (2022) Taiki Morinaga, “Fast neutrino flavor instability and neutrino flavor lepton number crossings,” Phys. Rev. D 105, L101301 (2022), arXiv:2103.15267 [hep-ph] .
  * Shalgar and Tamborra (2023a) Shashank Shalgar and Irene Tamborra, “Neutrino decoupling is altered by flavor conversion,” Phys. Rev. D 108, 043006 (2023a), arXiv:2206.00676 [astro-ph.HE] .
  * Shalgar and Tamborra (2023b) Shashank Shalgar and Irene Tamborra, “Neutrino flavor conversion, advection, and collisions: Toward the full solution,” Phys. Rev. D 107, 063025 (2023b), arXiv:2207.04058 [astro-ph.HE] .
  * Shalgar and Tamborra (2021) Shashank Shalgar and Irene Tamborra, “A change of direction in pairwise neutrino conversion physics: The effect of collisions,” Phys. Rev. D 103, 063002 (2021), arXiv:2011.00004 [astro-ph.HE] .
  * Johns (2023) Lucas Johns, “Collisional Flavor Instabilities of Supernova Neutrinos,” Phys. Rev. Lett. 130, 191001 (2023), arXiv:2104.11369 [hep-ph] .
  * Sasaki and Takiwaki (2022) Hirokazu Sasaki and Tomoya Takiwaki, “A detailed analysis of the dynamics of fast neutrino flavor conversions with scattering effects,” PTEP 2022, 073E01 (2022), arXiv:2109.14011 [hep-ph] .
  * Hansen  _et al._ (2022) Rasmus S. L. Hansen, Shashank Shalgar, and Irene Tamborra, “Enhancement or damping of fast neutrino flavor conversions due to collisions,” Phys. Rev. D 105, 123003 (2022), arXiv:2204.11873 [astro-ph.HE] .
  * Kato and Nagakura (2022) Chinami Kato and Hiroki Nagakura, “Effects of energy-dependent scatterings on fast neutrino flavor conversions,” Phys. Rev. D 106, 123013 (2022), arXiv:2207.09496 [astro-ph.HE] .
  * Johns and Xiong (2022) Lucas Johns and Zewei Xiong, “Collisional instabilities of neutrinos and their interplay with fast flavor conversion in compact objects,” Phys. Rev. D 106, 103029 (2022), arXiv:2208.11059 [hep-ph] .
  * Padilla-Gay  _et al._ (2022b) Ian Padilla-Gay, Irene Tamborra, and Georg G. Raffelt, “Neutrino fast flavor pendulum. II. Collisional damping,” Phys. Rev. D 106, 103031 (2022b), arXiv:2209.11235 [hep-ph] .
  * Lin and Duan (2023) Yu-Chia Lin and Huaiyu Duan, “Collision-induced flavor instability in dense neutrino gases with energy-dependent scattering,” Phys. Rev. D 107, 083034 (2023), arXiv:2210.09218 [hep-ph] .
  * Liu  _et al._ (2023) Jiabao Liu, Masamichi Zaizen, and Shoichi Yamada, “Systematic study of the resonancelike structure in the collisional flavor instability of neutrinos,” Phys. Rev. D 107, 123011 (2023), arXiv:2302.06263 [hep-ph] .
  * Kato  _et al._ (2023) Chinami Kato, Hiroki Nagakura, and Masamichi Zaizen, “Flavor conversions with energy-dependent neutrino emission and absorption,” Phys. Rev. D 108, 023006 (2023), arXiv:2303.16453 [astro-ph.HE] .
  * Burrows  _et al._ (2006) Adam Burrows, Sanjay Reddy, and Todd A. Thompson, “Neutrino opacities in nuclear matter,” Nucl. Phys. A 777, 356–394 (2006), arXiv:astro-ph/0404432 .
  * Janka (2012) Hans-Thomas Janka, “Explosion Mechanisms of Core-Collapse Supernovae,” Ann. Rev. Nucl. Part. Sci. 62, 407–451 (2012), arXiv:1206.2503 [astro-ph.SR] .
  * Mezzacappa  _et al._ (2020) Anthony Mezzacappa, Eirik Endeve, O. E. Bronson Messer, and Stephen W. Bruenn, “Physical, numerical, and computational challenges of modeling neutrino transport in core-collapse supernovae,” Liv. Rev. Comput. Astrophys. 6, 4 (2020), arXiv:2010.09013 [astro-ph.HE] .
  * Shalgar and Tamborra (2019) Shashank Shalgar and Irene Tamborra, “On the Occurrence of Crossings Between the Angular Distributions of Electron Neutrinos and Antineutrinos in the Supernova Core,” Astrophys. J. 883, 80 (2019), arXiv:1904.07236 [astro-ph.HE] .
  * Nagakura (2023) Hiroki Nagakura, “Roles of Fast Neutrino-Flavor Conversion on the Neutrino-Heating Mechanism of Core-Collapse Supernova,” Phys. Rev. Lett. 130, 211401 (2023), arXiv:2301.10785 [astro-ph.HE] .
  * Xiong  _et al._ (2023a) Zewei Xiong, Meng-Ru Wu, Gabriel Martínez-Pinedo, Tobias Fischer, Manu George, Chun-Yu Lin, and Lucas Johns, “Evolution of collisional neutrino flavor instabilities in spherically symmetric supernova models,” Phys. Rev. D 107, 083016 (2023a), arXiv:2210.08254 [astro-ph.HE] .
  * Martin  _et al._ (2021) Joshua D. Martin, Joseph A. Carlson, Vincenzo Cirigliano, and Huaiyu Duan, “Fast flavor oscillations in dense neutrino media with collisions,” Phys. Rev. D 103, 063001 (2021), arXiv:2101.01278 [hep-ph] .
  * Xiong  _et al._ (2023b) Zewei Xiong, Lucas Johns, Meng-Ru Wu, and Huaiyu Duan, “Collisional flavor instability in dense neutrino gases,” Phys. Rev. D 108, 083002 (2023b), arXiv:2212.03750 [hep-ph] .
  * (33) MPA Supernova Archive, https://wwwmpa.mpa-garching.mpg.de/ccsnarchive https://wwwmpa.mpa-garching.mpg.de/ccsnarchive/data/Bollig2016/.
  * Richers  _et al._ (2019) Sherwood A. Richers, Gail C. McLaughlin, James P. Kneller, and Alexey Vlasenko, “Neutrino Quantum Kinetics in Compact Objects,” Phys. Rev. D 99, 123014 (2019), arXiv:1903.00022 [astro-ph.HE] .
  * Tamborra  _et al._ (2017) Irene Tamborra, Lorenz Hüdepohl, Georg G. Raffelt, and Hans-Thomas Janka, “Flavor-dependent neutrino angular distribution in core-collapse supernovae,” Astrophys. J. 839, 132 (2017), arXiv:1702.00060 [astro-ph.HE] .
  * Morinaga  _et al._ (2020) Taiki Morinaga, Hiroki Nagakura, Chinami Kato, and Shoichi Yamada, “Fast neutrino-flavor conversion in the preshock region of core-collapse supernovae,” Phys. Rev. Res. 2, 012046 (2020), arXiv:1909.13131 [astro-ph.HE] .
  * Abbar  _et al._ (2020) Sajad Abbar, Huaiyu Duan, Kohsuke Sumiyoshi, Tomoya Takiwaki, and Maria Cristina Volpe, “Fast Neutrino Flavor Conversion Modes in Multidimensional Core-collapse Supernova Models: the Role of the Asymmetric Neutrino Distributions,” Phys. Rev. D 101, 043016 (2020), arXiv:1911.01983 [astro-ph.HE] .
  * Abbar  _et al._ (2019) Sajad Abbar, Huaiyu Duan, Kohsuke Sumiyoshi, Tomoya Takiwaki, and Maria Cristina Volpe, “On the occurrence of fast neutrino flavor conversions in multidimensional supernova models,” Phys. Rev. D 100, 043004 (2019), arXiv:1812.06883 [astro-ph.HE] .
  * Delfan Azari  _et al._ (2020) Milad Delfan Azari, Shoichi Yamada, Taiki Morinaga, Hiroki Nagakura, Shun Furusawa, Akira Harada, Hirotada Okawa, Wakana Iwakami, and Kohsuke Sumiyoshi, “Fast collective neutrino oscillations inside the neutrino sphere in core-collapse supernovae,” Phys. Rev. D 101, 023018 (2020), arXiv:1910.06176 [astro-ph.HE] .
  * Nagakura  _et al._ (2019) Hiroki Nagakura, Taiki Morinaga, Chinami Kato, and Shoichi Yamada, “Fast-pairwise Collective Neutrino Oscillations Associated with Asymmetric Neutrino Emissions in Core-collapse Supernovae,” Astrophys. J. 886, 139 (2019), arXiv:1910.04288 [astro-ph.HE] .
  * Abbar  _et al._ (2021) Sajad Abbar, Francesco Capozzi, Robert Glas, Hans-Thomas Janka, and Irene Tamborra, “On the characteristics of fast neutrino flavor instabilities in three-dimensional core-collapse supernova models,” Phys. Rev. D 103, 063033 (2021), arXiv:2012.06594 [astro-ph.HE] .
  * (42) Shashank Shalgar et al., in preparation.
  * Weiss  _et al._ (2004) A. Weiss, W. Hillebrandt, H. C. Thomas, and H. Ritter,  _Cox and Giuli’s Principles of Stellar Structure_ (2004).
  * Banerjee  _et al._ (2011) Arka Banerjee, Amol Dighe, and Georg G. Raffelt, “Linearized flavor-stability analysis of dense neutrino streams,” Phys. Rev. D84, 053013 (2011), arXiv:1107.2308 [hep-ph] .
  * Airen  _et al._ (2018) Sagar Airen, Francesco Capozzi, Sovan Chakraborty, Basudeb Dasgupta, Georg Raffelt, and Tobias Stirner, “Normal-mode Analysis for Collective Neutrino Oscillations,” JCAP 12, 019 (2018), arXiv:1809.09137 [hep-ph] .
  * Morinaga and Yamada (2018) Taiki Morinaga and Shoichi Yamada, “Linear stability analysis of collective neutrino oscillations without spurious modes,” Phys. Rev. D 97, 023024 (2018), arXiv:1803.05913 [hep-ph] .
  * Shalgar  _et al._ (2020) Shashank Shalgar, Ian Padilla-Gay, and Irene Tamborra, “Neutrino propagation hinders fast pairwise flavor conversions,” JCAP 06, 048 (2020), arXiv:1911.09110 [astro-ph.HE] .
  * Wu  _et al._ (2021) Meng-Ru Wu, Manu George, Chun-Yu Lin, and Zewei Xiong, “Collective fast neutrino flavor conversions in a 1D box: Initial conditions and long-term evolution,” Phys. Rev. D 104, 103003 (2021), arXiv:2108.09886 [hep-ph] .
  * Nagakura and Zaizen (2023) Hiroki Nagakura and Masamichi Zaizen, “Connecting small-scale to large-scale structures of fast neutrino-flavor conversion,” Phys. Rev. D 107, 063033 (2023), arXiv:2211.01398 [astro-ph.HE] .
  * Duan and Shalgar (2015) Huaiyu Duan and Shashank Shalgar, “Flavor instabilities in the neutrino line model,” Phys. Lett. B747, 139–143 (2015), arXiv:1412.7097 [hep-ph] .
  * Abbar  _et al._ (2015) Sajad Abbar, Huaiyu Duan, and Shashank Shalgar, “Flavor instabilities in the multiangle neutrino line model,” Phys. Rev. D92, 065019 (2015), arXiv:1507.08992 [hep-ph] .
  * Nagakura and Zaizen (2022) Hiroki Nagakura and Masamichi Zaizen, “Time-Dependent and Quasisteady Features of Fast Neutrino-Flavor Conversion,” Phys. Rev. Lett. 129, 261101 (2022), arXiv:2206.04097 [astro-ph.HE] .
  * O’Connor (2015) Evan O’Connor, “An Open-Source Neutrino Radiation Hydrodynamics Code for Core-Collapse Supernovae,” Astrophys. J. Suppl. 219, 24 (2015), arXiv:1411.7058 [astro-ph.HE] .
  * Itoh  _et al._ (1996) Naoki Itoh, Hiroshi Hayashi, Akinori Nishikawa, and Yasuharu Kohyama, “Neutrino Energy Loss in Stellar Interiors. VII. Pair, Photo-, Plasma, Bremsstrahlung, and Recombination Neutrino Processes,” Astrophys. J. Suppl. 102, 411 (1996).
  * Bollig  _et al._ (2017) Robert Bollig, Hans-Thomas Janka, Andreas Lohs, Gabriel Martinez-Pinedo, Charles J. Horowitz, and Tobias Melson, “Muon Creation in Supernova Matter Facilitates Neutrino-driven Explosions,” Phys. Rev. Lett. 119, 242702 (2017), arXiv:1706.04630 [astro-ph.HE] .
  * Thompson (2002) Todd A. Thompson,  _Topics in the theory of core-collapse supernovae_ , Ph.D. thesis, The University of Arizona. (2002).
  * Bruenn (1985) Stephen W. Bruenn, “Stellar core collapse: Numerical model and infall epoch,” Astrophys. J. Suppl. 58, 771–841 (1985).
  * Vogel (1984) Petr Vogel, “Analysis of the Anti-neutrino Capture on Protons,” Phys. Rev. D 29, 1918 (1984).
  * (59) Georg G. Raffelt,  _Stars as laboratories for fundamental physics: The astrophysics of neutrinos, axions, and other weakly interacting particles_ (The University of Chicago Press).
  * Betranhandy (2022) Aurore Betranhandy,  _Neutrino interactions and axion emission impact on core-collapse supernova simulations_ , Ph.D. thesis, Stockholm U. (2022).
  * (61) Evan O’Connor, private communication.