# Fast neutrino flavor conversions in a supernova: emergence, evolution, and
effects

Zewei Xiong  0000-0002-2385-6771<EMAIL_ADDRESS>GSI Helmholtzzentrum für
Schwerionenforschung, Planckstraße 1, 64291 Darmstadt, Germany  Meng-Ru Wu
0000-0003-4960-8706 Institute of Physics, Academia Sinica, Taipei 11529,
Taiwan  Institute of Astronomy and Astrophysics, Academia Sinica, Taipei
10617, Taiwan  Physics Division, National Center for Theoretical Sciences,
Taipei 10617, Taiwan  Manu George  0000-0002-8974-5986 Institute of Physics,
Academia Sinica, Taipei 11529, Taiwan  Chun-Yu Lin  0000-0002-7489-7418
National Center for High-performance Computing, Hsinchu 30076, Taiwan  Noshad
Khosravi Largani  0000-0003-1551-0508 Institute of Theoretical Physics,
University of Wroclaw, Pl. M. Borna 9, 50-204 Wroclaw, Poland  Tobias Fischer
0000-0003-2479-344X Institute of Theoretical Physics, University of Wroclaw,
Pl. M. Borna 9, 50-204 Wroclaw, Poland  Gabriel Martínez-Pinedo
0000-0002-3825-0131 GSI Helmholtzzentrum für Schwerionenforschung,
Planckstraße 1, 64291 Darmstadt, Germany  Institut für Kernphysik
(Theoriezentrum), Fachbereich Physik, Technische Universität Darmstadt,
Schlossgartenstraße 2, 64289 Darmstadt, Germany

###### Abstract

Fast flavor conversions (FFCs) of neutrinos, which can occur in core-collapse
supernovae (CCSNe), are multiangle effects. They depend on the angular
distribution of the neutrino’s electron lepton number (ELN). In this work, we
present a comprehensive study of the FFCs by solving the multienergy and
multiangle quantum kinetic equations with an extended set of collisional weak
processes based on a static and spherically symmetric CCSN matter background
profile. We investigate the emergence and evolution of FFCs in models
featuring different ELN angular distributions, considering scenarios with two
and three neutrino flavors. The spectrogram method is utilized to illustrate
the small-scale spatial structure, and we show that this structure of neutrino
flavor coherence and number densities in the nonlinear regime is qualitatively
consistent with the dispersion relation analysis. On the coarse-grained level,
we find that different asymptotic states can be achieved following the FFCs
depending on the locations and shapes of the ELN distributions, despite
sharing a common feature of the elimination of the ELN angular crossing. While
equilibration among different neutrino flavors may be achieved immediately
after the prompt FFCs, it is not a general outcome of the asymptotic state, as
subsequent feedback effects from collisional neutrino-matter interactions come
into play, particularly for cases where FFCs occur inside the neutrinosphere.
The impacts of FFCs and the feedback effect on the net neutrino heating rates,
the equilibrium electron fraction of CCSN matter, and the free-streaming
neutrino energy spectra are quantitatively assessed. Other aspects including
the impact of the vacuum term and the coexistence with other type of flavor
instabilities are also discussed.

##  I Introduction

Core-collapse supernovae (CCSNe) are cataclysmic events when massive stars
reach the end of their lives. During the evolution of a CCSN, a nascent proto-
neutron star (PNS) forms at the center as a prolific source of neutrinos.
These neutrinos play pivotal roles on the CCSN dynamics and the evolution of
chemical composition. They interact with the medium through both charged- and
neutral-current weak interactions in the proximity of the PNS and deposit
energy, facilitating the shock revival via reheating of material in the post-
shocked layer and leading to the eventual mass ejection. Particularly, the
charged-current interactions determine the proton-to-baryon ratio, denoted by
the electron fraction $Y_{e}$, which is a crucial quantity for the
nucleosynthesis results of CCSN explosions. A better knowledge of the flux
intensities and flavor content of neutrinos is necessary to robustly model the
inner dynamics of CCSNe, to predict the elemental compositions in CCSN ejecta,
and to determine the neutrino signals for the detection of the next galactic
supernova event (see e.g., Refs. [1, 2] for recent reviews).

The flavor oscillations of neutrinos among $\nu_{e}$, $\nu_{\mu}$, and
$\nu_{\tau}$ in vacuum and in medium have been well studied and confirmed by
various ground-based neutrino experiments [3]. When neutrino fluxes are
sufficiently high in CCSNe, the forward scattering among neutrinos themselves
leads to various collective phenomena of flavor instability. Particularly the
fast flavor conversion (FFC) associated with the fast flavor instability (FFI;
see e.g., Refs. [4, 5, 6, 7] for reviews) attracts great interest in recent
years owing to the vastly rapid conversion rate within nanoseconds and over a
distance shorter than a coin. Studies based on results from the multi-
dimensional CCSN simulations have shown a general occurrence of FFIs in
certain regions ahead of the shock wave, near the neutrinosphere, or even deep
inside the PNS [8, 9, 10, 11, 12, 13, 14, 15, 16, 17].

Since the first proposal of FFI in Ref. [18], the advance of both theory and
methodology throughout the past decade has improved our understanding of FFCs
and guided the line of research toward the ultimate goal of implementing
neutrino oscillations in CCSNe. A generic framework governing the coherent
flavor evolution and collisional neutrino-matter interactions is prescribed by
the neutrino quantum kinetic equation ($\nu$QKE) [19, 20, 21, 22]. The linear
stability analysis (LSA) of the $\nu$QKE provides a powerful tool to diagnose
the existence of flavor instabilities [23, 24, 25, 26]. It has been proved
based on the LSA that the FFI requires the presence of angular zero-crossings
where the angular distribution of the neutrino lepton number between any two
distinct flavors transitions from positive to negative values [27]. The
dispersion relation based on LSA implies the break of spatial homogeneity,
which results in small-scale spatial structure confirmed by several dynamical
simulations that numerically solved $\nu$QKE within tiny boxes [28, 29, 30,
31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42]. Quasi-steady states are found
to be achieved asymptotically on the coarse-grained level after neutrinos
undergo the kinematic decoherence [43, 44, 45, 29, 46]. It was reported by
these simulations that the FFCs lead to complete flavor equilibration or
partial flavor conversion depending on the ratio of neutrino fluxes and the
characteristics of the angular zero-crossings. Particularly when a complete
flavor equilibration or even overconversion is achieved, the neutrino flavor
content can be significantly affected, which was demonstrated to be able to
exert substantial impact on the neutrino-driven explosion mechanism and the
nucleosynthesis of CCSNe [47, 48, 49, 50, 51, 52], or even in the context of
neutron-star mergers [53, 54, 55, 56, 57, 58, 59, 60].

Recently, a growing amount of studies solving $\nu$QKE in a large spatial
domain ranging from the interior of the PNS to the region around the shock
wave were performed [61, 62, 63, 64, 65, 66]. These studies intended to
capture the flavor instabilities and to evolve both collective flavor
conversions as well as the neutrino transport simultaneously in a self-
consistent manner. Despite all these triumphs in the progress of FFCs, there
remain several unaddressed questions. (1) Does the small-scale structure
observed in the local simulations survive when the advection and neutrino-
matter interactions are considered? (2) What is the general outcome of FFCs as
well as its impact on the neutrino flux intensities and their ratios? Is
flavor equilibration the ultimate fate of FFCs? (3) What consequence could be
anticipated on CCSNe and the associated neutrino signals?

In addition to the FFI, there are other flavor instabilities such as the slow
(see Refs. [67, 1] for reviews), or collisional type (see Refs. [68, 69, 70,
71, 72, 73, 74, 75, 76, 77, 78, 79, 80]) that might coexist with the FFI.
Particularly, the collisional flavor instability (CFI), which was discovered
recently, can potentially trigger flavor conversions in regions close to or
inside the neutrinosphere where the FFI is known to play a role. Certainly,
there are still open issues regarding these topics. It is important to
understand whether and how different flavor instabilities interplay with each
other in more sophisticated supernova models.

In this work, we aim to provide answers to the above questions (1)–(3), and to
investigate the coexistence as well as impact of other flavor instabilities.
We extend our previous multi-energy and multi-angle neutrino flavor evolution
code [81, 71] to solve the $\nu$QKE by including a more complete set of
collisional weak interactions in static snapshots of background environments
obtained with spherically symmetric CCSN hydro simulations. Notice that the
multi-energy spherically symmetric neutrino transport in standard CCSNe
typically finds no angular crossing [82, 83, 77] because the imposed spherical
symmetry restricts the capability of capturing the multi-dimensional CCSN
features such as the fluid convection [83] and lepton-emission self-sustained
asymmetry (LESA) [84, 85]. These multi-dimensional features can lead to
attenuated $Y_{e}$ compared to the spherically symmetric model and create
conditions for FFIs. Therefore, we follow the spirit of Ref. [49] and
introduce radial-dependent $Y_{e}$ attenuation schemes in order to mimic the
suitable conditions for FFCs commonly observed in multi-dimensional
simulations. We use schemes with different attenuation strengths to probe a
variety of scenarios for the FFI close to and inside the neutrinosphere.

We adopt a similar two-step approach as in Ref. [71]. We first relax neutrino
profiles into stationary states up to $\simeq 1$ ms without including coherent
flavor oscillations. Those stationary states are then used as the initial
condition to simulate the evolution of neutrinos including diffusion,
collisions, and flavor oscillations. We then present a comprehensive analysis
for the FFCs covering three main aspects: the small-scale structure after the
FFI emerges from the linear to nonlinear regimes, the dynamical evolution of
the neutrino flavor content, and the effects on CCSN physics as well as the
neutrino signals.

This paper is organized as follows. We provide the descriptions and parameters
of our models for the FFCs in Sec. II. We use the LSA and the spectrogram
method introduced in Sec. III to characterize the small-scale features of FFI
in Sec. IV. We present a detailed study on the evolution of FFCs in Sec. V and
discuss their implications in the energy transport, chemical composition of
the matter, and free-streaming neutrino spectra in Sec. VI. Further
discussions and conclusions are given in Sec. VII. We adopt natural units with
$\hbar=c=k_{B}=1$ throughout the paper.

##  II Models

In this section we introduce our theoretical model for neutrino flavor
conversions in CCSNe. We begin with a review about the supernova model and the
simulations based on which the flavor study is conducted. We then revisit the
neutrino quantum kinetic approach and introduce our computational
implementation.

###  II.1 Supernova model

The CCSN simulation considered in this work is launched from the stellar
progenitor star with zero-age main sequence mass of 25 $M_{\odot}$ and solar
metalicity, from the stellar evolution series of Ref. [86], denoted as s25a28.
Our agile-boltztran supernova model is based on general relativistic spherical
symmetric neutrino radiation hydrodynamics in comoving coordinates [87, 88,
89, 90, 91], including six-species Boltzmann neutrino transport assuming
ultra-relativistic particles [92, 93]. The equations of radiation
hydrodynamics are solved on an adaptive baryon mass mesh [94, 95], for which
102 radial grid points are used for the current simulations, including about 5
$M_{\odot}$ of the stellar core, which covers parts of the extended silicon-
sulfur layer for this progenitor model. The neutrino distributions are solved
by means of the Boltzmann neutrino transport equation in the discrete-ordinate
method. They are discretized in terms of 6 propagation angles and 36 energy
bins (from 0.5–300 MeV) following the setup of Refs. [89, 96].

For the collision integral of the Boltzmann equation, a complete set of weak
processes is employed (see Table I in Ref. [97]). This includes the charged-
current neutrino emission and absorption (EA) involving heavy nuclei [98] and
unbound nucleons [97], isoenergetic neutrino scattering on nuclei and nucleons
(IS) [96, 89], inelastic neutrino scattering on electrons and positrons (NES)
[99, 88], and neutrino pair reactions (PR), including electron-positron
annihilation [96], nucleon-nucleon bremsstrahlung [100, 101] and the
annihilation of electron neutrino pairs [102, 103]. EA are implemented in the
full kinematics approach at the mean-field level, taking weak-magnetism
contributions self-consistently into account [104]. Note that muonic weak
processes are omitted here, such that muon and tauon flavors are not
distinguished for neutrinos and antineutrinos. For the equation of state (EOS)
used in the present supernova simulation, see details in Appendix A.

During the early post-bounce evolution, on the order of 10–25 ms after the
core bounce, the $\nu_{e}$ deleptonization burst is released. It is a standard
feature of CCSN phenomenology that results in the rapid drop of the electron
fraction $Y_{e}$, accordingly, from around $Y_{e}=0.3$–$0.4$ at core bounce to
$Y_{e}\simeq 0.1$, in the region where the neutrinos decouple from matter.
Even though neutrino decoupling is an energy and momentum dependent process
(c.f. Refs. [105, 106]), often the neutrino-phase-space-averaged
neutrinospheres of last elastic and inelastic scatterings are used as
references [107], inside of which the bulk part of the neutrino spectrum is
trapped and outside of which it is freely streaming. The post-bounce evolution
of this s25a28 progenitor model corresponds in spherically symmetric supernova
simulations to the failed branch, leading eventually to the PNS collapse and
the formation of a black hole [108, 103, 109]. Before that happens, the post-
bounce phase is determined by a steady state mass accretion from the
gravitationally unstable layers above the stellar core. The infalling material
falls onto the standing bounce shock. The latter is located at a radius of
about 80–100 km, depending on the progenitor model and the nuclear EOS.
Thereby, the infalling heavy nuclei are dissociated into bulk nuclear matter,
composed of unbound nucleons and light nuclei, that accumulate at the PNS
surface. This layer, illustrated in Fig. 1 via the shock at $r\simeq 85$ km
and the PNS surface at $r\simeq 40$ km depicting a shallow rise in density
(solid curve) and temperature (dotted curve), is subject to the neutrino
decoupling, spanning a range of temperatures from up to $T\simeq 10$ MeV, at
densities of a few times $10^{13}$ g cm-3, down to $T\simeq 3$ MeV which
corresponds to the shock location with densities around $10^{9}$ g cm-3 where
the temperature increases sharply due to the shock heating from $T\simeq 0.5$
MeV to $T\simeq 3$ MeV.

Figure 1:  Radial profiles of baryonic density, temperature $T/10$, and
$Y_{e}$ (before and after attenuations) for a snapshot at the post-bounce time
$t_{\rm pb}\approx 252$ ms. Orange, blue, and green vertical dashed lines
indicate the gain radius, the radii of $\nu_{e}$\- and $\bar{\nu}_{e}$-spheres
at $r\approx 63$ km, 40 km, and 35 km, respectively.

The timescale for the compression of the accumulated layer, settling at the
PNS surface into the steepening gravitational potential, is determined by the
high-density nuclear EOS, its temperature and isospin asymmetry dependence,
and the enclosed PNS mass. As a consequence, neutrino luminosities $L_{\nu}$
are obtained on the order of several $10^{52}$ erg s-1, after the $\nu_{e}$
deleptonization burst has been launched, and the average neutrino energies
have values of $\langle E_{\nu}\rangle\simeq 10$–20 MeV. For a recent review
of the role of neutrinos in the post-bounce phase of core-collapse supernovae,
see Refs. [110, 2] and references therein.

Instead of surveying over different time stages of supernova, we will focus on
one snapshot at a post-bounce time $t_{\rm pb}\approx 252$ ms in Fig. 1. Based
on this snapshot, we explore different scenarios for FFCs by applying a
radial-dependent attenuation to the original $Y_{e}$ profile to mimic suitable
conditions for FFCs observed in multi-dimensional supernova simulations. A
radial-dependent factor $b_{Y_{e}}=1-b_{1}/[1+\exp(r-b_{2}[\mathrm{km}])]$ is
multiplied to the original $Y_{e}$ profile with $b_{1}$ denoting the fraction
of attenuation and $b_{2}$ characterizing the radial range of attenuation.
With the attenuated electron fraction $Y_{e}$ and smaller electron chemical
potential $\mu_{e^{-}}$, the equilibrium neutrino chemical potential
$\mu_{\nu_{e}}^{\rm eq}=\mu_{e^{-}}+\mu_{p}-\mu_{n}$ is reduced to create a
condition less dominated by $\nu_{e}$’s, or even a condition where
$\bar{\nu}_{e}$’s dominate over $\nu_{e}$’s. We vary $b_{1}$ from 10% to 25%,
allowing different angular crossings dominated by $\nu_{e}$ or $\bar{\nu}_{e}$
to appear, which will be discussed later in Sec III.1. We show in Fig. 1 three
attenuated $Y_{e}$ profiles with different values of $(b_{1},b_{2}[{\rm
km}])=(0.1,42),(0.25,42)$, and $(0.25,30)$, respectively. Note that taking
$b_{2}=30$ km restricts the FFI solely within the PNS, but not extending to
the free-streaming regime. Also noted is that we only attenuate $Y_{e}$
without modifying the neutron and proton abundances for simplification. When
the latter abundances are modified accordingly based on the EOS, $\mu_{p}$
decreases, and $\mu_{n}$ increases. As a result, the attenuation factor,
$b_{1}$, required to achieve a similar scenario for each model is reduced by
approximately half. In addition, although the $Y_{e}$ profile with the largest
attenuation factor adopted here might be hard to be obtained in typical
supernovae, it may mimic other physical conditions such as the second collapse
in a hadron-quark phase transition supernova featuring a unique burst-like
neutrino signature dominated by $\bar{\nu}_{e}$ [111, 112, 113, 114, 115].

###  II.2 Neutrino quantum kinetic equation

We solve the neutrino quantum kinetic equation ($\nu$QKE) in both two- and
three-flavor schemes. The general equations governing the spatial and temporal
evolution of density matrices $\varrho$ (for neutrinos) and $\bar{\varrho}$
(for anti-neutrinos) are given by

$\displaystyle(\partial_{t}+v_{r}\partial_{r}+\frac{1-v_{r}^{2}}{r}\partial_{v_{r}})\varrho$
$\displaystyle=-i[\mathbf{H},\varrho]+\mathbf{C},$ (1)

and

$\displaystyle(\partial_{t}+v_{r}\partial_{r}+\frac{1-v_{r}^{2}}{r}\partial_{v_{r}})\bar{\varrho}$
$\displaystyle=-i[\bar{\mathbf{H}},\bar{\varrho}]+\bar{\mathbf{C}},$ (2)

with

$\varrho=\begin{bmatrix}\varrho_{ee}&\varrho_{e\mu}&\varrho_{e\tau}\\\
\varrho_{e\mu}^{*}&\varrho_{\mu\mu}&\varrho_{\mu\tau}\\\
\varrho_{e\tau}^{*}&\varrho_{\mu\tau}^{*}&\varrho_{\tau\tau}\end{bmatrix}\text{~{}and~{}}\bar{\varrho}=\begin{bmatrix}\bar{\varrho}_{ee}&\bar{\varrho}_{e\mu}&\bar{\varrho}_{e\tau}\\\
\bar{\varrho}_{e\mu}^{*}&\bar{\varrho}_{\mu\mu}&\bar{\varrho}_{\mu\tau}\\\
\bar{\varrho}_{e\tau}^{*}&\bar{\varrho}_{\mu\tau}^{*}&\bar{\varrho}_{\tau\tau}\end{bmatrix}$
(3)

in three flavor bases or

$\varrho=\begin{bmatrix}\varrho_{ee}&\varrho_{e\mu}\\\
\varrho_{e\mu}^{*}&\varrho_{\mu\mu}\end{bmatrix}\text{~{}and~{}}\bar{\varrho}=\begin{bmatrix}\bar{\varrho}_{ee}&\bar{\varrho}_{e\mu}\\\
\bar{\varrho}_{e\mu}^{*}&\bar{\varrho}_{\mu\mu}\end{bmatrix},$ (4)

in a simplified neutrino system of two flavors. Their diagonal elements are
normalized to the neutrino number density $n_{\nu_{i}}(r,t)=\int
dE\,dv_{r}\,\varrho_{ii}$ where $i=e,\,\mu,\,\tau$.

In general, the coherent propagation Hamiltonian $\mathbf{H}$ on the right-
hand side has three contributions. They are the vacuum mixing term

$\mathbf{H}_{\rm vac}(E)=\mathbf{U}\mathbf{M}\mathbf{U}^{\dagger},$ (5)

with the diagonal matrix in mass basis $\mathbf{M}=\mathrm{diag}[0,\delta
m_{\rm L}^{2},\delta m^{2}]/(2E)$ and the unitary mixing matrix $\mathbf{U}$
parametrized as in Eq. (14.34) of Ref. [3], the matter term

$\mathbf{H}_{\mathrm{mat}}=\mathrm{diag}[V_{\rm mat},0,0],$ (6)

with the effective potential $V_{\rm mat}$ corresponding to neutrino forward
scattering on ordinary matter, and the neutrino self-induced term

$\mathbf{H}_{\nu\nu}(v_{r})=\sqrt{2}G_{F}\int
dE^{\prime}\,dv_{r}^{\prime}(1-v_{r}v_{r}^{\prime})(\varrho-\bar{\varrho}^{*}),$
(7)

corresponding to neutrino forward scattering on other neutrinos, respectively.
In the two-flavor scheme, the vacuum term is explicitly

$\mathbf{H}_{\rm vac}(E)=\frac{\delta m^{2}}{4E}\begin{bmatrix}-\cos
2\theta_{V}&\sin 2\theta_{V}\\\ \sin 2\theta_{V}&\cos
2\theta_{V}\end{bmatrix}$ (8)

with one mass-square difference $\delta m^{2}$ and one mixing angle
$\theta_{V}$. Although the matter term may affect the behaviors of flavor
oscillations [116, 71], we will not include it explicitly in the following
discussion of this work.

The collisional term $\mathbf{C}$ includes all four types of weak processes in
Table 1. We use the same scheme for the EA and IS processes as in Ref. [71]
except that muonic EA reactions are not considered.

Table 1:  Set of weak processes considered in cose$\nu$, where $\nu$ and $\bar{\nu}$ are for all neutrino flavors and $N=n,p$. Label | Weak process | Abbreviation
---|---|---
(1a) | $\nu_{e}+n\leftrightarrows p+e^{-}$ | EA
(1b) | $\bar{\nu}_{e}+p\leftrightarrows n+e^{+}$ | EA
(1c) | $\bar{\nu}_{e}+p+e^{-}\leftrightarrows n$ | EA
(2a) | $\nu+N\leftrightarrows\nu+N$ | IS
(2b) | $\bar{\nu}+N\leftrightarrows\bar{\nu}+N$ | IS
(3a) | $\nu+e^{\pm}\leftrightarrows\nu+e^{\pm}$ | NES
(3b) | $\bar{\nu}+e^{\pm}\leftrightarrows\bar{\nu}+e^{\pm}$ | NES
(4a) | $\nu+\bar{\nu}\leftrightarrows e^{-}+e^{+}$ | PR
(4b) | $\nu+\bar{\nu}+N+N\leftrightarrows N+N$ | PR

Low-energy electron neutrinos and heavy-lepton neutrinos are thermalized by
NES and PR. Because of the high computational cost of the scattering kernels
associated with them, approximated schemes are employed. For PR, we integrate
the scattering kernels in Refs. [96, 117] into the emissivity $j_{\rm PR}$ and
opacity $\chi_{\rm PR}$ for all flavors in a similar way as in Ref. [118]
except assuming that the accompanied neutrinos or antineutrinos in the pair
processes obey the Fermi-Dirac distribution determined by the equilibrium
neutrino chemical potential (e.g., $\mu_{\nu_{e}}^{\rm eq}$ for $\nu_{e}$), in
order to take the blocking effect into account. For heavy-lepton neutrinos,
the equilibrium neutrino chemical potential is assumed to be zero. The
collisional term for neutrino can then be constructed in an effective way
similar to the EA processes as follows,

$\displaystyle\mathbf{C}_{\mathrm{PR}}(E)=$
$\displaystyle\frac{1}{2}\left\\{\mathrm{diag}[j_{e,\rm PR},j_{\mu,\rm
PR},j_{\tau,\rm PR}],\varrho_{\rm FO}-\varrho\right\\}$
$\displaystyle-\frac{1}{2}\left\\{\mathrm{diag}[\chi_{e,\rm PR},\chi_{\mu,\rm
PR},\chi_{\tau,\rm PR}],\varrho\right\\}$ (9)

with the fully-occupied differential number density $\varrho_{\rm FO}$ and the
curly bracket being the anti-commutator. The full term of NES in $\nu$QKE
involves channels of different Mandelstam variables [22], which increases the
complication of the treatment. For simplicity, we take the classical treatment
particularly for NES by only keeping the diagonal elements of $\mathbf{C}_{\rm
NES}$ [119, 120, 99, 88]. The thermalization of the energy spectra of heavy-
lepton neutrinos is guaranteed by the Kirchhoff-Planck relation (“detailed
balance”) between the in- and out-scattering kernels

$R_{\rm NES}^{\rm
in}(E,E^{\prime},v_{r},v_{r}^{\prime})=e^{\frac{E^{\prime}-E}{T}}R_{\rm
NES}^{\rm out}(E,E^{\prime},v_{r},v_{r}^{\prime}).$ (10)

In addition, for each flavor the NES scattering kernel is truncated at the
zeroth moment following the method in [99, 121] (See Appendix B for more
details).

###  II.3 Numerical setup

The $\nu$QKE is solved by cose$\nu$ under spherical symmetry within a radial
range uniformly discretized by $N_{r}$ radial grids between $r_{\rm ib}=20$ km
and $r_{\rm ob}=80$ km. At the inner boundary neutrinos of all species are in
Fermi distribution functions with chemical potentials given by weak
equilibrium with the processes of NES and PR included. The outer boundary is
chosen to be around the position of the shock wave. We employ the free-
streaming boundary condition there for forward propagating neutrinos with
$v_{r}\geq 0$. Although neutrinos can be back-scattered by heavy nuclei
outside the shock through the neutrino-nucleus interactions, which leads to
the presence of FFI around the shock [12], a recent work suggests that their
effect on the dominant neutrino outflow is negligible due to the relatively
diluted flux [35]. Therefore, we inject no neutrinos in the backward direction
with $v_{r}<0$.

Because the vacuum term $\mathbf{H}_{\rm vac}$ with a non-zero mixing angle
generates flavor mixing automatically, we do not employ any artificial
perturbation as the seed to trigger flavor conversions in the initial
condition.

For consistency across all models, the neutrino spectral and angular
distributions are all discretized into $N_{E}=15$ energy grids and
$N_{v_{r}}=100$ angular grids for both $\varrho$ and $\bar{\varrho}$. The
energy grids are taken between 1 MeV to 100 MeV spaced uniformly in
logarithmic scale, and the angular grids are linearly uniform between
$v_{r}=-1$ and 1.

Because both EA and IS rates that generally dominate over others strongly
depend on the neutrino energy and nucleonic density, taking a smaller maximum
neutrino energy and a larger radius of the inner boundary compared to our
previous work [71] reduces the maximum collisional rate that can be well
resolved within our simulation domain. Therefore, unlike in [71], no
attenuation is made to any collisional rate for all models in this work.

In contrast to the collisional rates, a common challenge is the tremendous
size of the neutrino self-induced term $\mathbf{H}_{\nu\nu}$ associated with
subnanosecond timescales and subcentimeter length scales. They are much
shorter than the scales related to the hydrodynamics (order of kilometers and
milliseconds), which makes it computationally difficult for a simulation to
well resolve both time and length scales with the unadjusted
$\mathbf{H}_{\nu\nu}$. A strategical mitigation in recent works [63, 71] is to
attenuate the values of $\mathbf{H}_{\nu\nu}$ while generally keeping the same
ordering that the magnitude of $\mathbf{H}_{\nu\nu}$ dominates over
collisional rates and vacuum oscillation frequency.111There are also other
strategies, for example, assuming that neutrinos are homogeneous in radial
blocks of a few hundred meters to reduce the needed number of the radial grids
in computation and to allow for no attenuation on $\mathbf{H}_{\nu\nu}$ [64].
The oscillation term $\mathbf{H}$ on the right-hand side of Eq. (1) becomes
$\mathbf{H}_{\rm vac}+a_{\nu\nu}(r)\mathbf{H}_{\nu\nu}$ with a radial-
dependent factor $a_{\nu\nu}(r)=a_{1}/[1+e^{(a_{2}-r)/a_{3}}]$ where
$a_{2}=30~{}{\rm km}$ and $a_{3}=2.5~{}{\rm km}$ for all our models. This
factor keeps almost constant $a_{\nu\nu}(r)\approx a_{1}$ for most of the
radial range with $r>35$ km. At inner radii where the neutrino number
densities are much larger, more attenuation is provided by the denominator in
the formula by around one more order of magnitude, i.e., the attenuation
factor $a_{\nu\nu}$ takes a value of $0.1a_{1}$ at $r\approx 24.5$ km. For
example, when $a_{1}$ is taken to be $10^{-3}$, the attenuation factor is
always greater than $10^{-4}$ for $r>24.5$ km and close to $10^{-3}$ for
$r>35$ km, which is a reasonable choice compared to other works [58, 66]. We
utilize this attenuation strategy and explore different choices of the
attenuation parameter $a_{1}$, with the radial resolutions $N_{r}$ adapted
accordingly.

We summarize all models and their parameters of this paper in Table 2. There
is no flavor instability in Model I. In both two- and three-flavor models, all
mixing angles in $\mathbf{H}_{\rm vac}$ are set to be $10^{-6}$ accounting for
the effective suppression by the matter. In three-flavor models, the CP
violation phase is omitted, and the lower mass-square difference is taken as
$\delta m_{\rm L}^{2}=7.39\times 10^{-5}~{}\mathrm{eV}^{2}$. In order to
investigate the effects of other types of flavor instabilities on the FFCs, we
run simulations with different values of $\delta m^{2}$ in the vacuum term and
introduce an auxiliary parameter $c_{\rm off}$ that are multiplied to all off-
diagonal elements of $\mathbf{C}$. The default value of $c_{\rm off}$ is 1.
When $c_{\rm off}$ is set to be zero, the CFI is prohibited.

Table 2:  Parameters used in each model. For all models, $r_{\rm ib}=20$ km to $r_{\rm ob}=80$ km, $t_{\mathrm{pb}}$=252 ms, $a_{2}=30$ km, $a_{3}=2.5$ km, $N_{v_{r}}=100$, $N_{E}=15$. model | $b_{1}$ |  $b_{2}$ [km] | no. of flavors | $a_{1}$ | $N_{r}$ | $c_{\rm off}$ |  $\delta m^{2}$ [eV2]
---|---|---|---|---|---|---|---
I | 0 |  |  |  |  |  |
II | 0.1 | 42 | 2 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
III | 0.15 | 42 | 2 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
IV | 0.2 | 42 | 2 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
V | 0.25 | 42 | 2 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
VI | 0.25 | 32 | 2 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
VII | 0.25 | 30 | 2 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
Vc0 | 0.25 | 42 | 2 | $10^{-3}$ | 25000 | 0 | $8\times 10^{-5}$
IIf3 | 0.1 | 42 | 3 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
IVf3 | 0.2 | 42 | 3 | $10^{-3}$ | 25000 | 1 | $8\times 10^{-5}$
IIa2 | 0.1 | 42 | 2 | $2\times 10^{-3}$ | 50000 | 1 | $8\times 10^{-5}$
IIa4 | 0.1 | 42 | 2 | $4\times 10^{-3}$ | 50000 | 1 | $8\times 10^{-5}$
IIa10 | 0.1 | 42 | 2 | $10^{-2}$ | 50000 | 1 | $8\times 10^{-5}$
IIIa2 | 0.15 | 42 | 2 | $2\times 10^{-3}$ | 50000 | 1 | $8\times 10^{-5}$
IIIa4 | 0.15 | 42 | 2 | $4\times 10^{-3}$ | 50000 | 1 | $8\times 10^{-5}$
IVa2 | 0.2 | 42 | 2 | $2\times 10^{-3}$ | 50000 | 1 | $8\times 10^{-5}$
IVa4 | 0.2 | 42 | 2 | $4\times 10^{-3}$ | 50000 | 1 | $8\times 10^{-5}$
IVv2 | 0.2 | 42 | 2 | $10^{-3}$ | 25000 | 1 | $2.4\times 10^{-3}$

##  III Analyses

###  III.1 Illustrative quantities

For a better illustration of our models, we define several quantities as
follows. The energy- and angle-integrated density matrices are defined by
$\langle\varrho\rangle_{E}(v_{r})=\int dE\,\varrho(E,v_{r})$ and
$\langle\varrho\rangle_{A}(E)=\int dv_{r}\,\varrho(E,v_{r})$, respectively.
The neutrino mean energy of a given flavor is calculated as

$\langle E_{\nu_{i}}\rangle=\frac{\int
dE\,dv_{r}\,E\varrho_{ii}(E,v_{r})}{\int dE\,dv_{r}\,\varrho_{ii}(E,v_{r})},$
(11)

where $i=e,\,\mu,\,\tau$. A dimensionless ratio,
$s_{e\mu}=|\langle\varrho_{e\mu}\rangle_{E}|/(|\langle\varrho_{ee}\rangle_{E}-\langle\varrho_{\mu\mu}\rangle_{E}|^{2}/4+|\langle\varrho_{e\mu}\rangle_{E}|^{2})^{1/2}$,
is defined to quantify the nonlinearity of flavor conversions caused by FFI.
The corresponding phase (angle in the complex plane) is defined as
$\phi_{e\mu}=\arg[\langle\varrho_{e\mu}\rangle_{E}\,{\rm
sgn}(\langle\varrho_{ee}\rangle_{E}-\langle\varrho_{\mu\mu}\rangle_{E})]$,
where $\mathrm{sgn}$ denotes the sign function. The electron-minus-muon lepton
number of neutrinos (E-MuLN, or shortened as ELN in the following discussion)
is
$G=\langle\varrho_{ee}\rangle_{E}-\langle\bar{\varrho}_{ee}\rangle_{E}-\langle\varrho_{\mu\mu}\rangle_{E}+\langle\bar{\varrho}_{\mu\mu}\rangle_{E}$.
When there only exists one zero-crossing in the ELN angular distribution, we
define two integrals $I_{+}=\int_{-1}^{1}dv_{r}\,G\Theta(G)$ and
$I_{-}=\int_{-1}^{1}dv_{r}\,G\Theta(-G)$, where $\Theta$ is the Heaviside step
function. Those two quantities allow us to evaluate the “deepness ratio” of
the ELN angular crossing based on the value of ${\rm
min}(|I_{-}|,|I_{+}|)/{\rm max}(|I_{-}|,|I_{+}|)$.

Model II
Model III
Model IV
Model V
Model VI
Model VII

Figure 2:  Radial profiles of ELN distribution in Models II–VII that have
different attenuated $Y_{e}$ profiles. The angular crossings of ELN are
indicated by the location of curves with the white color. Model I (without
$Y_{e}$ attenuation) is not shown as it contains no ELN angular crossing.
Figure 3:  Deepness ratios of the initial ELN distribution, ${\rm
min}(|I_{-}|,|I_{+}|)/{\rm max}(|I_{-}|,|I_{+}|)$, in Models II–V (see text
for the definition of $I_{\pm}$).

Except the original Model I without $Y_{e}$ attenuation and angular crossings,
Fig. 2 shows the initial ELN distributions for Models II–VII in the stationary
state obtained by solving Eqs. (1) and (2) without including the commutators
on the right-hand side. Clearly, these six models show a significant variation
in the radial shapes of their ELN angular crossings due to different $Y_{e}$
attenuation, while each of them only contains at most one angular crossing at
any specific radius. In particular, the ELN crossings in Models II–V extend
all the way to the outer boundary where neutrinos freely stream outwards. With
the least $Y_{e}$ attenuation imposed in Model II (corresponding to the red
curve in Fig. 1), a relatively shallow angular crossing begins to appear with
crossing radial velocity $v_{r,c}\lesssim 1$ at $r\approx 38$ km, and extends
to smaller $v_{r,c}$ as low as $v_{r,c}\approx 0.8$ at $r\approx 41$ km. For
the entire domain, the deepness ratios are less than 0.06 as shown in Fig. 3.
When we apply more attenuation on $Y_{e}$, which results in more production of
$\bar{\nu}_{e}$’s, the starting radii where the angular crossing appears shift
to $r\approx 33$ km and $28$ km in Model III and IV, respectively. The minimum
values of $v_{r,c}$ are $v_{r,c}\approx 0.3$ at $r\approx 38$ km and
$v_{r,c}\approx-0.4$ at $r\approx 32$ km accordingly. For Model V
(corresponding to the cyan curve in Fig. 1), the inner boundary becomes so
$\bar{\nu}_{e}$-dominated that the crossing starts from $v_{r,c}=-1$ at
$\approx 35$ km and monotonically approaches $v_{r,c}\approx 0.85$ as the
radius increases. The deepness ratios in Models III–V are generally greater
than 0.3, particularly reaching $\approx 1$ at $r\simeq 30$ km and 45 km in
Model IV. In the last two models (VI and VII), ELN crossings are only confined
inside the neutrino spheres as they originate from the change of sign of
$\mu^{\rm eq}_{\nu_{e}}$ due to the $Y_{e}$ attenuation. The crossing in Model
VI spans from $r\approx 40$ km to $r\approx 42$ km, while it is narrower in
Model VII (corresponding to the green curve in Fig. 1), ranging from $r\approx
38$ km to $r\approx 39$ km. These different radial shapes of ELN crossings
lead to very different FFI and the evolution of FFC, which will be discussed
in the following sections.

###  III.2 Linear stability analysis

We perform the LSA for models with two flavors. The analyses on three-flavor
models are expected to be similar when neutrinos of muon and tauon flavors are
identical, although a nontrivial interplay among flavors may happen when
difference in $\mu$ and $\tau$ neutrinos exist [122, 123]. Assuming
$|\varrho_{e\mu}|/|\varrho_{ee}-\varrho_{\mu\mu}|\ll 1$ and
$|\bar{\varrho}_{e\mu}|/|\bar{\varrho}_{ee}-\bar{\varrho}_{\mu\mu}|\ll 1$
before any instability leads to flavor conversion, $\nu$QKE can be linearized
to find unstable normal modes if flavor instability exists [23, 24].d We
follow our previous treatment in [71], neglecting the terms associated with
$\partial_{v_{r}}$ and IS but including the vacuum term. We further assume
that a collective mode of the perturbation
$\varrho_{e\mu}=Q(\Omega,K_{r},E,r,v_{r})e^{-i[\Omega t-K_{r}(r^{\prime}-r)]}$
and $\bar{\varrho}_{e\mu}^{*}=\bar{Q}(\Omega,K_{r},E,r,v_{r})e^{-i[\Omega
t-K_{r}(r^{\prime}-r)]}$ can develop locally near $r$, so the off-diagonal
parts of Eqs. (1) and (2) become222The dependence of $Q$ and $\bar{Q}$ on
$\Omega$, $K_{r}$, and $r$ is not shown explicitly.

$\displaystyle\left[\Omega-K_{r}v_{r}-\Phi(v_{r})+iC_{e\mu}(E)+\delta
m^{2}/2E\right]Q(E,v_{r})$ $\displaystyle=$
$\displaystyle-\sqrt{2}G_{F}[\varrho_{ee}(E,v_{r})-\varrho_{\mu\mu}(E,v_{r})]\int
dE^{\prime}\,dv_{r}^{\prime}\,\times$
$\displaystyle(1-v_{r}v_{r}^{\prime})[Q(E^{\prime},v_{r}^{\prime})-\bar{Q}(E^{\prime},v_{r}^{\prime})],$
(12)

and

$\displaystyle\left[\Omega-K_{r}v_{r}-\Phi(v_{r})+i\bar{C}_{e\mu}(E)-\delta
m^{2}/2E\right]\bar{Q}(E,v_{r})$ $\displaystyle=$
$\displaystyle-\sqrt{2}G_{F}[\bar{\varrho}_{ee}(E,v_{r})-\bar{\varrho}_{\mu\mu}(E,v_{r})]\int
dE^{\prime}\,dv_{r}^{\prime}\,\times$
$\displaystyle(1-v_{r}v_{r}^{\prime})[Q(E^{\prime},v_{r}^{\prime})-\bar{Q}(E^{\prime},v_{r}^{\prime})],$
(13)

where $\Phi(v_{r})=\sqrt{2}G_{F}\int
dE^{\prime}\,dv_{r}^{\prime}\,(1-v_{r}v_{r}^{\prime})[\varrho_{ee}(E^{\prime},v_{r}^{\prime})-\varrho_{\mu\mu}(E^{\prime},v_{r}^{\prime})-\bar{\varrho}_{ee}(E^{\prime},v_{r}^{\prime})+\bar{\varrho}_{\mu\mu}(E^{\prime},v_{r}^{\prime})]$,
$C_{e\mu}(E)$ and $\bar{C}_{e\mu}(E)$ are the damping rates attributed to EA
and the approximate PR terms. For example, $C_{e\mu}(E)=(j_{e,\rm
EA}+\chi_{e,\rm EA}+j_{e,\rm PR}+\chi_{e,\rm PR}+j_{\mu,\rm PR}+\chi_{\mu,\rm
PR})/2$. To account for the attenuation scheme used for $\mathbf{H}_{\nu\nu}$,
we modify Eqs. (III.2) and (III.2) by multiplying $G_{F}$ with the attenuation
factor $a_{\nu\nu}$, which affects the right-hand side and the term $\Phi$ on
the left-hand side.

For a given wave number $K_{r}$, the eigenvalues of $\Omega$ and their
corresponding eigenvectors can be numerically solved. The contribution of
$\Phi$ is often absorbed into $\Omega$ and $K_{r}$ in literature [24, 26], but
we do not follow this procedure in order to keep the original meaning of
$K_{r}$ as the wave number in the local frame. When $K_{r}=0$, it represents
the homogeneous mode. We highlight that we search for non-zero $K_{r}$ modes
in additional to the homogeneous one to avoid missing the most unstable mode.

We note that there may exist some other aspects that cannot be captured by the
linearized equations (III.2) and (III.2). For instance, the behaviors of
flavor instability may be affected by the angular advection term associated
with $\partial_{v_{r}}$ and the fact that neutrino number density is not
homogeneous. Although the pure FFI without collisions and with vanishing
vacuum oscillation frequency is free from the contamination of spurious modes
[124], which are associated with the use of an insufficient number of angular
grids in the discretization, those spurious modes may appear in multiangle
slow instabilities when the vacuum term is taken into account [125, 126, 127]
or in the CFI when collisions are included [128, 74]. Their growth rates
should approach zero as the adopted angular grid number increases toward the
continuous limit.

###  III.3 Spectrogram

In addition to the LSA, the onset of specific unstable modes as well as the
small-scale spatial structure can be simply analyzed by the Fourier transform
of relevant quantities over the entire simulation domain in the periodic-box
setup. In contrast, the inhomogeneous neutrino and matter profiles in the
spherical symmetric supernova model do not allow such a simple approach.
Nonetheless, the short-time Fourier transform (or spectrogram analysis)
provides a powerful tool in signal processing to diagnose the characteristic
features in the frequency domain over a short-time window for a time-varying
system, which has been generically used in the engineering and science (see,
e.g., [129, 130, 131]). Based on the same spirit, we can perform similar
spectrogram analysis to the neutrino flavor evolution outcomes except that
here we calculate the spectra in terms of wave number $K_{r}$ over a short-
length window.

The spectrogram analysis can be performed in various ways. In order to
mitigate the spectral leakage, we use the Tukey window function in the
spectrogram analysis [132, 133]. The flat top of Tukey window also naturally
provides less amplitude attenuation than in others such as a Gaussian window.
We divide the whole radial range of our simulation domain into $N_{r}^{\rm
sp}$ blocks with each block centered at a radius $r$, and calculate
discretized Fourier transform for the neutrino off-diagonal mixing at a
snapshot of simulation in each block

$\displaystyle\mathcal{F}_{r,K_{r}}=\frac{1}{\Delta
r^{\mathrm{sp}}}\int_{r-\frac{5}{8}\Delta
r^{\mathrm{sp}}}^{r+\frac{5}{8}\Delta
r^{\mathrm{sp}}}dr^{\prime}\,e^{-iK_{r}r^{\prime}}W(r^{\prime}-r)\times$
$\displaystyle\int
dv_{r}\,\langle\varrho_{e\mu}\rangle_{E}(v_{r},r^{\prime})\mathrm{sgn}[\langle\varrho_{e\mu}\rangle_{E}(v_{r},r^{\prime})-\langle\varrho_{e\mu}\rangle_{E}(v_{r},r^{\prime})],$
(14)

where $\Delta r^{\mathrm{sp}}=(r_{\rm ob}-r_{\rm ib})/N_{r}^{\rm sp}$ is the
radial width of the nonoverlapping part of each block, $K_{r}$ takes
discretized values $z\pi/\Delta r^{\mathrm{sp}}~{}\mathrm{km}^{-1}$ with
integers $z$, and $W(r)$ is the Tukey window function with a cosine lobe

$W(r)=\begin{cases}1&\mathrm{if}~{}|r|\leq\Delta r^{\mathrm{sp}}/2,\\\
\cos^{2}\left(\frac{4\pi|r|}{\Delta
r^{\mathrm{sp}}}-2\pi\right)&\mathrm{if}~{}|r|>\Delta
r^{\mathrm{sp}}/2.\end{cases}$ (15)

Depending on $N_{r}$ used in the simulations, we take $N_{r}^{\rm sp}$ as 50
or 100 to ensure sufficient resolution in $K_{r}$ space.

##  IV Emergence of flavor instability

###  IV.1 Scalability of the dispersion relation

Perturbed by the vacuum flavor mixing (or potentially quantum fluctuation),
the flavor waves of the neutrino fields can be launched and lead to flavor
transformation driven by the dispersion relation and the unstable collective
modes in LSA. Ignoring the contributions from other aspects such as advection,
collisions and the vacuum term, the dispersion relation branches for the
axisymmetric unstable mode can have distinct structures for different ELN
distributions [26]. For a shallow crossing, the dispersion relation typically
yields two separated unstable branches with respect to the wave number of the
perturbed flavor wave $K_{r}$. For a deep crossing, only one unstable branch
exists, and the maximum growth rate of the branch is enhanced.

All LSA results shown below are based on the stationary neutrino profiles
obtained without including oscillations (see e.g., Fig. 2 for the energy-
integrated profiles). The labels v0, v1, v2, and v-2 represent results using
four different values of $\delta m^{2}=0~{}\mathrm{eV}^{2}$, $8\times
10^{-5}~{}\mathrm{eV}^{2}$, $2.4\times 10^{-3}~{}\mathrm{eV}^{2}$ and
$-2.4\times 10^{-3}~{}\mathrm{eV}^{2}$, respectively. The labels c0 and c1
denote settings with $c_{\rm off}=0$ and 1, respectively. Note that a positive
(negative) $\delta m^{2}$ value represents the normal (inverted) neutrino mass
ordering.

We begin with Model II that has the most shallow angular crossing compared to
other models. Figure 4(a–c) shows the growth rates Im$(\Omega)$ of the
unstable mode in Model II with different choices of $\delta m^{2}$, $c_{\rm
off}$, and $a_{\nu\nu}$. From Eqs. (III.2) and (III.2), the dispersion
relation of the FFI is clearly scalable when the vacuum and collision terms
are omitted. In this case, applying the attenuation factor $a_{\nu\nu}$ simply
reduces the magnitude of $\Omega$ and $K_{r}$ for the same amount but does not
affect the shape of the dispersion relation branches.

(a)
(b)
(c)

(d)
(e)
(f)

Figure 4:  Growth rates Im$(\Omega)$ as functions of $K_{r}$ at $r=41.6$ km in
Models II (a–c), and V (d–f). The values of Im$(\Omega)_{\rm max}$ for the
cases with $a_{1}=10^{-2}$ and $a_{1}=10^{-3}$ in (a) and (d) are multiplied
by $10^{2}$ and $10^{3}$, respectively, to show the scalability of the
dispersion relation discussed in the text.

(a) $a_{\nu\nu}=1$, c1v2
(b) $a_{1}=10^{-3}$, c1v1
(c) $a_{1}=10^{-3}$, c1v2

(d) $a_{\nu\nu}=1$, c1v2
(e) $a_{1}=10^{-3}$, c1v1
(f) $a_{1}=10^{-3}$, c1v2

Figure 5:  Growth rates Im$(\Omega)~{}[{\rm km}^{-1}]$ as functions of $K_{r}$
and $r$ in Models II (a–c) and V (d–f). The cyan dashed curve marks $K_{r}$
with maximum growth rates.

The scalability approximately holds even if the collisional terms are
included. Figure 4(a) compares the unstable branches of Im$(\Omega)$ at
$r\approx 41.6$ km in Model II from the case without attenuation
($a_{\nu\nu}=1$) and from cases with attenuation factors $a_{1}=10^{-2}$ and
$a_{1}=10^{-3}$. For the later two with attenuation, their unstable branches
are shown by multiplying factors of $10^{2}$ and $10^{3}$ to both $K_{r}$ and
Im$(\Omega)$, respectively. Among those three cases, there is good agreement
on the rescaled growth rate in the dominant parts of the two separated
unstable branches.

We further examine the effects of collisions and vacuum terms on the
dispersion relation in Model II for the case with $a_{\nu\nu}=1$ and
$a_{1}=10^{-3}$ in Fig. 4(b) and (c) separately. Figure 4(b) shows that
without attenuation, the dispersion relation of FFI is hardly changed in the
presence of collisions and the vacuum term. This suggests that at this radius,
the neutrino number density is so high that other types of flavor
instabilities are overwhelmed by the large growth rate from the fast mode.

When taking the attenuation factor $a_{1}=10^{-3}$ in Fig. 4(c), the shape of
the dispersion relation remains quantitatively similar when including the
collisions in the case of c1v0 or a vacuum term with smaller $\delta m^{2}$ in
c0v1. A slight difference is that small nonzero growth rates Im$(\Omega)\simeq
1~{}{\rm km}^{-1}$, which connect the two otherwise separated branches of pure
fast mode when only the attenuated $\mathbf{H}_{\nu\nu}$ is included (case
c0v0), become visible in comparison to Fig. 4(b). This is because with the
attenuation, $\mathbf{H}_{\nu\nu}$ becomes closer to the $\mathbf{H}_{\rm
vac}$ in magnitude.333These extended nonzero growth rates may also originate
from the spurious mode associated with insufficient number of grids. However,
here we do not further clarify its exact origin, because the flavor conversion
at this radius are dominated by the FFI. When we take a larger vacuum term
with $\delta m^{2}=2.4\times 10^{-3}~{}{\rm eV}^{2}$, the slow flavor
instability becomes so strong that it mixes with the FFI and drastically
affects the structure of the dispersion relation branches. The growth rate of
the left unstable branch is enhanced by a factor of $\simeq 2$, and that of
the right one is suppressed.

For Model V with a much deeper ELN angular crossing at the same radius than
that in Model II, we show the corresponding LSA results in Fig. 4(d-f). In
this case, the dispersion relation has a single unstable branch with larger
growth rates. With unattenuated $\mathbf{H}_{\nu\nu}$, the maximum growth rate
in Model V shown in Fig. 4(d) is $\approx 8.4\times 10^{4}~{}{\rm km}^{-1}$,
which is about seven times larger than that in Model II, and is nearly
independent of the presence of the collisional and vacuum terms [Fig. 4(e)].
Moreover, due to the enhanced growth rate, even when an attenuation with
$a_{1}=10^{-3}$ is taken, the shape of the dispersion relation branch remains
hardly affected by the inclusion of collisional term and different sizes of
vacuum terms as shown in Fig. 4(f).

The scalability of dispersion relation discussed above also holds for a wider
radial range. Figure 5(a) shows the maximum growth rate Im$(\Omega)$ as a
function of $K_{r}$ and $r$ in Model II without attenuation for the case c1v2.
Similaly, two separated unstable branches span over the whole radial region
where there are ELN angular crossings. The maximum growth rate sits on the
unstable branch with smaller (more negative) $K_{r}$ for all radii, with the
corresponding $K_{r}$ value increasing from $K_{r}\simeq-1.5\times
10^{5}~{}{\rm km}^{-1}$ at $r=38$ km to $K_{r}\simeq 0$ km-1 at the outer
boundary as the neutrino flux decreases. The other branch locates in the range
of $0~{}{\rm km}^{-1}\lesssim K_{r}\lesssim 2\times 10^{4}~{}{\rm km}^{-1}$.
For the case with parameters of c1v1 and with $a_{1}=10^{-3}$ shown in Fig.
5(b), one sees that the shapes of the unstable branches remain largely similar
to those shown in Fig. 5(a), when both the $K_{r}$ range and the color scheme
of Im$(\Omega)$ are also rescaled by the same factor $a_{1}$. For the case
c1v2 with a larger $\delta m^{2}$ and with the attenuation parameter
$a_{1}=10^{-3}$, Fig. 5(c) shows that the shape of the dispersion relation for
the whole radial range is affected by the vacuum term. The values of
Im$(\Omega)$ in the lower branch are enhanced while those in the upper one are
suppressed, consistent with the result discussed earlier at $r\approx 41.6$
km. Considering that this change is related to the attenuation of
$\mathbf{H}_{\nu\nu}$ and does not appear in the original dispersion relation
without attenuation, we mainly adopt $\delta m^{2}=8\times 10^{-5}~{}{\rm
eV}^{2}$ in our simulation models to avoid this artifact.

Figure 5(d–f) compares the radial profiles of dispersion relation in Model V
for the case without attenuation to the attenuated cases with c1v1 and c1v2.
For all cases, two unstable branches only appear at $r\simeq 37$ km where a
very shallow angular crossing exists. Those two branches merge into one at
larger radii. In the original case without attenuation in Fig. 5(d), the value
of $K_{r}$ corresponding to the maximum growth rate increases from
$K_{r}\simeq-2.5\times 10^{5}~{}{\rm km}^{-1}$ at $r=37$ km to close to zero
at the outer boundary. The other two cases in Fig. 5(e) and (f) also follow
this trend with a factor of $10^{3}$ smaller in terms of $K_{r}$ and the
growth rate, showing again very good scalability discussed earlier.

Quantitatively, we compare in Fig. 6(a) and (b) the maximum growth rates in
Model II and V for cases with different attenuation factors as well as with
different collision and vacuum terms. Once again, the values of
Im$(\Omega)_{\rm max}$ for the attenuated cases with $a_{1}=10^{-2}$ and
$a_{1}=10^{-3}$ are multiplied by $10^{2}$ and $10^{3}$, respectively. For
Model II, Fig. 6(a) clearly shows that a good scalability can be reached
without attenuation or when the attenuated cases take the smaller $\delta
m^{2}$ (case c1v1). When the larger $\delta m^{2}$ (case c1v2) is adopted for
the attenuated case with $a_{1}=10^{-3}$, it artificially enhances the maximum
growth rates for the whole region with FFI. For Model V shown in Fig. 6(b),
because the angular crossing is deeper and the maximum growth rate is
generally larger than in Model II, the impact due to the vacuum term is
relatively smaller. Even for the case with the larger $\delta m^{2}$ (c1v2)
and $a_{1}=10^{-3}$, the rescaled maximum growth rates closely follow the
original ones, with deviations only up to $\lesssim 20\%$ at $r\gtrsim 60$ km.

Also shown in Fig. 6 are the maximum growth rates estimated by two previously
proposed analytical formulas $\sqrt{2}G_{F}\sqrt{|I_{-}I_{+}|}$ and
$\sqrt{2}G_{F}(2|I_{-}I_{+}|)/(|I_{-}|+|I_{+}|)$, which are proportional to
the geometric and harmonic averages of $|I_{-}|$ and $|I_{+}|$, respectively
[11, 134]. These estimated rates qualitatively agree with each other, as well
as with the values computed by the LSA.

(a)
(b)

Figure 6:  Maximum growth rates Im$(\Omega)_{\rm max}$ as functions of radius
$r$ in Models II (a) and V (b). For the cases with $a_{1}=10^{-2}$ and
$a_{1}=10^{-3}$, the values of Im$(\Omega)_{\rm max}$ are multiplied by
$10^{2}$ and $10^{3}$, respectively, for easier comparison. Black and orange
dotted curves are obtained using the simple analytical formulas in Refs. [11]
and [134], respectively.

With all the discussions above, we conclude that the overall behavior of the
dispersion relation in the region where the ELN angular crossings exist is
predominantly determined by the FFI. The attenuated schemes can retrieve those
structures of unstable modes to a good approximation by taking appropriate
values of $\delta m^{2}$ for specific cases.

###  IV.2 Other type of flavor instability

(a)
(b)

Figure 7:  Maximum growth rates Im$(\Omega)_{\rm max}$ as functions of radius
$r$ for in Models II (a) and V (b). The gray shaded areas indicate the radial
regions where ELN angular crossings exist in the initial ELN distribution
without oscillations. Blue and red curves represent results with parameters in
the cases c0v2 and c1v0 respectively. Dotted and solid curves represent the
cases using $a_{\nu\nu}=1$ without and with excluding the spurious non-
collective modes respectively. The dashed curve in (b) is with $a_{1}=10^{-3}$
and with all the spurious modes excluded.

For regions without angular crossings, it is curious to see whether other type
of flavor instabilities such as the CFI or the slow type may exist. In Model
IV where the angular crossing extends to more inner radius, the FFI
essentially dominates all radii of concern. While for other models, we find
that CFI exists around the neutrinosphere in Model V with the strongest
$Y_{e}$ attenuation leading to the dominance of $\bar{\nu}_{e}$’s over
$\nu_{e}$’s, but it is suppressed in the $\nu_{e}$-dominated Models I–III. For
Models VI and VII with FFI only appearing in limited range inside the
neutrinosphere, their ELN angular profiles outside the neutrinosphere are
similar to Model I and thus not discussed in detail below. We do not find any
slow instabilities that are strong enough to be distinguished from spurious
modes in any of the models.

As two representative instances, we show the maximum growth rates calculated
from the LSA for Models II and V in Fig. 7(a) and (b) for $r<40$ km by the
dotted lines. Although in radii without angular crossings (nonshaded region),
nonzero small growth rates are found for cases c0v2 and c1v0 in both Models II
and V, some of them turn out to be spurious modes that contain noncollective
features. When applying a practical selection rule detailed in Appendix C to
preclude the spurious modes, none of those in Model II survive [see the solid
lines in Fig. 7(a)]. In Model V, the unstable mode due to the nonzero vacuum
term is also spurious, whereas the mode due to CFI persists [solid red line in
Fig. 7(b)]. When the $\mathbf{H}_{\nu\nu}$ attenuation with $a_{1}=10^{-3}$ is
applied, the growth rates due to CFI are nearly the same as the original case
without attenuation in regions without ELN crossings [dashed red line in Fig.
7(b)]. Notice that this CFI mode also appears visibly around $r\simeq 30$ km
over a wide range of $K_{r}$ in Fig. 5(e) and (f).

Some further comments regarding the spurious modes and the CFI are in order.
First, the spurious modes in principle should also exist in simulations since
the LSA results presented here are based on the same number of energy and
angular grids as in simulations. However, we find that these noncollective
spurious modes do not grow in the linear regime but are suppressed by the
angular derivative in the advection term
$(1-v_{r}^{2})r^{-1}\partial_{v_{r}}$. Thus, their presence does not affect
the numerical outcome discussed in the next section.

Second, in contrast to our previous work [71] where CFI were found to be
present in snapshots that are similar to the unattenuated Model I here, the
CFI now only exists in Model V in which $\bar{\nu}_{e}$’s dominate over
$\nu_{e}$’s due to the significant $Y_{e}$ attenuation. We find that the main
distinction is related to the newly added NES process in this work. Without
the NES process, the heavy-lepton neutrinos are thermalized less efficiently
and have higher mean energy than the electron flavor neutrinos (see the
comparison shown in Fig. 22 in Appendix B). Given that all flavors of
neutrinos are still approximately isotropic at this radius, the distinct
behaviors can be understood by examining the quantity

$\gamma_{\rm coll}=\frac{\int
dE\,dv_{r}\,[C_{e\mu}(\varrho_{ee}-\varrho_{\mu\mu})-\bar{C}_{e\mu}(\bar{\varrho}_{ee}-\bar{\varrho}_{\mu\mu})]}{\int
dE\,dv_{r}\,(\varrho_{ee}-\varrho_{\mu\mu}-\bar{\varrho}_{ee}+\bar{\varrho}_{\mu\mu})},$
(16)

which is the imaginary part of the simple formula, Eq. (13) obtained in Ref.
[73]. The denominator, the total ELN, is positive in Models I and II, as well
as in all cases of the Ref. [71], but is negative in Model V. For the
numerator, it is mainly dominated by the first term
$C_{e\mu}(\varrho_{ee}-\varrho_{\mu\mu})$ due to the larger EA rate associated
with $\nu_{e}$. In Ref. [71], the heavy-lepton neutrinos dominate over
$\nu_{e}$’s at higher energy in the region where the CFI can develop. Weighted
by the collisional rate $C_{e\mu}$, which is in general proportional to
$E^{2}$, the high-energy part of the neutrinos contributes predominantly to
the total integral of the numerator and leads to a negative $\gamma_{\rm
coll}$. With a positive denominator, the associated CFI corresponds to the
minus-type defined in Ref. [73], whose growth rate can be well approximated as
$|\gamma_{\rm coll}|$ and as large as a few km-1, mainly determined by the
collisional rates of the corresponding high-energy part. On the other hand, in
Models I and II considered here, the NES process thermalizes the spectra of
heavy-lepton neutrinos so that the mean energy of $\nu_{\mu}$ is close to
$\nu_{e}$. As a result, $\gamma_{\rm coll}$ is positive, which implies that
the corresponding CFI would be the plus-type of Ref. [73] if it exists.
However, unlike the minus-type solution, the CFI may not be assured for the
plus-type. The actual growth rate can be less than $\gamma_{\rm coll}$ or even
completely vanishes, depending on the neutrino spectral properties and the
total ELN (see, e.g., Fig. 2 in Ref. [73]). In Model V, $\bar{\nu}_{e}$’s
dominate over $\nu_{e}$’s with a negative total ELN. With the numerator being
positive as in Models I and II, $\gamma_{\rm coll}$ is negative. As a result,
the CFI is of the minus-type, for which $|\gamma_{\rm coll}|$ is a good
approximation to indicate the corresponding growth rate.

The above discussion suggests the following intriguing fact. Although the EA
processes in Table 1 only involve the electron flavor neutrinos from the
classical viewpoint, it can damp the coherent state between $\nu_{e}$ and
$\nu_{\mu}$ in $\nu$QKE so that the CFI does depend on the properties of
heavy-lepton neutrinos. We also note that although CFI is suppressed in our
models dominated by $\nu_{e}$’s, it does not imply that CFI is generally not
present in all $\nu_{e}$-dominated conditions. Multidimensional CCSN
simulations with convection may potentially result in conditions favoring the
minus-type CFI when the dynamical evolving heavy-lepton neutrinos are less
efficiently thermalized by the NES process.

###  IV.3 Spectrogram

In this subsection we perform the analysis based on the simulation results
using the spectrogram defined in Sec. III.3 when flavor conversions reach the
nonlinear regime, and compare the structures to the dispersion relation
obtained by the LSA in the previous subsections. Note that here the model
parameters are given in Table 2, and the $\mathbf{H}_{\nu\nu}$ attenuation is
applied to all the models.

Before presenting the spectrogram analysis, we first illustrate the evolution
of FFC from the onset of FFI in Fig. 8 for Model II as an example. Since the
maximum growth rate of FFI, Im$(\Omega)_{\rm max}$, in this model peaks at
$r\simeq 42$ km, FFI grows faster to reach the nonlinear regime around that
radius first. The left panel of Fig. 8 shows that $s_{e\mu}$, the
dimensionless ratio representing the off-diagonal flavor mixing introduced in
Sec. III.1, reaches the nonlinear regime over a wide range from $r\simeq 40$
km to $r\simeq 58$ km by the time of $t=32~{}\mu$s. The range over which the
nonlinear FFC takes place further widens to $r\simeq 70$ km at $t=64~{}\mu$s
and extends to the outer boundary at a much later time $t=160~{}\mu$s.

Figure 8:  Evolution of the dimensionless ratio of off-diagonal mixing
$s_{e\mu}$ in $30$ km $<r<80$ km (left panel), and its associated phase angle
$\phi_{e\mu}$ within a smaller radial range $50$ km $<r<55$ km (right panel)
for Model II.

(a) Model II, $32~{}\mu$s
(b) Model II, $64~{}\mu$s
(c) Model II, $160~{}\mu$s
(d) Model IIa2, $32~{}\mu$s
(e) Model IIa4, $32~{}\mu$s
(f) Model IIa10, $32~{}\mu$s
(g) Model V, $1.6~{}\mu$s
(h) Model V, $6.4~{}\mu$s
(i) Model V, $32~{}\mu$s
(j) Model Vc0, $32~{}\mu$s
(k) Model III, $160~{}\mu$s
(l) Model IV, $160~{}\mu$s

Figure 9:  Spectrograms with respect to $K_{r}$ and $r$ for different models
at different times. The color bar is for the Fourier transformed magnitude
$\log_{10}\big{|}\mathcal{F}_{r,K_{r}}[{\rm cm}^{-3}]\big{|}$. The cyan dashed
curve shows the $K_{r}$ with maximum growth rates calculated from LSA.

On the right panel in Fig. 8, we also show the radial shape of the complex
quantum phase $\phi_{e\mu}$ defined in Sec. III.1 in $50$ km $<r<55$ km to
illustrate the coherent feature of FFC. The phase angle $\phi_{e\mu}$ not only
oscillates rapidly in time during FFCs, but also develops coherent spatial-
dependent oscillatory patterns. Moreover, the pattern is highly correlated for
a wide range of $v_{r}$ due to the collective nature of FFC. The typical
oscillation length scale of the synchronized phase angle is $\lesssim 0.2$ km,
broadly consistent with the scale determined by the maximally unstable $K_{r}$
shown in Fig. 5(b) in the same radial range. This simple comparison suggests
that the characteristic length scale determined by the dispersion relation in
the linear regime correlates with the coherent length scale at late times when
FFC are in the nonlinear or asymptotic regime [135, 33], which will be further
discussed below with the spectrograms. We also note that this spatial
variation can lead to a substantial contribution to the flavor evolution
equation through the advection term $v_{r}\partial_{r}$, and thus participate
in determining the outcome of flavor evolution.

The spectrograms shown in Fig. 9(a–c) compare the Fourier transformed
magnitudes $\mathcal{F}_{r,K_{r}}$ at the same three time frames in Fig. 8
with the curve of $K_{r}$ corresponding to the maximum growth rate from LSA.
At $t=32~{}\mu$s, the dominant wave number $K_{r}$ of the off-diagonal flavor
mixing in $\mathcal{F}_{r,K_{r}}$ generally aligns well with the maximally
unstable $K_{r}$ values derived with LSA. Relatively small differences of
$\simeq 20$–30 km-1 that exist are mainly related to the change of ELN
distribution due to flavor conversions in simulation, while the LSA results
are based on unoscillated neutrino profiles. We note that near $r=40$ km, the
range of dominant $K_{r}$ widens over time to encloses the homogeneous modes
with $K_{r}\approx 0~{}{\rm km}^{-1}$. The amplitude of the off-diagonal
flavor mixing between $r\approx 42$ km and 50 km decreases at a later time
$t=160~{}\mu$s, with $\mathcal{F}_{r,K_{r}}\lesssim 10^{30}~{}{\rm cm}^{-3}$.
For $r\gtrsim 50$ km, similar to the case with a shallow angular crossing
studied under the periodic boundary [36], our result suggests that a steady
wave-like pattern along the radial direction can also survive in a more
realistic setting without adopting the periodic boundary condition.

The radial-dependent coherent pattern of the flavor mixing remains
qualitatively similar when different attenuation factors are used. Figure
9(d–f) compares the spectrograms for Models IIa2, IIa4, and IIa10 at the same
time $t=32~{}\mu$s. The features of the spectrograms in all three models again
agree well with the results derived from the LSA. The patterns also show a
good scalability as discussed before. Because models with less attenuation
evolve faster, they reach the asymptotic pattern at earlier times as shown in
Fig. 9(e) and (f).

In Model V, because the maximum growth rates are generally larger than those
in Model II, flavor conversions occur faster. Figure 9(g) shows that a
significant flavor mixing with $\mathcal{F}_{r,K_{r}}>10^{31}~{}{\rm cm}^{-3}$
is attained for $r\approx 42$ km at $t=1.6~{}\mu$s. At $t=6.4~{}\mu$s shown in
Fig. 9(h), the FFC in the region from $\approx 36$ km to $60$ km evolves to
the nonlinear regime. The dominant wave numbers in $r\lesssim 40$ km and
$r\gtrsim 50$ km still qualitatively follow the curve of maximum-growth-rate
$K_{r}$ from the LSA. However, in the region between 40 and 50 km, the
spectrogram covers a much wider wave number modes, spreading from
$K_{r}\simeq-300~{}{\rm km}^{-1}$ to $300~{}{\rm km}^{-1}$. This becomes even
more evident at $t=32~{}\mu$s shown in Fig. 9(i). Also at $t=32~{}\mu$s, the
emergence of CFI leads to flavor mixing between $r\approx 27$ km to 33 km with
the unstable modes centered on $K_{r}\simeq 0~{}{\rm km}^{-1}$.

We use Model Vc0 to disentangle the effect due to the interplay between the
CFI and the FFI. Figure 9(j) shows that the asymptotic spectrogram of Model
Vc0 at the same time $t=32~{}\mu$s generally follows the trend of $K_{r}$ from
LSA. The dominant $K_{r}$ of the spectrogram spans from
$K_{r}\simeq-250~{}{\rm km}^{-1}$ at $r\approx 36$ km to
$K_{r}\simeq-30~{}{\rm km}^{-1}$ at $r\approx 46$ km. Between $r\simeq 40$ km
and $60$ km, the range of $K_{r}$ in the spectrogram becomes wider, consistent
with similar features shown in Fig 9(h) and (i) in the same radial range.

Finally, we also show in Fig. 9(k) and (l) the asymptotic spectrograms in
Models III and IV, respectively. For both models, the radial shape of the
dominant $K_{r}$ modes again agree qualitatively with the LSA results,
although quantitative differences do exist. In any case, our spectrogram
analysis reveals that asymmetric patterns with respect to $K_{r}$ generally
exist, which suggests that the outcome of FFCs and the flavor waves are not
necessarily homogeneous. The overall development of small-scale structures
reported here also qualitatively agrees with the findings in local simulations
adopting periodic boundary conditions.

Besides the off-diagonal flavor mixing, the diagonal elements of the density
matrices associated with the neutrino number density also contain the same
small-scale structures. In Fig. 10, we show the neutrino number densities of
Models III and IIIa2 between $r=38$ km and 40 km for illustration. Clearly,
oscillatory patterns exist in general for all three different neutrino
species. In Model III, its typical oscillatory wavelength is $\simeq 0.07$ km,
corresponding to a $K_{r}\simeq 90~{}{\rm km}^{-1}$, which is consistent with
the dominant mode of spectrogram in Fig. 9(k). With a less attenuation in
Model IIIa2, the oscillatory length scale become correspondingly smaller. The
presence of small-scale structure may lead to implications associated with the
refractive energy and the potential self-thermalization of the neutrino gas
[136, 137].

Figure 10:  Radial profiles of neutrino number densities between $r=38$ km and
40 km at $t=64~{}\mu$s in Models III and IIIa2.

##  V Evolution of fast flavor conversions

(a)
(b)
(c)
(d)
(e)
(f)

Figure 11:  Radial profiles of neutrino number densities for $\nu_{e}$,
$\bar{\nu}_{e}$, and $\nu_{\mu}$ in Models II–VII (a–f).

(a)
(b)
(c)
(d)
(e)
(f)

Figure 12:  Radial profiles of neutrino mean energies for $\nu_{e}$,
$\bar{\nu}_{e}$, $\nu_{\mu}$, and $\bar{\nu}_{\mu}$ in Models II–VII (a–f).

In this section, we will present more details regarding the evolution of the
FFC and its impact on key properties of neutrinos. Figures 11 and 12 show the
neutrino number density and mean energy of different species for Models II–VII
at different times. For the number densities, because the initial values of
$\nu_{\mu}$ and $\bar{\nu}_{\mu}$ are rather close to each other, and given
that the fast pairwise flavor conversion produces equal amount of $\nu_{\mu}$
and $\bar{\nu}_{\mu}$, only the $\nu_{\mu}$ ones are shown in Fig. 11. For the
average energies, all fours species are shown explicitly in Fig. 12. Also
noted is that the values shown in these two figures are coarse-grained ones,
computed by averaging over a size of 0.6 km from the original simulation
outputs, for the purpose of capturing the global trend. Below, we will first
compare the similarities and differences between those models. Different
adopted schemes such as the number of neutrino flavors, the amount of
attenuation on $\mathbf{H}_{\nu\nu}$, and the size of the vacuum term will
then be discussed.

###  V.1 Model II

We start with Model II in which the angular crossing is the shallowest. Prior
to the occurrence of the FFC, the number densities of all neutrino species
decrease smoothly with radius, following a hierarchy of $n_{\nu_{e}}\gtrsim
n_{\bar{\nu}_{e}}>n_{\nu_{\mu}}$ in the radial range shown in Fig. 11(a). Due
to the FFI discussed earlier, the FFC takes place and leads to the pairwise
conversion from $\nu_{e}$ and $\bar{\nu}_{e}$ to muon flavor neutrinos. At
$t=32~{}\mu$s, the region between $r\approx 40$ km and $r\approx 57~{}$km is
influenced by FFC. The impact of FFC extends to $r\approx 70$ km at
$t=64~{}\mu$s and eventually reaches the outer boundary at a later time.
Notice that there are “spiky” features just behind the outermost range
affected by FFC, e.g., at $r\approx 56$ km at $t=32~{}\mu$s and $r\approx 68$
km at $t=64~{}\mu$s. Those are consistent with what reported in local
simulations adopting periodic boundaries where a larger amount of flavor
conversions can exist as a transient phenomenon before the system settles down
to the asymptotic state.

In Model II, $n_{\nu_{e}}$ and $n_{\bar{\nu}_{e}}$ in the region affected by
FFC are reduced by $\simeq 10\%$ to $20\%$, while $n_{\nu_{\mu}}$ increases by
$\simeq 50\%$–$70\%$. Once the asymptotic state in a certain radial range is
reached, the number densities do not vary much anymore, illustrated by
comparing the curves at $t=64~{}\mu$s and at $t=320~{}\mu$s in $40$
km$\lesssim r\lesssim 68$ km. Clearly, flavor equilibration is not achieved in
Model II. This is because the FFC here is restricted within a limited angular
range (see e.g., Fig. 8 and later in Fig. 14).

For the mean energies of $\nu_{\mu}$ and $\bar{\nu}_{\mu}$, differently from
the similar evolution of $n_{\nu_{\mu}}$ and $n_{\bar{\nu}_{\mu}}$, they can
be affected very differently by FFCs. Figure 12(a) shows that $\langle
E_{\nu_{\mu}}\rangle$ is reduced substantially to $\approx 16$–17 MeV for
$r>40$ km, while $\langle E_{\bar{\nu}_{\mu}}\rangle$ remains around $18$–20
MeV. For $\nu_{e}$ and $\bar{\nu}_{e}$, their mean energies increase by less
than 1 MeV.

###  V.2 Models III and IV

For Models III and IV that possess deeper angular crossings in terms of
$v_{r}$, the ELN crossings also appear at radii deeper inside the
neutrinospheres (see Fig. 2). We show in Fig. 11(b) and (c) the evolution of
the neutrino number densities. At $t=32~{}\mu$s, FFCs take place between
$r\approx 32$ km to $r\approx 76$ km in Model III, and from $r\approx 28$ km
till the outer boundary in Model IV. The deeper angular crossings in both
models lead to an increased amount of the converted neutrinos from the
electron to the muon flavor than that in Model II. Particularly, the nearly
equal number densities in 50 km$\lesssim r\lesssim 70$ km of all different
flavors indicate that flavor equilibration is approximately achieved in that
radial range in both models. The zig-zag pattern around $r=70$–72 km in Model
III is associated with the transient phenomenon from the linear to nonlinear
regime aforementioned in Model II.

Differently from Model II, the existence of the ELN angular crossings at
deeper radii inside the neutrino spheres of $\nu_{e}$ enables an intriguing
feedback effect between the FFC and the collisional EA processes after the
initial phase of FFC in both models. In Model III, the number densities of
$\nu_{e}$ and $\bar{\nu}_{e}$ around $r\simeq 45$ km at $t=64~{}\mu$s increase
by a similar amount from those at $t=32~{}\mu$s as shown in Fig. 11(b). This
is because the EA collisions around the neutrinospheres of both $\nu_{e}$ and
$\bar{\nu}_{e}$ act to restore their number densities to the stationary values
prior to the occurrence of FFC. The increase of $n_{\nu_{e}}$ and
$n_{\bar{\nu}_{e}}$ continues and also affects regions outside the
neutrinosphere. For instance, at $t=160~{}\mu$s, $n_{\nu_{e}}$ and
$n_{\bar{\nu}_{e}}$ become closer to their initial values for $r\gtrsim 40$
km. Similar behaviors also appear in Model IV shown in Fig. 11(c). Nearly
flavor equilibration is achieved at $t=32~{}\mu$s initially, followed by the
increase of $n_{\nu_{e}}$ and $n_{\bar{\nu}_{e}}$ locally in $40$ km$\lesssim
r\lesssim 66$ km at $t=128~{}\mu$s, and later at $t=320~{}\mu$s the
restoration of $n_{\nu_{e}}$ and $n_{\bar{\nu}_{e}}$ back to the initial
values for $r\gtrsim 40$ km. The restoration is stronger than in Model III,
particularly for $n_{\nu_{e}}$, which almost overlaps with the unoscillated
level. We note that because this feedback involves mostly the collisional
processes and advection, the associated timescale of the number density
increase is not related to the conversion rate of FFC. For instance, in Model
IV, the range affected by this feedback effect extends from $r\simeq 40$ km at
$t=32~{}\mu$s to $r\approx 66$ km at $t=128~{}\mu$s, consistent with the
corresponding neutrino propagation time, which is a few tens of km/c.

In addition to the restoration of $\nu_{e}$ and $\bar{\nu}_{e}$ number
densities discussed above, the feedback effect contains an interesting twist
in reshaping the ELN angular distribution, which then leads to different local
ELN angular distributions compared to those after the prompt flavor
conversions. One major difference is related to the conservation of the
neutrino lepton number. In the prompt FFCs, in addition to the elimination of
the ELN angular crossings, the neutrino and antineutrino lepton numbers, as
well as the total ELN are locally conserved to a good approximation. However,
all these three numbers are no longer conserved quantities on a longer
timescale when the collision and neutrino diffusion play a role. Although the
FFC still works to erase the ELN angular crossings, the collisional weak
processes such as EA can change those lepton numbers as well as the total ELN.
Consequently, the angular distributions of the ELN and each neutrino species
keep evolving on timescales determined by the collision rates, while the FFC
functions to ensure that no new ELN angular crossings appear during the
evolution. The system can then settles down to a quasi-steady state determined
by the FFC, collisions, and advection all together.

To illustrate this effect, we first show in Fig. 13 the neutrino angular
distributions for different times at $r\approx 36$ km, which is inside the
$\nu_{e}$-sphere in both Model III and IV. In Model III, the initial ELN
distribution transitions from positive to negative at $v_{r,c}\approx 0.37$
[Fig. 13(a)]. With the condition $|I_{+}|>|I_{-}|$, the prompt FFC between the
electron and muon flavors erases the ELN crossing by making it positive in the
entire range of $v_{r}$, accompanied with significant flavor swap for
$v_{r}>0$ at $t=3.2~{}\mu$s [Fig. 13(b)]. Since this radius is inside the
neutrinosphere, the EA processes of $\nu_{e}$ and $\bar{\nu}_{e}$ with matter
tend to restore their angular distributions back to the initial state prior to
FFC. Due to the higher opacity for $\nu_{e}$’s than $\bar{\nu}_{e}$’s, the
angular distribution of $\nu_{e}$’s increases faster, which leads to a more
positive ELN distributions at $t=16~{}\mu$s [Fig. 13(a)]. On a longer
timescale, the restoration of $\nu_{e}$’s slows down as its angular
distribution becomes closer to the initial state. Due to the restoration of
$\bar{\nu}_{e}$’s, the ELN distribution starts to decrease, which potentially
leads to a new angular crossing. However, the FFC can smear any newly formed
angular crossing and keep the ELN angular distribution positive during this
process. Finally, the neutrino distributions settle into a quasi-steady state
at $t\approx 64~{}\mu$s where the ELN angular distribution is nearly flat for
$v_{r}\gtrsim 0.5$ [Fig. 13(a)]. Meanwhile, the $\nu_{e}$ and $\bar{\nu}_{e}$
distributions are restored to values closer to their initial states for
$v_{r}\lesssim 0.5$ [Fig. 13(b)].

(a)
(c)
(b)
(d)

Figure 13:  Evolution of the angular distributions of the ELN, $G(v_{r})$
(upper panels), and that of each neutrino species at $r\approx 36$ km (lower
panels) in Models III (a–b) and IV (c–d).

(a)
(c)
(e)
(b)
(d)
(f)

Figure 14:  Evolution of the angular distributions of the ELN, $G(v_{r})$
(upper panels), and that of each neutrino species at $r\approx 50$ km (lower
panels) in Models II (a–b), III (c–d), and IV (e–f). Note that
$\langle\varrho_{\nu}\rangle_{E}$ are shown in $-0.2\leq v_{r}\leq 1$ because
most of neutrinos propagate outwards.

For Model IV, the neutrino distribution is dominated by $\bar{\nu}_{e}$’s with
$|I_{+}|<|I_{-}|$ at the same radius $r\simeq 36$ km with an initial crossing
at $v_{r}\simeq-0.2$. Thus, flavor equilibration takes place mainly for
negative $v_{r}$ at $t=3.2~{}\mu$s [Fig. 13(c)], and the ELN angular
distribution becomes essentially negative for all $v_{r}$’s. However, the
faster EA rates for $\nu_{e}$’s increase the ELN distribution and turn it into
positive values at $t=16~{}\mu$s [Fig. 13(c)]. Later, the effect of slower
$\bar{\nu}_{e}$ rates take effects, and lowers the values of ELN distribution
for a major $v_{r}$ range. At $t=64~{}\mu$s, the distributions also settle
into an approximately steady state, with nearly flat ELN distribution around
zero for $v_{r}\gtrsim 0$ [Fig. 13(c)], as well as the
$\langle\varrho_{\nu_{e}}\rangle_{E}\simeq\langle\varrho_{\bar{\nu}_{e}}\rangle_{E}$
and
$\langle\varrho_{\nu_{\mu}}\rangle_{E}\simeq\langle\varrho_{\bar{\nu}_{\mu}}\rangle_{E}$
[Fig. 13(d)].

The collisional feedback effect inside the neutrinosphere also affects the
flavor evolution in the outer part. Figure 14(c–f) shows that in both Model
III and IV, the prompt FFC at $r=50$ km and $t=32~{}\mu$s turns the ELN
distribution to the same sign, with significant amount of $\nu_{e}$’s and
$\bar{\nu}_{e}$’s being converted. However, as neutrinos continuously stream
out, what are observed at $r=50$ km at later times reflect the state of
neutrinos escaping from inside the neutrinosphere at earlier times. At
$t=160~{}\mu$s, the angular distributions are restored back to $\simeq
80\%$–90% of the initial values for $\nu_{e}$’s and $\bar{\nu}_{e}$’s in
$v_{r}\gtrsim 0.8$ in Model III, or even beyond for $\nu_{e}$’s in
$v_{r}\gtrsim 0.9$ in Model IV. The asymptotic ELN distributions are almost
zero for $v_{r}\gtrsim 0.7$ in both Models, showing very different shapes from
the prompt FFC.

For comparison, we also show in Fig. 14(a–b) the corresponding distributions
in Model II. From $t=32~{}\mu$s to $t=160~{}\mu$s, the angular distributions
of $\nu_{e}$’s and $\bar{\nu}_{e}$’s are only slightly changed. This is
because in Model II, the FFCs only start at $r\approx 40$ km located just
outside the $\nu_{e}$-sphere. No significant amount of $\nu_{e}$’s and
$\bar{\nu}_{e}$’s are produced additionally from within the neutrinosphere.

In spite of the intriguing feedback effect on the neutrino distributions
discussed above, the main impact of FFC on the radial profiles of the neutrino
mean energies in Model III and IV are relatively similar to that in Model II,
as shown in Fig. 12(b) and (c). FFC slightly enhances $\langle
E_{\nu_{e}}\rangle$ and $\langle E_{\nu_{e}}\rangle$ by $\lesssim 1$ MeV in
all radii, while it largely reduces $\langle E_{\nu_{\mu}}\rangle$ by $\simeq
2$ MeV in regions for $r\gtrsim 35$ km. One visible difference from Model II
here is that for $r\lesssim 35$ km, both $\langle E_{\nu_{\mu}}\rangle$ and
$\langle E_{\bar{\nu}_{\mu}}\rangle$ are reduced due to the backward-
propagating neutrinos from outer radii.

###  V.3 Model V

Model V is special because of the coexistence of CFI and FFI. Before
discussing the combined effect originated from both flavor instabilities, we
first discuss the results in the corresponding Model Vc0 where the the CFI is
switched off. Figure 15 shows the radial profiles of the neutrino number
densities and their mean energies in Model Vc0. The impact of FFCs here is
qualitatively similar to those discussed earlier. At $t=32~{}\mu$s, the prompt
FFC affects the radial range from $r\approx 35$ km to the outer boundary and
leads to approximate flavor equilibration in $40$ km$\lesssim r\lesssim 54$ km
where the deepness ratio of the initial ELN, rewritten as $|I_{+}/I_{-}|$, is
larger (see Fig. 3). As neutrinos continue propagating outwards, their flavor
content is replaced by the emission and oscillation properties close to the
neutrinosphere. Together with the feedback mechanism related to the EA
processes, the number densities of $\nu_{e}$ and $\bar{\nu}_{e}$ are also
enhanced for $r\gtrsim 40$ km at $t=320~{}\mu$s as in Models III and IV. For
$r\gtrsim 50$ km, the mean energies of $\nu_{\mu}$ and $\bar{\nu}_{\mu}$
become much closer to those of $\nu_{e}$ and $\bar{\nu}_{e}$, respectively.
Because the angular crossings for radii inside the $\nu_{e}$ neutrinosphere in
Model Vc0 is located at negative $v_{r}$, the FFC affects more the backward
propagating $\nu_{\mu}$’s in that region, which leads to a minor enhancement
of $n_{\nu_{\mu}}$ near $r=32$ km as well as moderate reductions of the mean
energies for heavy-lepton neutrinos there.

(a)
(b)

Figure 15:  Radial profiles of neutrino number densities (a) and mean energies
(b) in Model Vc0.

For Model V where the collisional flavor conversions is also present, Fig.
11(d) shows that in addition to the same FFCs in $r\gtrsim 35$ km, the CFI
also leads to flavor transformations and thus the enhanced $n_{\nu_{\mu}}$ in
$r\simeq 28$–35 km inside the neutrinosphere. As a result, $n_{\nu_{\mu}}$
becomes almost the same as $n_{\nu_{e}}$ already at $t=32~{}\mu$s and remains
till $t=320~{}\mu$s. For the mean energies, once again, the muon flavor
neutrinos become closer to those of electron flavor as shown in Fig. 12(f).

We notice that in both Models V and Vc0 at $t=320~{}\mu$s, the neutrino number
densities in regions for $r\gtrsim 50$ km become $n_{\nu_{e}}\simeq
n_{\bar{\nu}_{e}}$, different from the initial hierarchy of
$n_{\nu_{e}}<n_{\bar{\nu}_{e}}$. Because the flavor content is still dominated
by $\bar{\nu}_{e}$’s after the prompt conversions, this effect cannot be
attributed to the pairwise FFC. It suggests that the collisional feedback
effect inside the neutrinosphere leads to a larger net change of flavor
content in $\bar{\nu}_{e}$’s than $\nu_{e}$’s.

###  V.4 Models VI and VII

For Models VI and VII, the initial ELN angular crossings only reside in a
limited radial range at 29 km $\lesssim r\lesssim 33$ km and 27 km $\lesssim
r\lesssim 29$ km (see Fig. 2), respectively. Figure 11(e) and (f) shows that
the FFC affects mostly the number density of heavy-lepton flavors outside the
regions where the crossings exist initially. For $\nu_{e}$ and
$\bar{\nu}_{e}$, their number densities remain nearly unchanged for all radii
due to the following reasons. First of all, since these crossings are inside
the neutrinospheres, the neutrino-matter interaction outside the FFI region
plays an important role in resetting the properties of $\nu_{e}$ and
$\bar{\nu}_{e}$. Secondly, as the neutrino angular distributions are more
isotropic inside the neutrinosphere, neutrinos not affected by FFCs can also
diffuse back into the region of FFC to mitigate the impact of flavor
conversions. Moreover, even inside the region of FFC, the dynamical feedback
from the collisions mentioned above in Models III and IV also acts to restore
the distributions of $\nu_{e}$ and $\bar{\nu}_{e}$ back to the level close to
those before oscillations. For the heavy-lepton neutrinos, their spectra at
inner radii become closer to those of $\nu_{e}$ and $\bar{\nu}_{e}$ when the
system reaches the quasi-steady state, which limits further flavor
conversions. Since the FFI region in Model VI is closer to the neutrinosphere
than in Model VII, the overall impact of flavor conversions is larger in Model
VI.

Similar behaviors are found in terms of neutrino mean energies as well. In
both models, the mean energy for heavy-lepton neutrinos are reduced
significantly, while the mean energies of $\nu_{e}$ and $\bar{\nu}_{e}$ only
increase very little. The impact in Model VI is larger than in Model VII.

###  V.5 Two flavors v.s. three flavors

To illustrate the effects of adopting schemes with different numbers of
neutrino flavors in the simulations, we choose Models II and IV as the
representative examples for cases without and with the dynamical collisional
feedback. The results from Models II and IV are compared with Models IIf3 and
IVf3 in Fig. 16. For Models II and IIf3, the radial profile of the number
density ratios $n_{\nu}(t=320~{}\mu{\rm s})/n_{\nu}(t=0)$ for all flavors,
including $\nu_{\tau}$ and $\bar{\nu}_{e}$ in the three-flavor scenario, are
shown in Fig. 16(a). Note that we use different scales for the electron
flavors (left axis) and the heavy-lepton flavors (right axis) for better
readability. For Models IV and IVf3, the number density ratios at
$t=32~{}\mu$s and $320~{}\mu$s are shown in Fig. 16(b) and (c), respectively,
to illustrate the impact after the prompt FFC and that at the asymptotic state
after the collisional feedback takes place.

These plots show that with three flavors, the enhancement of $\nu_{\mu}$ and
$\bar{\nu}_{\mu}$ number densities is reduced, due to the additional flavor
conversion channel from electron to tauon flavors. The enhancement in Model
IIf3 is $\simeq 50\%$ of that in Model II. For Models IVf3 and Model IV, the
enhancement of heavy-lepton flavor number densities in Model IVf3 is reduced
by less than $50\%$, particularly at $t=320~{}\mu$s. This implies that the
flavor conversions result in a larger amount of the total heavy-lepton flavor
number density in Model IVf3 than in Model IV.

(a)
(b)
(c)

Figure 16:  Comparison of the radial profiles of neutrino number density
ratios $n_{\nu}(t)/n_{\nu}(t=0)$ between two-flavor (dashed) and three-flavor
(solid) schemes for Model II at $320~{}\mu$s (a), Model IV at $32~{}\mu$s (b)
and at $320~{}\mu$s (c). Note that the scales for heavy-lepton neutrinos are
labeled on the right axes.

(a)
(b)
(c)

Figure 17:  Comparison of the radial profiles of neutrino number density
ratios $n_{\nu}(t)/n_{\nu}(t=0)$ (blue) in Models II (a), III (b), and IV (c)
at $64~{}\mu$s with different attenuation factors $a_{1}=10^{-3}$ (blue),
$2\times 10^{-3}$ (green), $4\times 10^{-3}$ (red), and $1\times 10^{-2}$
(black in Model II). Note that the scales for the ratio of $\nu_{\mu}$ are
labeled on the right axes.  Figure 18:  Comparison of the radial profiles of
neutrino number density ratios $n_{\nu}(t)/n_{\nu}(t=0)$ between Models IV and
IVv2 at $64~{}\mu$s. The scale for the ratio of $\nu_{\mu}$ is labeled on the
right axis.

For $\nu_{e}$ and $\bar{\nu}_{e}$, their number densities are generally more
reduced by flavor conversions, while the amount of reduction in Models IIf3
and IVf3 from the corresponding two-flavor cases differs. In Model IIf3, the
reduction is relatively minor and only appears at $r\gtrsim 50$ km. In Model
IVf3, the reduction at $t=32~{}\mu$s is much more significant than in Model
IIf3, reflecting the shift of $\nu_{e}$ and $\bar{\nu}_{e}$ survival
probability from $\simeq 1/2$ to $\simeq 1/3$ when reaching flavor
equilibration (see, e.g. [33]). At $t=320~{}\mu$s, although the $\nu_{e}$’s
and $\bar{\nu}_{e}$’s number densities are restored back to values closer to
those prior to FFCs due to the collisional feedback effect, the three-flavor
case still shows more reduction than the two-flavor case. Overall, including
the tauon flavor allows more flavor conversion of $\nu_{e}$’s and
$\bar{\nu}_{e}$’s, which is consistent with the observed increase of the total
heavy-lepton flavor content.

###  V.6 Attenuation on the neutrino self-induced term

We discussed in Sec. IV the impact of $\mathbf{H}_{\nu\nu}$ attenuation on the
dispersion relation and the small-scale structure. In this subsection, we
examine its impact on the flavor evolution in macroscopic scale. Figure
17(a)–(c) shows the neutrino number density ratios in Models II–IV at
$t=64~{}\mu$s with different attenuation factors $a_{1}$. At this time, only
the case in Model II with $a_{1}=10^{-3}$ shows substantial differences in the
number density ratios at $r\gtrsim 70$ km. The reason is that the evolution in
Model II is mainly determined by the prompt FFCs. Thus, cases with more
attenuation (smaller $a_{1}$) evolve slower than those with less attenuation.
Another slight difference is that the starting radius of FFCs with less
attenuation is slightly shifted leftwards from $r\approx 39$ km to $\approx
38$ km.

For Models III and IV, the number density ratios in Fig. 17(b) and (c) exhibit
a transition at $r\simeq 50$ km, within which the impact of collisional
feedback appears. Since the timescale of the collisional feedback effect does
not depend on the attenuation factor, provided that the FFC timescale remains
short enough, the results shown in these two panels are therefore hardly
affected by the choice of different attenuation factors.

###  V.7 Vacuum term

We also discussed in Sec. IV that the general trend of the full dispersion
relation is almost independent of the vacuum term when no
$\mathbf{H}_{\nu\nu}$ attenuation is applied. When taking $a_{1}=10^{-3}$, the
dispersion relation in Model II, which has the smallest instability growth
rate, can be affected by the adopted value of $\delta m^{2}$, while the other
models such as IV and V remain nearly unaffected. Here, we examine whether the
simulation results at late times are affected by the strength of the vacuum
term by comparing Model IV to Model IVv2 at $t=64~{}\mu$s in Fig. 18. We find
that taking a different $\delta m^{2}$ indeed only results in negligible
effects, which is consistent with the findings reported in Refs. [128, 66,
138]. For both the radial regions affected by the prompt flavor conversion and
that affected by the collisions feedback, a larger $\delta m^{2}$ only leads
to insignificant deviations.

##  VI Effects and observables

After discussing the evolution of FFCs, we examine in this section their
implications and feedback to key quantities relevant to CCSN physics,
including the net neutrino heating rates and the equilibrium electron
fraction. In addition, the impact of FFCs on the free-streaming neutrino
energy spectra, which are important inputs for the expected supernova neutrino
signals, will be discussed.

###  VI.1 Heating-cooling rates

The specific net heating rate dominated by the EA process of $\nu_{e}$ and
$\bar{\nu}_{e}$, processes (1a–1c) in Table 1, is given by

$\displaystyle\dot{q}=\frac{1}{\rho}\int dE\,$ $\displaystyle E[\chi_{e,\rm
EA}\varrho_{ee}-j_{e,\rm EA}(\varrho_{ee,\rm FO}-\varrho_{ee})$
$\displaystyle+\bar{\chi}_{e,\rm EA}\bar{\varrho}_{ee}-\bar{j}_{e,\rm
EA}(\bar{\varrho}_{ee,\rm FO}-\bar{\varrho}_{ee})].$ (17)

We show the above net heating rate in Models II–V for cases with and without
FFCs in Fig. 19. Notice that the rates before FFCs (black curve) are
practically indistinguishable on the plot for all models despite different
attenuation strengths adopted to $Y_{e}$ profiles. This is because the
reduction of $\nu_{e}$’s and enhancement of $\bar{\nu}_{e}$’s are
complementary to each other. Without FFCs, the regions with net cooling
($\dot{q}<0$) and heating ($\dot{q}>0$) are separated by the gain radius at
$r\approx 63$ km. In the cooling region, the maximum net cooling rate is as
low as $\dot{q}\approx-8.4\times 10^{21}~{}{\rm erg}\,{\rm g}^{-1}\,{\rm
s}^{-1}$ at $r\simeq 46$ km.

(a)
(b)
(c)

Figure 19:  Comparison of the specific heating-cooling rate after the prompt
FFCs at $t=32~{}\mu$s (a) and in the asymptotic state at $t=160~{}\mu$s (b)
for Models II (blue), III (green), IV (cyan), and V (red). The original
heating rates without flavor conversions are shown by the black curves and are
indistinguishable for all models. The blue dashed curve in (a) is for Model IV
at $t=3.2~{}\mu$s. The magenta curve in (c) is for Model Vc0 at $t=32~{}\mu$s,
compared to Model V (red).

With FFCs, other solid curves in Fig. 19(a) show the radial profiles of
$\dot{q}$ at $t=32~{}\mu$s for Models II–V. Since the pairwise FFCs from
electron to heavy-lepton flavors lead to the reduction of heating rates in
these models, $\dot{q}$ are reduced in the entire domain. As a result, the
gain radius shifts to larger radii around $r\simeq 69$ km, consistent with
what found in [49, 51]. At this time, $\dot{q}$ in models with more $Y_{e}$
attenuation that exhibit more flavor conversions of $\nu_{e}$ to $\nu_{\mu}$
shown in Fig. 11 are more negative, with Model V showing the most negative
$\dot{q}\simeq-1.3\times 10^{22}~{}{\rm erg}\,{\rm g}^{-1}\,{\rm s}^{-1}$.

It should also be highlighted that because the prompt flavor conversion takes
place in different timescales among models and at different radii, the heating
rate can be reduced more at an earlier time. For example, the cyan dashed
curve in Fig. 19(a) demonstrates that the net cooling from $r\simeq 28$ km to
$\simeq 46$ km are drastically enhanced at $t=3.2~{}\mu$s in Model IV. The
reduction of $\dot{q}$ at this time is particularly severe around $r\simeq 30$
km and $\simeq 44$ km, where the initial ELN deepness ratios are $\approx 1$
(Fig. 3), allowing nearly complete flavor equilibration or even overconversion
after the prompt FFC. However, for models in which the feedback between FFCs
and collisional weak processes inside neutrinospheres provides supply of
additional $\nu_{e}$’s and $\bar{\nu}_{e}$’s, the reduction of $\dot{q}$ can
be alleviated at later times. This can lead to the change of ordering in terms
of $\dot{q}$ reduction from what was discussed above at $t=32~{}\mu$s. Figure
19(b) shows that the reductions of $\dot{q}$ in Models III and IV are less
than that in Model II around $r\simeq 46$ km at $t=160~{}\mu$s when the
systems approach toward asymptotic states. For Model V, although it remains to
be the model affected by FFCs the most, the impact also becomes slightly
smaller at this time. As for Models VI and VII, because the impact of FFC on
$\nu_{e}$’s and $\bar{\nu}_{e}$’s flavor content are small, $\dot{q}$ in these
two models are only affected insignificantly, and are thus not shown in Fig.
19.

Besides the FFCs, the CFI that is present in Model V also contributes to the
reduction of $\dot{q}$. Figure 19(c) shows that an additional $\dot{q}$
reduction appears at $r\simeq 32$ km and $t=32~{}\mu$s in Model V, which is
otherwise not present in the corresponding model without the CFI (Model Vc0).
The impact of this component also largely diminishes at later times as shown
in Fig. 19(b).

We would like to caution that the impact of FFC on the net heating rate
discussed here is under the assumption that background matter profiles are
static. In realistic situations, even if the prompt flavor conversions can be
reset by the collisional weak processes, its imprint on the matter background
may introduce further feedback to neutrinos and lead to an outcome that sits
in between the state after the prompt conversions and the asymptotic one shown
above. This additional feedback effect may also depend on involved
hydrodynamic timescales and will need to be evaluated in future work.

###  VI.2 Equilibrium electron fraction

Another pertinent implication of FFC on supernova physics is related to the
nucleosynthesis of elements. Although our background matter profile is taken
from a supernova model that does not explode, we expect that the analyses of
the impact of FFC on the chemical composition in our models can still shed
light on what may happen in successfully exploded supernovae with different
progenitor masses as well as the neutrino-driven wind at a later epoch.

For this purpose, we examine the effect of FFC on the equilibrium $Y_{e}$,
which is a key quantity relevant to supernova nucleosynthesis. For $T\lesssim
1$ MeV, the equilibrium $Y_{e}$ is dominated by the neutrino absorption
reactions and can be approximated by

$Y_{e}^{\nu,\mathrm{eq}}=\left(1+\frac{L_{\bar{\nu}_{e}}}{L_{\nu_{e}}}\frac{\langle
E_{\bar{\nu}_{e}}\rangle-2\Delta_{np}+1.2\Delta_{np}^{2}/\langle
E_{\bar{\nu}_{e}}\rangle}{\langle
E_{\bar{\nu}_{e}}\rangle+2\Delta_{np}+1.2\Delta_{np}^{2}/\langle
E_{\nu_{e}}\rangle}\right)^{-1},$ (18)

where $\Delta_{np}=m_{n}-m_{p}$ is the neutron–to-proton mass difference
[139]. Notice that although all neutrinos propagate radially outwards at the
outer boundary of our simulation domain, there remains a fraction of neutrinos
with non-negligible transverse velocities, i.e., $v_{r}\neq 1$. Assuming that
the neutrino flux is conserved and no other types of flavor oscillations occur
outside the outer boundary, the neutrino spectra beyond the simulation domain
for $r>r_{\rm ob}$ can then be approximated by

$\displaystyle\langle\varrho\rangle_{A}(E)=$
$\displaystyle\int_{0}^{1}dv_{r}\,\varrho(E,v_{r},r_{\rm ob})\times$
$\displaystyle\left(\frac{r_{\rm
ob}}{r}\right)^{2}\frac{v_{r}}{\sqrt{1-(\frac{r_{\rm
ob}}{r})^{2}(1-v_{r}^{2})}},$ (19)

which can be further used to calculate the neutrino luminosity and mean energy
relevant to nucleosynthesis.

(a)
(b)

Figure 20:  Comparison of the equilibrium electron fraction $Y_{e}^{\nu,\rm
eq}$ after the prompt FFCs (a) and at a later time $t=160~{}\mu$s (b) for
Models II (blue), III (green), IV (cyan), and V (red). The transparent dashed
curves show the corresponding profiles without flavor conversions for
comparisons.

The equilibrium $Y_{e}^{\nu,\mathrm{eq}}$ in Models II–V are shown from $r=60$
km to 200 km in Fig. 20. The values asymptotically approach constants at
larger radii where neutrinos propagate predominately along the radial
direction. Unlike the heating rate, which is hardly affected by the $Y_{e}$
attenuation for cases without any flavor conversions, the
$Y_{e}^{\nu,\mathrm{eq}}$ profiles depend on the ratio between $\nu_{e}$’s and
$\bar{\nu}_{e}$’s and hence vary from Models II–IV without FFCs (transparent
dashed curves). At $r=200$ km, $Y_{e}^{\nu,\mathrm{eq}}$ are $\approx 0.527$,
0.511, 0.493, and 0.476 in Models II–V without FFCs, respectively, which
reflects the increasing dominance of $\bar{\nu}_{e}$ capture rates that result
in more neutron-rich conditions.

The impact of the prompt FFCs on $Y_{e}^{\nu,{\rm eq}}$ is shown in Fig.
20(a).444Because the timescale for prompt FFCs at $r_{\rm ob}$ varies from
model to model, we show the profiles at $t=96~{}\mu$s for Model II, at
$t=64~{}\mu$s for Model III, and at $t=32~{}\mu$s for Models IV and V,
respectively. Since the prompt FFCs transform an approximately same amount of
$\nu_{e}$ and $\bar{\nu}_{e}$ into the heavy-lepton flavors, it implies that
for cases with an initial asymmetry between $\nu_{e}$ and $\bar{\nu}_{e}$, the
asymmetry gets further enhanced by FFCs. As a result, a $\nu_{e}$-dominated
condition becomes even more dominated by $\nu_{e}$’s, and vice versa. In
addition, given that the mean energy of $\nu_{e}$ is generally lower than
those of other species, FFC enhances $\langle E_{\nu_{e}}\rangle$ relatively
more than $\langle E_{\bar{\nu}_{e}}\rangle$, which also leads to the
increased $Y_{e}^{\nu,\mathrm{eq}}$. Combining these two effects, Fig. 20(a)
shows that $Y_{e}^{\nu,\mathrm{eq}}$ at $r=200$ km are increased by $\approx
0.013$, 0.013, 0.008, and -0.001 in Models II–V, respectively. The impact is
generally larger for models whose $Y_{e}^{\nu,\mathrm{eq}}$ deviate more from
0.5. Note that the impact in Model II is similar to that in Model III, due to
the shallower ELN crossing in Model II that only allows a smaller amount of
FFCs.

Similar to the impact on the heating rate, the impact of flavor conversions on
$Y_{e}^{\nu,\mathrm{eq}}$ at the asymptotic times are reshaped by the feedback
effect from the collisional weak processes in Models III–V. As discussed in
Sec. V.2, the feedback effect tends to bring the flavor content of $\nu_{e}$
back closer to the level prior to oscillations than for $\bar{\nu}_{e}$, this
increases the ratio of $\nu_{e}$ number density over $\bar{\nu}_{e}$ by more
than 10% compared to the ratio before conversions, which leads to a large
change on $Y_{e}$. Figure 20(b) shows that for Model II without feedback
effect, $Y_{e}^{\nu,\mathrm{eq}}$ at $t=160~{}\mu$s is indeed similar to that
shown in Fig. 20(a). In contrast, $Y_{e}^{\nu,\mathrm{eq}}$ in Models III–V
are substantially increased by $\simeq 0.02$–0.04 at the asymptotic times
compared to the values immediately after the FFCs, due to the increased
$\nu_{e}$ to $\bar{\nu}_{e}$ number density ratios. These results suggest that
with the feedback effect, a neutron-rich condition may even be turned into a
proton-rich condition (Model V). Moreover, the feedback effect may further
enhance the proton richness than the outcome of prompt FFCs. All these may
facilitate the production of light $p$-nuclei through the $\nu p$-process
[140, 47].

###  VI.3 Free-streaming neutrino spectra

Finally, we examine the impact of the flavor conversions on the energy spectra
of neutrinos at sufficiently large radii outside the simulation boundary. We
follow the same extrapolation method in Sec. VI.2 to obtain the differential
number luminosity

$\frac{dL_{\rm num}}{dE}(E)=4\pi\int_{0}^{1}dv_{r}\,\varrho(E,v_{r},r_{\rm
ob})r_{\rm ob}^{2}v_{r}$ (20)

at sufficiently large radii. Note that no additional flavor conversions are
included in order to focus on the impact of FFCs on the free-streaming
neutrinos.

Figure 21(a–e) shows the free-streaming neutrino energy spectra at
$t=320~{}\mu$s in Models II–VI for cases with and without flavor conversions.
In addition, Fig. 21(f) shows the spectra after the prompt FFCs at
$t=96~{}\mu$s in Model IV. A common feature is that the degeneracy between
$\nu_{\mu}$’s and $\bar{\nu}_{\mu}$’s are all broken after the flavor
conversions given that $\nu_{e}$’s and $\bar{\nu}_{e}$’s have different energy
spectra prior to flavor conversions.

(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)

Figure 21:  Free-streaming neutrino energy spectra at $t=320~{}\mu$s in Models
II (a), III (b), IV (c), V (d), VI (e), IIf3 (g), IIIf3 (h), and IVf3 (i).
Note that the panel (f) shows the spectrum immediately after the prompt FFC at
$t=96~{}\mu$s in Model IV for illustration.

However, these plots also clearly show that flavor equilibration is not a
general outcome of FFCs. Only the Model V [Fig. 21(d)] as well as the earlier
phase after FFCs in Model III (not shown) and IV [Fig. 21(f)] obtain outcomes
that are close to flavor equilibration. For the rest, equilibration is not
achieved due to the initial shallow ELN angular crossings in Model II, the
dynamical feedback from collisions after the FFCs in Model III and IV, and a
complete reset from the collisional weak processes when the angular crossing
is inside neutrinospheres in Model VI. For cases that include three flavors
shown in Fig. 21(g–i) for Model II–IV, the results are qualitatively similar
to that in the corresponding two-flavor cases.

We also emphasize again that for models where the collisional feedback effect
is important, it results in a net total amount of neutrino luminosity summing
all flavors compared to the values without flavor conversions, as discussed
earlier in Sec. V.2. Even though the flavor equilibration is achieved in Model
V, the energy spectra of all neutrino species lie above the average of the
unoscillated level. Based on this observation, we expect that the FFC can lead
to a significant enhancement for the energy loss rate due to neutrino
emission, which will further affect the CCSN dynamics. This effect, however,
needs to be fully evaluated via a CCSN hydrodynamic simulation that
incorporates the FFCs.

##  VII Discussion and conclusions

We have utilized our multi-group and discrete-ordinate collective neutrino
oscillation simulation code based on cose$\nu$ to comprehensively study the
emergence, evolution, and effects of FFCs in a spherically symmetric supernova
background obtained with agile-boltztran. The collisional weak processes in
this study have been extended from Ref. [71] to include neutrino-electron
scatterings and neutrino pair reactions. We adopted radial dependent $Y_{e}$
attenuation schemes that allow us to probe various scenarios with the ELN
angular crossings residing outside, inside, and across the neutrinosphere, as
well as conditions with neutrino number densities dominated by $\nu_{e}$’s or
$\bar{\nu}_{e}$’s, prior to the FFCs. We have also investigated the impact of
assumed number of neutrino flavors, the choice of the $\mathbf{H}_{\nu\nu}$
attenuation factor, the size of the vacuum term, and the coexistent CFI. All
models examined in this work are listed in Table 2.

To analyze our simulation results, we introduced the spectrogram method that
is suitable for the inhomogeneous neutrino gas in Sec. III.3 to inspect the
small-scale structure of the flavor coherence. The spectrograms are compared
to the dispersion relation obtained from the LSA in Sec. IV. We examined in
detail the evolution history in both two- and three-flavor models as well as
the impact on FFCs from other flavor instabilities in Sec. V. Finally, we
discussed the general effects of FFCs on CCSN physics and the free streaming
neutrino energy spectra in Sec. VI.

We confirmed several characteristics of FFCs reported in a recent study [66].
After the FFCs take place, neutrino profiles eventually reach an asymptotic
state on the coarse-grained level, with the angular crossings in the ELN
distributions eliminated in general. The asymptotic state of neutrinos shows a
strong dependence on the collisional neutrino-matter interactions,
particularly when the ELN angular crossings appear inside the neutrinospheres.

On top of that, we summarize our new findings as follows:

  1. 1.

The comparison of our LSA and spectrogram results shows that the emergence of
FFI and the early evolution of FFCs qualitatively follow the prediction of
dispersion relation from the LSA (Fig. 9). The unstable modes with larger
growth rates develop predominantly, which corresponds to the emergence of
small-scale structure in the off-diagonal flavor coherence. Although the
dominant modes of the small-scale structure may be shifted to slightly
different $K_{r}$ during the later evolution or even disappear for some radial
ranges, they generally survive and stay rather close to the initially most
unstable $K_{r}$ mode from LSA in Models II–IV and Vc0. In addition, we also
observed that the small-scale structure can appear in the diagonal elements,
i.e., neutrino number density, with a length scale consistent with that in the
off-diagonal flavor coherence (Fig. 10).

  2. 2.

We showed that a good scalability of the dispersion relation exists even when
the attenuation on $\mathbf{H}_{\nu\nu}$ is employed (Figs. 4–6). Although
this implies the typical length of the small-scale structure changes with the
attenuation factors, the neutrino flavor evolution has rather weak dependence
on the chosen attenuation factors so long as the dominant $K_{r}$ modes of the
dispersion relation are well resolved (Fig. 17).

  3. 3.

The evolution of the neutrino flavor content exhibits a large variation among
different models (Fig. 11). The prompt flavor conversions immediately after
the emergence of the local FFI can lead to nearly complete flavor
equilibration for models with deep ELN angular crossings (Models III–V).
Models with shallow crossings result in incomplete flavor conversions,
demonstrated by Model II.

  4. 4.

We found that after the prompt FFC phase in Models III–V, an interesting
feedback effect due to the neutrino-matter interaction just below the
neutrinosphere continues to restore the $\nu_{e}$’s and $\bar{\nu}_{e}$’s, and
reshapes the local ELNs as well as the angular distribution of each neutrino
species (Figs. 13 and 14). This effect can largely alter the flavor state of
neutrinos around and outside the neutrinosphere without incurring an
additional ELN angular crossing, leading to an asymptotic state much different
from that obtained right after the prompt FFCs. Although heavy-lepton flavor
neutrinos do not directly experience collisional restoration, their flavor
contents are also enhanced through further flavor conversions from the
restored $\nu_{e}$’s and $\bar{\nu}_{e}$’s. Detailed explanations on this
dynamical feedback were discussed in Sec. V.2.

  5. 5.

When the ELN angular crossings only reside inside the neutrinosphere but not
outside (Models VI and VII), the impact of FFCs can affect the heavy-lepton
flavor neutrinos in the free-streaming regime, as they suffer less collisional
processes on diffusing outwards. However, negligible effects are observed for
$\nu_{e}$’s and $\bar{\nu}_{e}$’s owing to the large opacity that they
experience (Fig. 11). Overall, the asymptotic outcome of FFCs in the free-
streaming regime depends on the how deep the FFI region is inside the
neutrinosphere.

  6. 6.

Regarding the implication of FFCs on the CCSN dynamics, we found that with the
prompt FFCs, it reduces the net heating rate behind the shock due to the
reduced $\nu_{e}$ and $\bar{\nu}_{e}$ absorption rates, resulting in more
cooling below as well as less heating above the gain radius, consistent with
[66]. However, we also found that for Models III–V in which the collisional
feedback plays an important role after the prompt FFCs, the restoration of
$\nu_{e}$’s and $\bar{\nu}_{e}$’s at later times reduces the impact of FFCs on
the enhanced cooling in the asymptotic state (Fig. 19).

  7. 7.

The feedback mechanism also leads to another interesting consequence on
$Y_{e}^{\nu,\rm eq}$. Right after the prompt FFCs, because similar amount of
flavors is converted, the asymmetry between $\nu_{e}$ and $\bar{\nu}_{e}$
number densities is increased. Hence, it enhances the proton or neutron
richness depending on the initial value of $Y_{e}$. However, the collisional
feedback effect restores more $\nu_{e}$’s than $\bar{\nu}_{e}$’s due to the
disparity between their collisional rates. It systematically leads to higher
$Y_{e}^{\nu,\rm eq}$ and makes material more proton-rich in the asymptotic
state (Fig. 20).

  8. 8.

For the emerging neutrino energy spectra after the FFCs relevant to the
neutrino detection, we found that flavor equilibration is not a general
outcome. Significant differences between neutrinos of electron and heavy-
lepton flavors can exist (Fig. 21). This finding indicates that the subsequent
slow flavor instability and Mikheyev–Smirnov–Wolfenstein matter effect can be
relevant and will need to be considered for studying the neutrino signals.

  9. 9.

Although the above conclusions are based on the two-flavor scenario, we found
that the general features of flavor conversions remain qualitatively similar
to that with three flavors. Nevertheless, including a third flavor can result
in additional flavor conversions into the additional sector, which leads to
less amount of $\nu_{e}$ and $\bar{\nu}_{e}$. For the heavy-lepton sectors,
the enhanced amounts for muon and tauon flavors in the three-flavor scenarios
are less than that for muon flavors in the two-flavor cases. However, the
total amount of heavy-lepton flavors in three-flavor cases is in fact larger
than that in the two-flavor cases (Fig. 16).

  10. 10.

We showed that the dispersion relation in our simulation domain is barely
affected by the inclusion of the vacuum term when no attenuation is applied to
$\mathbf{H}_{\nu\nu}$. When attenuation is applied, one may pragmatically
adopt a reduced $\delta m^{2}$ to maintain the scalability of the dispersion
relation. In this case, the simulation results are largely insensitive to the
presence of the vacuum term (Fig. 18).

  11. 11.

We found that the occurrence and the strength of CFI, which is intrinsically a
multienergy phenomenon due to the energy dependence of the EA collisional
rates, can depend on the energy spectra of heavy-lepton neutrinos. When the
NES processes are included, the condition for the CFI in this work becomes
different from that in Ref. [71] and can be understood using the formula
developed in Ref. [73] (Sec. IV.2). Among all models considered in this paper,
noticeable flavor conversion due to the CFI is only observed inside the FFC
region in Model V when the initial condition is dominated by
$\bar{\nu}_{e}$’s. While this affects less the properties of neutrinos than
the FFC does, it can lead to a transient enhancement of the net cooling rate
in the corresponding region (Fig. 19).

In summary, our comprehensive study of neutrino flavor conversions in the
post-shocked region supports that it is important to include FFCs in supernova
simulations, in order to fully understand the role of neutrinos in supernova
dynamics, nucleosynthesis, as well as to predict the CCSN neutrino signals.
Our work also highlights the need of using energy-dependent and realistic
neutrino collision rates in relevant simulations, because the condition and
evolution of FFC and CFI can be affected by those collision processes either
directly or through feedback. Certain assumptions that were introduced in our
work such as the spherical symmetry, static supernova background, and the
missing coherent forward scattering with the matter remain to be further
examined and will be studied in future.

###### Acknowledgements.

We thank Ninoy Rahman and Gang Guo for discussions on collisional weak
processes, Yi-Siou Wu on the comparison of two- and three-flavor results, and
Sherwood Richers on the effective method of treating neutrino pair reactions.
ZX and GMP acknowledge support of the European Research Council (ERC) under
the European Union’s Horizon 2020 research and innovation programme (ERC
Advanced Grant KILONOVA No. 885281), the Deutsche Forschungsgemeinschaft (DFG,
German Research Foundation) – Project-ID 279384907 – SFB 1245, and MA
4248/3-1. MRW and MG acknowledge supports from the National Science and
Techonology Council, Taiwan under Grant No. 111-2628-M-001-003-MY4, and the
Academia Sinica (Project No. AS-CDA-109-M11). MRW also acknowledges supports
from the Physics Division of the National Center for Theoretical Sciences,
Taiwan. CYL acknowledges support from the National Center for High-performance
Computing (NCHC). TF and NKL acknowledge support from the Polish National
Science Center (NCN) under Grant No. 2020/37/B/ST9/00691. The core-collapse
supernova simulations were performed at the Wroclaw Center for Scientific
Computing and Networking (WCSS). The neutrino quantum kinetic simulations were
performed at the GSI Helmholtzzentrum für Schwerionenforschung. We acknowledge
the following software: matplotlib [141], numpy [142], and scipy [143].

## References

  * Mirizzi  _et al._ [2016] A. Mirizzi, I. Tamborra, H.-T. Janka, N. Saviano, K. Scholberg, R. Bollig, L. Hüdepohl, and S. Chakraborty, Nuovo Cimento Rivista Serie 39, 1 (2016), arXiv:1508.00785 [astro-ph.HE] .
  * Fischer  _et al._ [2023] T. Fischer, G. Guo, K. Langanke, G. Martinez-Pinedo, Y.-Z. Qian, and M.-R. Wu, (2023), arXiv:2308.03962 [astro-ph.HE] .
  * Workman  _et al._ [2022] R. L. Workman, V. D. Burkert, V. Crede, E. Klempt, and et. al., Prog. Theor. Exp. Phys. 2022, 083C01 (2022).
  * Tamborra and Shalgar [2021] I. Tamborra and S. Shalgar, Annu. Rev. Nucl. Part. Sci. 71, 165 (2021), arXiv:2011.01948 .
  * Richers and Sen [2023] S. Richers and M. Sen, in  _Handb. Nucl. Phys._ (Springer, 2023) pp. 3771–3787, arXiv:2207.03561 .
  * Capozzi and Saviano [2022] F. Capozzi and N. Saviano, Universe 8, 94 (2022), arXiv:2202.02494 .
  * Volpe [2023] M. C. Volpe, arXiv e-prints (2023), arXiv:2301.11814 .
  * Abbar  _et al._ [2019] S. Abbar, H. Duan, K. Sumiyoshi, T. Takiwaki, and M. C. Volpe, Phys. Rev. D 100, 043004 (2019), arXiv:1812.06883 [astro-ph.HE] .
  * Delfan Azari  _et al._ [2019] M. Delfan Azari, S. Yamada, T. Morinaga, W. Iwakami, H. Okawa, H. Nagakura, and K. Sumiyoshi, Phys. Rev. D 99, 103011 (2019), arXiv:1902.07467 .
  * Delfan Azari  _et al._ [2020] M. Delfan Azari, S. Yamada, T. Morinaga, H. Nagakura, S. Furusawa, A. Harada, H. Okawa, W. Iwakami, and K. Sumiyoshi, Phys. Rev. D 101, 023018 (2020), arXiv:1910.06176 [astro-ph.HE] .
  * Nagakura  _et al._ [2019] H. Nagakura, T. Morinaga, C. Kato, and S. Yamada, Astrophys. J. 886, 139 (2019), arXiv:1910.04288 .
  * Morinaga  _et al._ [2020] T. Morinaga, H. Nagakura, C. Kato, and S. Yamada, Phys. Rev. Res. 2, 012046(R) (2020), arXiv:1909.13131 [astro-ph.HE] .
  * Abbar  _et al._ [2020] S. Abbar, H. Duan, K. Sumiyoshi, T. Takiwaki, and M. C. Volpe, Phys. Rev. D 101, 043016 (2020), arXiv:1911.01983 .
  * Glas  _et al._ [2020] R. Glas, H.-T. Janka, F. Capozzi, M. Sen, B. Dasgupta, A. Mirizzi, and G. Sigl, Phys. Rev. D 101, 063001 (2020), arXiv:1912.00274 .
  * Nagakura  _et al._ [2021] H. Nagakura, A. Burrows, L. Johns, and G. M. Fuller, Phys. Rev. D 104, 083025 (2021), arXiv:2108.07281 [astro-ph.HE] .
  * Abbar  _et al._ [2021] S. Abbar, F. Capozzi, R. Glas, H. T. Janka, and I. Tamborra, Phys. Rev. D 103, 063033 (2021), arXiv:2012.06594 .
  * Harada and Nagakura [2022] A. Harada and H. Nagakura, Astrophys. J. 924, 109 (2022), arXiv:2110.08291 [astro-ph.HE] .
  * Sawyer [2009] R. F. Sawyer, Phys. Rev. D 79, 105003 (2009), arXiv:arXiv:0803.4319 .
  * Sigl and Raffelt [1993] G. Sigl and G. Raffelt, Nuclear Physics B 406, 423 (1993).
  * Vlasenko  _et al._ [2014] A. Vlasenko, G. M. Fuller, and V. Cirigliano, Phys. Rev. D 89, 105004 (2014), arXiv:1309.2628 .
  * Volpe [2015] C. Volpe, Int. J. Mod. Phys. E 24, 1541009 (2015), arXiv:1506.06222 .
  * Blaschke and Cirigliano [2016] D. N. Blaschke and V. Cirigliano, Phys. Rev. D 94, 033009 (2016), arXiv:1605.09383 [hep-ph] .
  * Banerjee  _et al._ [2011] A. Banerjee, A. Dighe, and G. Raffelt, Phys. Rev. D 84, 053013 (2011), arXiv:1107.2308 .
  * Izaguirre  _et al._ [2017] I. Izaguirre, G. Raffelt, and I. Tamborra, Phys. Rev. Lett. 118, 021101 (2017), arXiv:1610.01612 [hep-ph] .
  * Capozzi  _et al._ [2017] F. Capozzi, B. Dasgupta, E. Lisi, A. Marrone, and A. Mirizzi, Phys. Rev. D 96, 043016 (2017), arXiv:1706.03360 .
  * Yi  _et al._ [2019] C. Yi, L. Ma, J. D. Martin, and H. Duan, Phys. Rev. D 99, 063005 (2019), arXiv:1901.01546 [hep-ph] .
  * Morinaga [2022] T. Morinaga, Phys. Rev. D 105, L101301 (2022), arXiv:2103.15267 .
  * Martin  _et al._ [2020] J. D. Martin, C. Yi, and H. Duan, Phys. Lett. B 800, 135088 (2020), arXiv:1909.05225 [hep-ph] .
  * Bhattacharyya and Dasgupta [2021] S. Bhattacharyya and B. Dasgupta, Phys. Rev. Lett. 126, 061302 (2021), arXiv:2009.03337 [hep-ph] .
  * Bhattacharyya and Dasgupta [2020] S. Bhattacharyya and B. Dasgupta, Phys. Rev. D 102, 063018 (2020), arXiv:2005.00459 .
  * Wu  _et al._ [2021] M.-R. Wu, M. George, C.-Y. Lin, and Z. Xiong, Phys. Rev. D 104, 103003 (2021), arXiv:2108.09886 [hep-ph] .
  * Richers  _et al._ [2021a] S. Richers, D. E. Willcox, N. M. Ford, and A. Myers, Phys. Rev. D 103, 083013 (2021a), arXiv:2101.02745 .
  * Richers  _et al._ [2021b] S. Richers, D. Willcox, and N. Ford, Phys. Rev. D 104, 103023 (2021b), arXiv:2109.08631 .
  * Zaizen and Morinaga [2021] M. Zaizen and T. Morinaga, Phys. Rev. D 104, 083035 (2021), arXiv:2104.10532 [hep-ph] .
  * Abbar and Capozzi [2022] S. Abbar and F. Capozzi, J. Cosmol. Astropart. Phys. 2022, 051 (2022), arXiv:2111.14880 [astro-ph.HE] .
  * Richers  _et al._ [2022] S. Richers, H. Duan, M.-R. Wu, S. Bhattacharyya, M. Zaizen, M. George, C.-Y. Lin, and Z. Xiong, Phys. Rev. D 106, 043011 (2022), arXiv:2205.06282 .
  * Bhattacharyya and Dasgupta [2022] S. Bhattacharyya and B. Dasgupta, Phys. Rev. D 106, 103039 (2022), arXiv:2205.05129 .
  * Grohs  _et al._ [2023a] E. Grohs, S. Richers, S. M. Couch, F. Foucart, J. P. Kneller, and G. C. McLaughlin, Phys. Lett. B 846, 138210 (2023a), arXiv:2207.02214 .
  * Zaizen and Nagakura [2023a] M. Zaizen and H. Nagakura, Phys. Rev. D 107, 103022 (2023a), arXiv:2211.09343 [astro-ph.HE] .
  * Zaizen and Nagakura [2023b] M. Zaizen and H. Nagakura, Phys. Rev. D 107, 123021 (2023b), arXiv:2304.05044 [astro-ph.HE] .
  * Xiong  _et al._ [2023a] Z. Xiong, M.-R. Wu, S. Abbar, S. Bhattacharyya, M. George, and C.-Y. Lin, Phys. Rev. D 108, 063003 (2023a).
  * Nagakura  _et al._ [2023] H. Nagakura, L. Johns, and M. Zaizen, (2023), arXiv:2312.16285v1 .
  * Raffelt and Sigl [2007] G. G. Raffelt and G. Sigl, Phys. Rev. D 75, 083002 (2007), arXiv:hep-ph/0701182 [hep-ph] .
  * Abbar and Volpe [2019] S. Abbar and M. C. Volpe, Phys. Lett. B 790, 545 (2019), arXiv:1811.04215 .
  * Johns  _et al._ [2020] L. Johns, H. Nagakura, G. M. Fuller, and A. Burrows, Phys. Rev. D 102, 103017 (2020).
  * Xiong  _et al._ [2023b] Z. Xiong, M.-R. Wu, and Y.-Z. Qian, Phys. Rev. D 108, 043007 (2023b), arXiv:2303.05906 .
  * Xiong  _et al._ [2020] Z. Xiong, A. Sieverding, M. Sen, and Y.-Z. Qian, Astrophys. J. 900, 144 (2020), arXiv:2006.11414 [astro-ph.HE] .
  * Fujimoto and Nagakura [2023] S.-I. Fujimoto and H. Nagakura, Mon. Not. R. Astron. Soc. 519, 2623 (2023), arXiv:2210.02106 .
  * Nagakura [2023] H. Nagakura, Phys. Rev. Lett. 130, 211401 (2023), arXiv:2301.10785 [astro-ph.HE] .
  * Ehring  _et al._ [2023a] J. Ehring, S. Abbar, H. T. Janka, G. Raffelt, and I. Tamborra, Phys. Rev. D 107, 103034 (2023a), arXiv:2301.11938 .
  * Ehring  _et al._ [2023b] J. Ehring, S. Abbar, H. T. Janka, G. Raffelt, and I. Tamborra, Phys. Rev. Lett. 131, 061401 (2023b), arXiv:2305.11207 .
  * Balantekin  _et al._ [2023] A. B. Balantekin, M. J. Cervia, A. V. Patwardhan, R. Surman, and X. Wang, (2023), arXiv:2311.02562 .
  * Wu  _et al._ [2017] M.-R. Wu, I. Tamborra, O. Just, and H.-T. Janka, Phys. Rev. D 96, 123015 (2017), arXiv:1711.00477 [astro-ph.HE] .
  * George  _et al._ [2020] M. George, M.-R. Wu, I. Tamborra, R. Ardevol-Pulpillo, and H.-T. Janka, Phys. Rev. D 102, 103015 (2020), arXiv:2009.04046 [astro-ph.HE] .
  * Li and Siegel [2021] X. Li and D. M. Siegel, Phys. Rev. Lett. 126, 251101 (2021), arXiv:2103.02616 [astro-ph.HE] .
  * Just  _et al._ [2022] O. Just, S. Abbar, M. R. Wu, I. Tamborra, H. T. Janka, and F. Capozzi, Phys. Rev. D 105, 083024 (2022), arXiv:2203.16559 .
  * Fernández  _et al._ [2022] R. Fernández, S. Richers, N. Mulyk, and S. Fahlman, Phys. Rev. D 106, 103003 (2022), arXiv:2207.10680 [astro-ph.HE] .
  * Nagakura [2023] H. Nagakura, Phys. Rev. D 108, 103014 (2023), arXiv:2306.10108 .
  * Froustey  _et al._ [2023] J. Froustey, S. Richers, E. Grohs, S. Flynn, F. Foucart, J. P. Kneller, and G. C. McLaughlin, (2023), arXiv:2311.11968 .
  * Grohs  _et al._ [2023b] E. Grohs, S. Richers, S. M. Couch, F. Foucart, J. Froustey, J. Kneller, and G. McLaughlin, (2023b), arXiv:2309.00972 .
  * Stapleford  _et al._ [2020] C. J. Stapleford, C. Fröhlich, and J. P. Kneller, Phys. Rev. D 102, 081301(R) (2020), arXiv:1910.04172 [astro-ph.HE] .
  * Shalgar and Tamborra [2023] S. Shalgar and I. Tamborra, Phys. Rev. D 108, 043006 (2023), arXiv:2206.00676 .
  * Nagakura and Zaizen [2022] H. Nagakura and M. Zaizen, Phys. Rev. Lett. 129, 261101 (2022), arXiv:2206.04097v1 .
  * Shalgar and Tamborra [2023] S. Shalgar and I. Tamborra, Phys. Rev. D 107, 063025 (2023), arXiv:2207.04058 [astro-ph.HE] .
  * Nagakura and Zaizen [2023] H. Nagakura and M. Zaizen, Phys. Rev. D 107, 063033 (2023), arXiv:2211.01398 [astro-ph.HE] .
  * Nagakura and Zaizen [2023] H. Nagakura and M. Zaizen, (2023), arXiv:2308.14800 .
  * Duan  _et al._ [2010] H. Duan, G. M. Fuller, and Y.-Z. Qian, Annu. Rev. Nucl. Part. 60, 569 (2010), arXiv:1001.2799 [hep-ph] .
  * Johns [2023a] L. Johns, Phys. Rev. Lett. 130, 191001 (2023a), arXiv:2104.11369 .
  * Padilla-Gay  _et al._ [2022] I. Padilla-Gay, I. Tamborra, and G. G. Raffelt, Phys. Rev. D 106, 103031 (2022), arXiv:2209.11235 .
  * Johns and Xiong [2022] L. Johns and Z. Xiong, Phys. Rev. D 106, 103029 (2022), arXiv:2208.11059 [hep-ph] .
  * Xiong  _et al._ [2023] Z. Xiong, M.-R. Wu, G. Martínez-Pinedo, T. Fischer, M. George, C.-Y. Lin, and L. Johns, Phys. Rev. D 107, 083016 (2023), arXiv:2210.08254 [astro-ph.HE] .
  * Lin and Duan [2023] Y.-C. Lin and H. Duan, Phys. Rev. D 107, 083034 (2023), arXiv:2210.09218 .
  * Xiong  _et al._ [2023] Z. Xiong, L. Johns, M.-R. Wu, and H. Duan, Phys. Rev. D 108, 083002 (2023), arXiv:2212.03750 .
  * Kato  _et al._ [2023] C. Kato, H. Nagakura, and M. Zaizen, Phys. Rev. D 108, 023006 (2023), arXiv:2303.16453 [astro-ph.HE] .
  * Liu  _et al._ [2023a] J. Liu, M. Zaizen, and S. Yamada, Phys. Rev. D 107, 123011 (2023a), arXiv:2302.06263 .
  * Akaho  _et al._ [2024] R. Akaho, J. Liu, H. Nagakura, M. Zaizen, and S. Yamada,  109, 023012 (2024), arXiv:2311.11272 .
  * Liu  _et al._ [2023b] J. Liu, R. Akaho, A. Ito, H. Nagakura, M. Zaizen, and S. Yamada, Phys. Rev. D 108, 123024 (2023b), arXiv:2310.05050 .
  * Kato  _et al._ [2023] C. Kato, H. Nagakura, and L. Johns, (2023), arXiv:2309.02619 .
  * Fiorillo  _et al._ [2023] D. F. G. Fiorillo, I. Padilla-Gay, and G. G. Raffelt, (2023), arXiv:2312.07612 .
  * Shalgar and Tamborra [2023] S. Shalgar and I. Tamborra, (2023), arXiv:2307.10366 .
  * George  _et al._ [2023] M. George, C.-Y. Lin, M.-R. Wu, T. G. Liu, and Z. Xiong, Comput. Phys. Commun. 283, 108588 (2023), arXiv:2203.12866 .
  * Tamborra  _et al._ [2017] I. Tamborra, L. Hüdepohl, G. G. Raffelt, and H.-T. Janka, Astrophys. J. 839, 132 (2017), arXiv:1702.00060 [astro-ph.HE] .
  * Nagakura  _et al._ [2020] H. Nagakura, A. Burrows, D. Radice, and D. Vartanyan, Mon. Not. R. Astron. Soc. 492, 5764 (2020), arXiv:1912.07615 .
  * Tamborra  _et al._ [2014] I. Tamborra, F. Hanke, H. T. Janka, B. Müller, G. G. Raffelt, and A. Marek, Astrophys. J. 792, 96 (2014), arXiv:1402.5418 .
  * Glas  _et al._ [2019] R. Glas, H.-T. Janka, T. Melson, G. Stockinger, and O. Just, Astrophys. J. 881, 36 (2019), arXiv:1809.10150 .
  * Rauscher  _et al._ [2002] T. Rauscher, A. Heger, R. D. Hoffman, and S. E. Woosley, Astrophys. J. 576, 323 (2002), arXiv:0112478 [astro-ph] .
  * Mezzacappa and Bruenn [1993a] A. Mezzacappa and S. W. Bruenn, Astrophys. J. 405, 637 (1993a).
  * Mezzacappa and Bruenn [1993b] A. Mezzacappa and S. W. Bruenn, Astrophys. J. 410, 740 (1993b).
  * Mezzacappa and Bruenn [1993c] A. Mezzacappa and S. W. Bruenn, Astrophys. J. 405, 669 (1993c).
  * Liebendörfer  _et al._ [2001] M. Liebendörfer, A. Mezzacappa, and F.-K. Thielemann, Phys. Rev. D 63, 104003 (2001), arXiv:astro-ph/0012201 [astro-ph] .
  * Liebendörfer  _et al._ [2004] M. Liebendörfer, O. E. B. Messer, A. Mezzacappa, S. W. Bruenn, C. Y. Cardall, and F. K. Thielemann, Astrophys. J., Suppl. 150, 263 (2004), arXiv:astro-ph/0207036 [astro-ph] .
  * Lindquist [1966] R. W. Lindquist, Ann. Physics1 37, 487 (1966).
  * Fischer  _et al._ [2020] T. Fischer, G. Guo, G. Martínez-Pinedo, M. Liebendörfer, and A. Mezzacappa, Phys. Rev. D 102, 123001 (2020), arXiv:2008.13628 [astro-ph.HE] .