Detectability of Late-time Supernova Neutrinos with Fallback Accretion onto Protoneutron Star

We investigate the late-time neutrino emission powered by fallback mass accretion onto a protoneutron star (PNS), using neutrino radiation-hydrodynamic simulations with full Boltzmann neutrino transport. We follow the time evolution of the accretion flow onto the PNS until the system reaches a quasi-steady state. A standing shock wave is commonly formed in the accretion flow, whereas the shock radius varies depending on the mass accretion rate and the PNS mass. A sharp increase in temperature emerges in the vicinity of the PNS (∼10 km), which characterizes neutrino emission. Both the neutrino luminosity and the average energy become higher with increasing mass accretion rate and PNS mass. The mean energy of the emitted neutrinos is in the range of 10 ≲ ϵ ≲ 20 MeV, which is higher than that estimated from PNS cooling models (≲10 MeV). Assuming a distance to core-collapse supernova of 10 kpc, we quantify neutrino event rates for Super-Kamiokande (Super-K) and Deep Underground Neutrino Experiment (DUNE). The estimated detection rates are well above the background, and their energy-dependent features are qualitatively different from those expected from PNS cooling models. Another notable feature is that the neutrino emission is strongly flavor dependent, exhibiting that the neutrino event rate hinges on the neutrino oscillation model. We estimate them in the case with the adiabatic Mikheev–Smirnov–Wolfenstein model, and show that the normal and inverted mass hierarchy offer a large number of neutrino detections in Super-K and DUNE, respectively. Hence the simultaneous observation with Super-K and DUNE of fallback neutrinos will provide a strong constraint on the neutrino mass hierarchy.


INTRODUCTION
Most massive stars with zero-age main sequence mass ≳ 8M ⊙ end their lives as core-collapse supernovae (CC-SNe).In the central region, a protoneutron star (PNS) is formed as a consequence of gravitational collapse of its iron core.A huge gravitational energy of the PNS is released via neutrinos.The neutrinos absorbed behind the stalled shock wave can foster shock expansion.
After shock revival, the evolution of the PNS proceeds on a Kelvin-Helmholtz timescale of neutrino cooling (of the order of ten seconds) (Burrows & Lattimer 1986), that is much longer than the timescale on which shock revival is initiated by the neutrino-heating mechanism (≲ 1s) (see for a recent review, e.g., Burrows & Vartanyan 2021).The neutrino signal thus carries imprints of not only the explosion mechanism but also the subsequent evolution of the remnant system.This argument is in line with neutrino data of SN 1987A (Hirata et al. 1987;Bionta et al. 1987).On the other hand, the low-statistics neutrino data from SN 1987A provided lit-tle constraint on the detailed features of neutrinos such as flavor-dependent features and their time structure.With current neutrino detectors the next nearby CCSN will provide high-statistics signals placing constraints on the CCSN dynamics.
Significant progress has been made in the last decades in neutrino detection techniques.Various types of neutrino detectors such as water Cherenkov (Ikeda et al. 2007;Hyper-Kamiokande Proto-Collaboration et al. 2018;Abbasi et al. 2011), liquid argon (Abi et al. 2021), liquid scintillators (Asakura et al. 2016;An et al. 2016), and dark-matter detectors (Lang et al. 2016) are currently operating or planned, offering the means to distinguish neutrino flavors (Horiuchi & Kneller 2018).The sensitivity of detectors has steadily improved, and the detector size has also increased by more than an order of magnitude compared to those used in the 1980s, suggesting that the long-term neutrino signal will be detectable for nearby CCSNe (see, e.g., Suwa et al. 2019;Li et al. 2021).
This has motivated the CCSN community to develop theoretical models for long-term evolution of PNS cooling.Numerical simulations have also been performed (see, e.g., Fischer et al. 2010;Roberts et al. 2012;Roberts 2012), and various physical quantities relevant to PNS cooling have also been investigated: the nuclear matter equation of state (Nakazato & Suzuki 2019;Nakazato & Suzuki 2020;Nakazato et al. 2022;Sumiyoshi et al. 2023), progenitor dependence (Nakazato et al. 2013), and neutrino matter interactions (Fischer et al. 2012;Martínez-Pinedo et al. 2012;Fischer et al. 2020;Pascal et al. 2022;Sugiura et al. 2022).These works will play pivotal roles to place a constraint on microphysical parameters in supranuclear density of PNS in real observations.
In this paper, we discuss the late time neutrino emission of CCSN from a different perspective: fallback accretion (FBA) onto the PNS.Most previous studies have a priori assumed that FBA has no influence on the neutrino signal.One thing we do notice here is, however, that large amounts of FBA have been observed rather commonly in recent multi-dimensional (multi-D) CCSN simulations (see, e.g., Burrows & Vartanyan 2021;Bollig et al. 2021;Nagakura et al. 2021a).More interestingly, they may last a very long time (≫ 10s) (see, e.g., Fig. 2 in Janka et al. 2022) due to the shock deceleration or reverse shock that occurs after the shock wave passes the CO/He-core interface (Fryxell et al. 1991) and He/H interface (Chevalier 1989).This suggests that the neutrino emission from FBA potentially overwhelms those radiated from PNS.
The impact of FBA in the late time neutrino emission was investigated by the pioneering work of Fryer (2009).This study showed that FBA has a large influence on the neutrino luminosity and their average energy.It should be noted, however, that there are potential systematic uncertainties in their models; for instances, the inner boundary of the computational domain in the simulations is located much outside the neutrino sphere (which will be shown in Sec.4.2 and see also Table 2 in Fryer (2009)), and the neutrino transport was handled with a gray flux-limited diffusion approximation (Herant et al. 1994).These simplifications prevented them from studying detailed features of neutrinos from FBA, and they may discard some important properties inherent in FBA.
In this paper, we investigate the neutrino emission driven by FBA onto the PNS by performing neutrino radiation-hydrodynamics simulations covering the optically thick region in the computational domain.We analyze the neutrino emission by systematically changing the accretion rate (at the outer boundary) and the PNS mass.Based on the numerical simulations, we provide some key features of the neutrino signal from FBA, and then the neutrino event rates in some representative CCSN neutrino detectors are estimated.
This paper is organized as follows.In Sec. 2, we start with providing an overview of FBA onto PNS in the post explosion phase, and then we describe our approach to study the neutrino signal powered by FBA.The numerical methods and the PNS models are described in Sec. 3.All results of our numerical simulations are encapsulated in Sec. 4. The detectability of neutrinos is discussed in Sec. 5. We summarize our findings in section 6.

FALLBACK ACCRETION IN CCSNE
Even after the shock wave begins its runaway expansion, a certain amount of post-shock matter is bound by the gravity of the PNS, and it eventually returns back to the PNS.Such FBA in CCSNe has been studied in the literature from the early 1970s.The importance of FBA was first pointed out by Colgate (1971).They suggested that FBA is necessary to explain the consistent amount of nucleosynthetic yields.From the observational point of view, some previous studies suggested that FBA has an influence on both electromagnetic- (Dexter & Kasen 2013) and neutrino emission (Fryer 2009) in the late phase.We also note that FBA potentially accounts for some peculiar energetic (Moriya et al. 2018(Moriya et al. , 2019) ) and weak CCSN explosions (Moriya et al. 2010).If FBA leads to an oversupply of mass onto the PNS, it may trigger a black hole formation (Zhang et al. 2008;Chan et al. 2018).If the core is rapidly rotating, gammaray burst would occur following the collapsar scenario (MacFadyen & Woosley 1999;MacFadyen et al. 2001;Perna et al. 2014).FBA can also affect the PNS spin (Barrère et al. 2022;Ronchi et al. 2022;Coleman & Burrows 2022) and its spin-kick alignment indicated by some pulsar observations (Johnston et al. 2005(Johnston et al. , 2007;;Ng & Romani 2007;Janka et al. 2022).
FBA in CCSNe can be categorized into several phases (Chevalier 1989).In the early post shock revival phase, it would be chaotic due to the turbulent accretion flows originated from multi-D fluid instabilities in the postshock region.It should also be mentioned that asymmetric shock revival can lead to large FBA from the angular region where the shock expansion is weaker (Nagakura et al. 2021a;Bollig et al. 2021).In the very later phase, which is referred to as the uniform expansion phase ≳ 10 3 s, the accretion rate simply scales as Ṁ ∝ t −5/3 (Chevalier 1989).This scaling is verified by various numerical studies (Zhang et al. 2008;Dexter & Kasen 2013;Janka et al. 2022).Note that the accretion rate on this free-expansion phase may be significantly enhanced by the arrival of the reverse shock onto PNS.It has also been suggested that strong FBA can be formed by the deceleration of the shock at the CO/He-core or He/H interfaces (Janka et al. 2022) (see also Fig. 2 of Zhang et al. (2008), in which the enhancement of FBA is clearly visible).
In this paper, we pay attention to the phase of ≳ 10 s after core bounce.In this phase, the PNS temperature at the surface becomes less than ∼ 3 MeV (Roberts 2012;Nakazato et al. 2013), and the neutrino emission gradually subsides in the Kelvin-Helmholtz timescale.This suggests that the neutrino emission can be dominated by FBA, inferred from the previous works (Fryer 2009;Nagakura et al. 2020Nagakura et al. , 2021a;;Bollig et al. 2021).It should also be noted that we develop a general discussion of FBA neutrino emission without specifying any late phases in this study, since our approach can be applied to different situations.Nevertheless, the increase of FBA by a reverse shock created at the CO/He-core or He/H interfaces is an intriguing phase, since a large FBA may happen at a very late phase of CCSNe (≳ 10 3 s) (Zhang et al. 2008).

Numerical Method
We perform neutrino radiation-hydrodynamics simulation of FBA onto the PNS.We employ a general relativistic Boltzmann radiation-hydrodynamics code, in which we solve the hydrodynamic equation and the Boltzmann equation for neutrino transport, both in general relativity.The numerical details can be found in a series of papers (Nagakura et al. 2014(Nagakura et al. , 2017(Nagakura et al. , 2019;;Akaho et al. 2021;Akaho et al. 2023).Although our code can treat multi-dimension in space, we assume spherically symmetry in this study.
The Boltzmann neutrino transport directly solves the neutrino distribution function under multi-species, multi-energy, and multi-angle treatments.This allows us to develop accurate models of the neutrino radiation field even in semi-transparent regions, whereas the accuracy is not guaranteed in other approximate methods (e.g., MGFLD and two-moment methods).The interesting issue of quantifying the error for each approximate transport method is beyond the scope of this paper.
Assuming spherical symmetry, the distribution function can be expressed with four variables; time (t), radius (r), neutrino energy (ϵ), and the zenith angle in the momentum space (θ ν ).We note that the spacetime metric is determined by solving the Tolman-Oppenheimer-Volkov (TOV) equations used to construct the NS, as we explain later in section 3.2.The explicit form of the Boltzmann equation in the conserved form can be writ-ten as (Shibata et al. 2014) where α, and γ ij denote the lapse function and the spatial metric, respectively.The factors ω (t) and ω (θν ) are defined as where p µ denotes the neutrino four-momentum.The tetrad bases are given as Although the FBA is a priori multi-D, the accretion energy converts to thermal energy in the vicinity of PNS, and eventually spreads all over the PNS surface.This suggests that the asymmetry of neutrino emission becomes milder than that of FBA, and numerical simulations also support this assumption (Vartanyan et al. 2019).Whether large asymmetries of neutrino emission can be created by non-radial FBA is an interesting question which we defer to future work.We also note that the spherically symmetric conditions artificially suppress the PNS convection.One may wonder if this may cause to underestimate the diffusion component of the neutrino luminosity.According to recent multi-D simulations, however, the PNS convection subsides by ∼ 5s after bounce (Nagakura et al. 2021a), suggesting that it does not affect FBA neutrinos in the late phase (t ≳ 10s).
We employ 512 radial grid points covering the range r ∈ [0 : 100] km.We note that the resolution is high around the PNS surface (where the minimum mesh width is ∼ 30 m).Such a high spatial resolution is mandatory in studying FBA, since the scale heights of matter-and neutrino-radiation field around the PNS are very small.The number of energy mesh is 20, covering the range ϵ ∈ [0 : 300] MeV and the mesh width is logarithmically distributed.The zenith angle θ ν ∈ [0 : π] has 10 grid points.
We employ the Furusawa-Togashi EOS (Furusawa et al. 2017) with some extension.It should be mentioned that, in the case with low mass accretion rate, the thermodynamical quantities can be outside of the range covered by the EOS table.To deal with this issue, we extended the EOS table in a pragmatic way; it is smoothly connected to the gamma-law EOS as the pressure given as P = (Γ − 1)ρϵ, where ρ and ϵ denote the density and the specific internal energy, respectively.The gamma law index Γ is obtained from the edge of the EOS table.We found that Γ is almost 4/3 for various input parameters.We note, however, that our prescription is rather pragmatic, and it does not have the ability to capture the realistic matter evolution; in particular for shock dynamics.For this reason, we stop the calculation if the shock wave reaches the position where thermodynamical quantities are out of the range of the Furusawa-Togashi EOS.On the other hand, these prescriptions do not compromise the present result, since neutrino emission occurs in high density regions, which are always covered by the Furusawa-Togashi EOS.
Regarding neutrino-matter interactions, we incorporate them based on the standard set (Bruenn 1985) with some extensions: energy-changing neutrino-electron scattering and the nucleon-nucleon bremsstrahlung.See Sumiyoshi et al. (2005); Nagakura et al. (2017) for the numerical implementations of these reactions.

Models
The accurate determination of the neutrino emission from FBA requires resolving the PNS surface where the accretion energy is converted to thermal energy.We also note that the neutrino opacity hinges on the energy, and the low energy neutrinos can escape from the very high density region (> 10 14 g/cm 3 ), exhibiting that the PNS structure needs to be determined to quantify the neutrino energy spectrum.Hence, we construct the PNS structure by assuming steady state, before running FBA simulations.
We assume an isotropic temperature of T = 2 MeV and the electron fraction of Y e = 0.05 inside the PNS as a reference model, which represents the matter state of the PNS in the cooling phase.For the sake of completeness, the temperature dependence in the neutrino signal is also checked in this study (see Sec. 5.3).Given T and Y e , we prepare two different PNS structures by solving TOV equations, by varying the central baryon mass density.The first one has the central density of ρ c = 8.5 × 10 14 g • cm −3 , leading to the total mass M PNS = 1.41M ⊙ .For the second one, we set ρ c = 1.2 × 10 15 g • cm −3 , that leads to M PNS = 1.98M ⊙ .The spacetime metric obtained from solving the TOV equations is used for the radiation-hydrodynamic simulations, and kept fixed in time.
Below, we describe our FBA model.One thing to note here is that the self-consistent treatment of FBA requires successful CCSN explosion models.It should be noted, however, that detailed features of FBA such as the mass accretion rate, the thermodynamical states, and their time evolution strongly depend on the progenitor, the timing of the shock revival, and the ejecta morphology.In this study, we are not interested in such details of the neutrino signal, but rather in generic features that can be applied to any types of FBA.To this end, we treat FBA in a simple manner capturing the essential features.In our models, we assume a mass inflow from the outer boundary of the computational domain.
The accretion rate is one of the control parameters, and we study four cases: Ṁ = 10 −2 , 5 × 10 −3 , 2 × 10 −3 , and  10 The choice of mass accretion rate is motivated by previous studies of FBA by (Chan et al. 2018;Moriya et al. 2019;Janka et al. 2022));.According to their results, strong FBA (10 −3 M ⊙ /s) can occur in the late phase (> 10s) for some progenitors.In this study, we increase the mass accretion rate to see its dependence on the neutrino luminosity.We note that it is necessary to check this dependence by decreasing the mass accretion rate for the sake of completeness, but these simulations are currently not available due to some technical problems associated to EOS tables.Addressing this issue is postponed to a future work.The matter temperature is set as T = 0.5 MeV.We note that this setup (cold FBA) leads to a conservative estimation of the neutrino signal by FBA.Y e is set to be 0.5.We run each model until the system reaches a quasi-steady state.
For computational reasons, the temperature inside 8 km is fixed in time.It is well inside the PNS; in fact the matter density is ρ > 1.5 × 10 14 g/cm 3 and its temperature is also low (∼ 2MeV), indicating that the boundary condition does not affect the neutrino emission.On the other hand, this prescription is necessary to prevent the PNS from over-cooling.In fact, if the temperature evolution is fully solved, its monotonic decrease with time may cause numerical crashes.As we shall show below, the PNS temperature does not affect the neutrino emission by FBA, indicating that this prescription does not compromise the present result.

Matter distributions
We first focus on the matter distributions of FBA after the system has settled to a quasi-steady state.The time it took to reach the steady state varies for different models, in the range of 50 ≲ t ≲ 120 ms.It took longer for lower accretion rate models.Figures 1, 2, and 3 show the density, temperature, and four-velocity of the fluid.In these figures, color distinguishes models with different mass accretion rates Ṁ and PNS masses for M PNS = 1.98M ⊙ , and three models As can be seen in these figures, an accretion shock wave is formed due to FBA.We note that similar phenomena are observed in recent multi-D CCSN simulations (see, e.g., Fig. 9 in Nagakura et al. (2021a)).According to these CCSN models, FBA tends to be cold or lower entropy (otherwise the thermal pressure hampers accretion), and the downflow onto the PNS becomes supersonic.At the surface of the PNS, the fluid needs to be subsonic, implying that a shock wave is inevitably formed.As displayed in these figures, the shock position is larger for the lower accretion rate and the lower PNS mass.This is attributed to the lower ram pressure in the preshock region.
Another notable feature displayed in Fig. 2 is that sharp peak profiles in the temperature distribution emerge in the vicinity of the PNS.To see the profile clearly, we magnify the corresponding region (10km ≤ r ≤ 13km), which is displayed in the same panel.The sharp increase of temperature is due to the hard wall of the PNS surface, in which the matter density changes by four orders of magnitude over a width ≲ 1 km.In this region, the kinetic energy of FBA can be efficiently converted into thermal energy, and therefore the temperature of FBA also increases rapidly with decreasing radius. 1 We also note that the thermal energy of matter is proportional to Ṁ v ff 2 , where v ff denotes the free-fall velocity where the kinetic energy is dissipated.Since the dissipation region is less sensitive to the PNS mass, v ff is proportional to M 0.5 PNS .This is the rationale behind the higher peak temperature for the higher mass accretion rate and the higher PNS mass (see Fig. 2).
Contrary to the trend which we have discussed so far, the temperature decreases rapidly with decreasing radius at ≲ 11.5km and ≲ 12km for M PNS = 1.98M ⊙ and 1.41M ⊙ , respectively (see the magnified figure of Fig. 2).This exhibits that neutrino cooling gives feedback on the matter distribution.On the other hand, the temperature profile is very complicated in the transition layer between the cold PNS envelope and the inner edge of FBA.As we shall show below, weak processes are responsible for the complex radial profile in the temperature distribution.It is also interesting to note that the matter profile in the region 10km ≲ r ≲ 11km does not 1 A smaller PNS radius means that a larger gravitational energy is converted into thermal energy.Since our focus is the late phase, the PNS radius is thought to have shrunk to a small radius due to cooling (in our setting, ∼ 11 km).This situation is actually advantageous for creating a high temperature peak and leads to a larger amount of neutrino emission.depend on the mass accretion rate.Although we postpone the detailed investigation to future work, this may be due to a self-regulation mechanism around the PNS surface.Since the matter pressure needs to be connected smoothly across the layer, the fluid element at the PNS surface undergoes shrinking.This implies that the gravitational energy is converted into thermal energy, which also accounts for the increase of neutrino luminosity, in particular for heavy-leptonic neutrinos (ν x ).

Neutrino distributions
Before going into details, let us first provide the information on species-dependent neutrino spheres.As a reference, we show them in the case of Ṁ = 10 −3 M ⊙ • s −1 and the PNS mass 1.98M ⊙ .The neutrino spheres for the energy of 23.4 MeV, roughly corresponding to the average energy of neutrinos, are 11.2 km 10.7 km, 3.66 km for ν e , νe , and ν x , respectively.The density at each neutrino sphere is 2.55 × 10 12 g • cm −3 , 1.88 × 10 14 g • cm −3 and 1.11 × 10 15 g • cm −3 , respectively.This exhibits that the neutrino sphere is located at higher matter density than in the early post-bounce phase (a few hundreds of milliseconds after core bounce); for instance, the neutrino sphere of ν e is located at ∼ 10 11 g/cm 3 in the early post-bounce phase.This difference can be understood as follows.The density gradient becomes so steep in the late phase, indicating that the scale height in this region becomes small.Since the neutrino optical depth is determined not only by the local reaction rate but also by the scale height, the optical depth tends to be smaller in the late phase for the region with the same matter density.It is also worthy to note that these neutrino spheres are located much deeper than the inner boundary adopted in the simulations of Fryer (2009).
To delve into the neutrino feedback on matter, we portray the radial profiles of the energy flux (F ν ) of each species of neutrinos in the case of the PNS mass 1.98M ⊙ and the accretion rate Ṁ = 10 −3 M ⊙ • s −1 (see Fig 4 .).In the figure, neutrino fluxes are multiplied by a factor r2 .We note that F ν r 2 is approximately constant in space, if there are no neutrino emission and absorp- tion 2 .This indicates that the information on neutrino cooling (or heating) is imprinted in the radial profile of F ν r 2 .As shown in the top panel of Fig 4, neutrino fluxes for ν e and νe increase with radius in the region 11km ≤ r ≤ 11.5km, indicating that these neutrinos are substantially produced there.The fact that ν x is approximately constant in space indicates that ν x is not produced in this region.It is also worthy to note that ν e absorption dominates over emission in the narrow region at ∼ 10.9km (see blue line).A similar profile is also observed for νe at smaller radius (see the red line).
It is also informative to see the temperature profile as a function of the matter density, which is displayed in the bottom panel of Fig 4 .For ν e (ν e ), strong neutrino production occurs at very high density 5 × 10 13 g/cm 3 ≲ ρ ≲ 2 × 10 14 g/cm 3 (10 14 g/cm 3 ≲ ρ ≲ 2 × 10 14 g/cm 3 ).When these neutrinos propagate out-wards in the lower density environment, neutrino absorption becomes dominant in the region 10 13 g/cm 3 ≲ ρ ≲ 5 × 10 13 g/cm 3 (5 × 10 13 g/cm 3 ≲ ρ ≲ 10 14 g/cm 3 ), but neutrino emission again dominates over absorption until ρ ∼ 10 9 g/cm 3 .These non-monotonic profiles of ν e and νe fluxes are clearly associated with the matter temperature profile, which shall be discussed later.Our result also suggests that there is a substantial amount of diffusion component for both ν e and νe in their energy fluxes, which are missing components in the simulations of Fryer (2009).We also find that ν x profile is much simpler than others; ν x is mainly produced at 5×10 13 g/cm 3 ≲ ρ ≲ 2×10 14 g/cm 3 and then they freely escape from the system.This suggests that ν x production mainly occurs in such a high density region.It is, hence, mandatory to cover the high density region in numerical simulations to quantify the neutrino signal from FBA.To see what weak process accounts for the neutrino emission and absorption, we display the radial profile of inverse mean free path of each weak process in Figs. 5  and 6, for the accretion rate of Ṁ = 10 −3 M ⊙ • s −1 and the PNS mass 1.98M ⊙ model.These figures show the inverse mean free path as a function of the ra-dius and the matter density, respectively.Except for the high density region (∼ 5 × 10 13 g/cm 3 for ν e and ∼ 2 × 10 14 g/cm 3 for νe ), electron-capture on free proton and positron-capture on free neutron dominate the ν e and νe emission, respectively.We also find that nucleon-nucleon bremsstrahlung becomes dominant in the high density regions, leading to a non-monotonic radial profile of neutrino opacity.This is responsible for the non-monotonic profile for both neutrino fluxes in ν e and ν.This leads to the complex radial profile of matter temperature.Indeed, the inverse mean free path peaks at 11 ≲ r ≲ 12 km for ν e and 10 ≲ r ≲ 11 km for νe , and these spatial positions are roughly the same as those at the temperature dips.This exhibits that temperature dips are caused by the neutrino cooling by ν e and νe .For ν x , on the other hand, the neutrino opacity is dominated by nucleon-nucleon bremsstrahlung at the emission region (5 × 10 13 g/cm 3 ≲ ρ ≲ 2 × 10 14 g/cm 3 ).Although the electron-positron pair production becomes dominant at ρ ≲ 10 12 g/cm 3 , the emissivity is very small.In fact, the radial profile of F ν r 2 for ν x is almost constant in space over the low density region (see Fig 4).

Neutrino luminosity and mean energy
Figure 7 summarizes the energy luminosity for all simulated models.As shown in the figure, larger accretion rates and larger PNS masses leads to higher luminosities.The energy luminosity for ν e and νe are of the order of o(10 51 ) erg • s −1 , and ν x luminosities are below 10 51 erg • s −1 .One thing we do notice here is that neutrino luminosities obtained in our simulations are systematically higher than those reported in Fryer (2009).This is again due to the fact that the simulations of Fryer (2009) did not cover the high density region, which results in underestimating neutrino luminosities.Our result suggests that it is mandatory to include the high density region in theoretical models to quantify the neutrino signal from FBA and to extract physical information from the neutrino signal in real observation (see Sec. 5 for more details).Luminosities for ν x hardly depend on the accretion rates.This is because ν x are mainly emitted from the inner PNS, as we saw in the radial profiles of the flux.The luminosities of ν e and νe are also not proportional to the accretion rates and shifted to higher values, due to the same reason.
In Fig. 8, we provide mean energies of the emitted neutrinos.The mean energy of ν e is ϵ ∼ 13 MeV for the highest accretion case ( Ṁ = 10 −2 M ⊙ •s −1 ), and it is still ϵ ≳ 10 MeV for other cases with lower mass accretion rate.The mean energy of νe is always higher than that of ν e , and reaches a maximum of ϵ ∼ 17 MeV.We also find that, similar to the luminosity, larger accretion rates and PNS masses lead to higher neutrino mean energies.This is due to the higher matter temperature in the neutrino emission region (see Fig. 2).
It is worthy to note that ν x has the lowest mean energy among three flavors, which is ϵ ∼ 10 MeV.This tendency is clearly different from the canonical hierarchy of neutrino mean energy in CCSNe.In general, the mean energy of ν x is the highest among all flavors of neutrinos in early post-bounce phase, and then the mean energy of all flavors becomes almost identical in the late phase.This exhibits that FBA leads to a qualitatively different neutrino emission from PNS cooling, and that the neutrino detection rate should depend on the neutrino oscillation model.The low luminosity and low mean energy of ν x are attributed to the temperature distribution in the ν x emission region.As shown in Fig. 4, most of the ν x are produced in the region of 5 × 10 13 g/cm 3 ≲ ρ ≲ 2 × 10 14 g/cm 3 .This region corresponds to the transition layer between the PNS surface and the inner edge of FBA.Although the matter temperature sharply increases with radius, it is still very low (≲ 4MeV).As a result, both the luminosity and the mean energy of ν x become much lower than those associated to ν e and νe .We note that the radius of the emission region for ν x is smaller than for other flavors.This causes a lower neutrino luminosity, although this effect is minor since the difference of emission region among all flavors of neutrinos is only ≲ 1km.
We remind the readers that our current focus is the late phase (t ≳ 10 s), where the typical mean energy of the diffusive neutrino component from PNS is ϵ ≲ 10 MeV (Suwa et al. 2019).This is much smaller than neutrinos from FBA.As discussed in Sec. 5, higher  luminosities and mean energies of neutrinos are more favorable for neutrino detection.Our results support the claim in Fryer (2009) that FBA can substantially increase the neutrino event rate, which is quantified in the next section.

Detectors and neutrino oscillation models
We evaluate the detectability of the FBA neutrinos by two representative terrestrial neutrino detectors, Super-Kamiokande (Hereafter Super-K) and Deep Underground Neutrino Experiment (DUNE).In this estimation, we employ the neutrino cross section data taken from the SNOwGLoBES (SNOwGLoBES : SuperNova Observatories with GLoBES).We ignore any smearing effects caused by the detector response and various noises just for simplicity.
Super-K is a water-Cherenkov detector using pure water (Fukuda et al. 2003) with gadolinium compound loaded recently (Abe et al. 2022).The main detection channel of Super-K is the inverse-beta interaction νe + p → e + + n. (4) We assume the fiducial volume of 32.5 kton for the estimation of the event rate.Its update version, Hyper-Kamiokande is also under construction (Hyper-Kamiokande Proto-Collaboration et al. 2018).Its fiducial volume will be 220 kton, and the detection rate can be easily scaled from the result of Super-K.We assume pure water for the evaluation of the event rates.It should be mentioned that the gadolinium-loading in SK plays an important role to decouple the FBA neutrino signal from the background (Li et al. 2022;Simpson et al. 2019).Unlike the strong neutrino burst in the early post-bounce phase, the luminosity is lower and the timescale is longer for FBA indicating that the reduction of the background is very important to identify the signal.
DUNE is a future-planned neutrino detector.It will use liquid argon as the neutrino detector medium.The main detection channel of DUNE is the neutrino-argon charged-current interaction ν e + 40 Ar → e − + 40 K * . (5) We assume a full volume of 40 kton for the estimation of the event rate.
For the estimation of the neutrino flux arriving on the earth, we take into account the neutrino oscillation effect in the same way as Dighe & Smirnov (2000); Nagakura et al. (2021a).Neutrino flavors are assumed to convert adiabatically by the Mikheyev-Smirnov-Wolfenstein (MSW) effect.Although it is a simple oscillation model, this provides an essential feature of how detectability of FBA neutrinos depends on flavor conversions.
Following Nagakura et al. (2021a), the neutrino fluxes arriving on earth F e , Fe , F x , Fx (corresponds to ν e , νe , ν x , νx , respectively) are calculated from the values of the fluxes without neutrino oscillation (F 0 e , F 0 e , F 0 x , F 0 x ) as: where p, p are survival probabilities.In the normal-mass hierarchy case, they are defined as On the other hand, in the inverted-mass hierarchy case, they are defined as The values of the neutrino mixing angles θ 12 , θ 13 are assumed to be sin 2 θ 12 = 2.97×10 −1 and sin 2 θ 13 = 2.15× 10 −2 , adopted from Capozzi et al. (2017).We assume F 0 x = F 0 x in this study.

Neutrino Event Rates
Figure 9 shows the neutrino event rates per unit time, in which we integrate over energy, while the energydependent feature is discussed later.The distance is assumed to be 10 kpc.As shown in Fig. 9, the event rate clearly depends on the mass hierarchy, where the difference is more than double.In the case with the normal(inverted)-mass hierarchy, p (p) becomes small, indicating that neutrinos (anti-neutrinos) undergo large flavor conversions.As shown in Sec.4.3, both the energy luminosities and the average energies of ν e and νe are higher than those of ν x at the source, indicating that the large flavor conversion results in reducing the ν e and νe number flux.As a result, the number of event rate at Super-K and DUNE becomes lower in the case of inverted-mass hierarchy and normal one, respectively.Hence, simultaneous observation of FBA neutrinos with Super-K and DUNE will provide a strong constraint on neutrino mass hierarchy.
The estimated event rate is found to be o(10) s −1 for the accretion rate of Ṁ ∼ 10 −3 M ⊙ •s −1 .This result also suggests that if we detect a large number of neutrinos in the very late phase, the detection will be an evidence for the occurrence of FBA neutrinos.It is also worthy to note that similar accretion rates were found in previous studies Chan et al. (2018); Moriya et al. (2019); Janka et al. (2022) in the late phase.
We also find that the dependence of the event rate on the mass accretion rate hinges on the neutrino oscillation model.In the case of normal (inverted) mass hierarchy, the detection rate at Super-K (DUNE) becomes remarkably higher for higher mass accretion rates, whereas it is less sensitive to the accretion rate in the case of inverted (normal) one.This trend can also be understood through the species-dependent feature of neutrino emission at the CCSN source.As shown in Figs.7 and 8, both the luminosity and the average energy of ν x weakly depend on the mass accretion rate, and therefore the large flavor conversion makes the detection count at each detector less sensitive to the mass accretion rate.Nevertheless, the number of event count at each detector is remarkably higher than that emitted from the PNS.As a reference, we show the case for the neutrino signal without FBA but only with isothermal PNS of T = 3MeV in Fig. 9 (see below for the details).This figure illustrates that the detection rate of neutrinos from FBA is remarkably higher than in the case with thermal neutrinos from the PNS.It is also worthy to note that flavor dependent features would be resolved by using other reaction channels or joint analysis with other detectors (see, e.g.Beacom et al. 2002;Dasgupta & Beacom 2011;Nagakura 2021), that would provide a key information to distinguish the neutrinos powered by FBA from those radiated only from the inner region of the PNS.
In Figs. 10 and 11, we show the energy-dependent event rate, in which the energy bin is set to be 1 MeV.As references, we also make a plot for the case with the purely thermal emission (Fermi-Dirac distribution with zero chemical potential) with a PNS temperature of 2 and 3MeV.The emission radius is assumed to be 11 km, and we also take into account the gravitational redshift in this estimation.For all simulated models, the event rates are orders of magnitude larger than the background event rates (see the latest experimental data of Super-K in Harada et al. (2023b).These figures illustrate that a large neutrino emission can be expected in the case of higher mass accretion rate.Another notable feature found in these figures is the high energy tail in each spectrum.Even in the case with the low mass accretion rate ( Ṁ = 10 −3 M ⊙ •s −1 ), neutrinos with ≳ 30MeV may be observed.It should be noted that these high energy neutrinos cannot be detected in the late phase (t ≫ 10s) by thermal emission of PNS, unless the source is extremely close (Nakazato et al. 2022).If we detect them in real observations in the late phase, these neutrinos would be generated by FBA.We note that Figs. 10 and 11 show the energy event rate per second, indicating that the actual event count may be a factor of > 10 larger than this value (since we are currently considering in the phase of > 10s after core bounce).

PNS Temperature Dependence
It interesting to see how the neutrino signal from FBA depends on the PNS temperature.The increase of the PNS temperature would lead to higher neutrino emission inside the PNS, which potentially alters the neutrino signal.In this test, we employ the same numerical setup as that used in our model with the PNS mass of 1.98M ⊙ and the accretion rate of Ṁ = 10 −3 M ⊙ • s −1 except for the PNS temperature.We consider two cases: 3 MeV and 4 MeV.We note that T = 4 MeV is too hot for PNS in the late phase which we consider in this paper (> 10s after core bounce), but the result is still informative.
In Fig. 12, we show the energy spectrum of the neutrino event rate at Super-K and DUNE in the case of normal-and inverted mass hierarchy, respectively.We note that each oscillation model corresponds to the case having the lower number of event rate than the other mass hierarchy.As shown in these figures, even in these pessimistic cases, the PNS temperature does not affect the neutrino event rate.This result supports the claim that neutrinos from FBA overwhelm the thermal neutrinos from the PNS, unless they are extremely hot (T ≫ 4 MeV).

SUMMARY AND CONCLUSION
In this paper we investigated neutrino emission from fallback mass accretion (FBA) onto PNS in the late phase of CCSN (> 10s) by using general relativistic neutrino radiation-hydrodynamic simulations with full Boltzmann neutrino transport.In our numerical simulations, we covered the very high density region (> 10 14 g/cm 3 ) where we set a quasi-steady PNS structure as initial conditions.We changed the mass accretion rate in a parametric manner, and ran each simulation until the system settled to a quasi-steady state.
We found that a higher accretion rate and a higher PNS mass leads to a higher temperature in the transition layer from the PNS surface to the inner edge of FBA, where most of the neutrinos are radiated.As a result, both luminosities and mean energies of neutrinos tend to be higher with increasing mass accretion rates.On the other hand, the sensitivity of neutrino emission to the mass accretion rate hinge on neutrino species.Although ν e and νe emission strongly vary with the mass accretion rate, ν x is less sensitive.This is due to the fact that ν x is produced in the highest density region (ρ ≳ 5 × 10 13 g/cm 3 ), indicating that the impact of FBA on the temperature distribution tends to be weak.Nevertheless, both the luminosity and the mean energy of ν x are remarkably higher than those estimated by standard PNS cooling models.The present study supports the claim by Fryer (2009) that FBA can substantially change the neutrino emission in the late phase of CCSN.On the other hand, we also find that most of the neutrinos by FBA are produced in the high density region which the simulations of Fryer (2009) did not cover.As a result, the neutrino luminosities in his estimation are underestimated by a factor of ≳ 5, and this systematic error has a non-negligible effect to extract physical information from neutrino signal in real observations.We also find that the dominant weak processes for neutrino emission depends on species: electron-capture by free proton, positron-capture by free neutron, and nucleon-nucleon bremsstrahlung for ν e , νe , and ν x , respectively.Although the electron-positron pair can be a dominant emission process for ν x in the low density region, the emissivity is too low to change the neutrino flux.
Based on the numerical results, we estimate the expected event rate for Super-K and DUNE with the adiabatic MSW oscillation model.One thing we need to stress is that neutrino emission from FBA has a rich flavor-dependent structure, indicating that the neutrino observation should depend on the neutrino oscillation model.Indeed, the difference of event rate between normal-and inverted mass hierarchy at each detector becomes more than double.In short, the detection rate tends to be smaller if the flavor conversion is strong.This is attributed to the fact that ν x luminosity and mean energy are systematically lower than those of other species.Nevertheless, the event rate is the order of o(100) s −1 for the optical case with the highest accretion rate in both detectors, and still o(10) s −1 for the least optimal setting, which is much larger than the canonical PNS cooling model.We also provide energy-dependent features in the neutrino signal.We find that the peak energy of neutrino detection is remarkably higher than the thermal emission of PNS with ≤ 3MeV.Our result suggests that high energy neutrinos (≳ 30MeV) may be observed in the late phase, which will be evidence that neutrinos are emitted by FBA.
As a final remark, we point out a couple of limitations in our study.First, we assumed spherical symmetry.In the multi-D case, the accretion shock wave may be unstable to non-radial perturbations (Blondin et al. 2003;Yamasaki & Yamada 2005, 2006, 2007;Foglizzo et al. 2007), and FBA is usually accompanied by turbulence (Vartanyan et al. 2022), which potentially leads to temporal variations in the neutrino signal.On the other hand, it would be hard to resolve the temporal variation by the current-and even future-planned neutrino detectors, unless the CCSN source is very close (see, e.g., Nagakura et al. 2021b).This is because the neutrino luminosity is very low in the late phase, and the temporal variation would be smeared out by noise.It should be mentioned, however, that the thermodynamical properties in the post-shock flow may be influenced by the shock instability, which may change the neutrino signal.We postpone this detailed study to future work.Second, the number of models simulated in this study is limited due to the computational cost.It should be stressed that high spatial resolution is required to resolve both matter and neutrino distributions around the surface of the PNS, implying that the timestep is severely limited by the Courant condition.To prepare for future observations, however, we need a systematic study by covering wider ranges of PNS masses and mass accretion rates than those studied in the present paper.The EOS dependence is also worthy of investigation.Numerical simulations are not suitable to carry out such a systematic study, and therefore we are planning to take a semi-analytic approach to address this issue.If we cover full parameter space, we may be able to infer the EOS parameters or accretion rates from the future neutrino detection.An analysis pipeline based on a Bayesian approach has been already developed for thermal neutrino detection (Harada et al. 2023a).The results with similar approach will be reported elsewhere.

Figure 1 .
Figure 1.Radial profiles of the density for the 1.98M⊙ model (top) and the 1.41M⊙ (bottom).The different colors indicate different accretion rates.

Figure 2 .
Figure 2. Same as figure 1, but for the temperature.

Figure 3 .
Figure 3. Same as figure 1, but for the four-velocity.

Figure 4 .
Figure 4. Energy flux times the square of radius for the PNS mass 1.98M⊙ and the accretion rate Ṁ = 10 −3 M⊙•s −1 .The horizontal axis are the radius and density, for top and bottom panels, respectively.

Figure 5 .
Figure 5. Radial profile of the inverse mean free path for νe (top), νe (middle) and νx (bottom) with the energy of 23.4 MeV, for the model with the accretion rate Ṁ = 10 −3 M⊙ • s −1 and the PNS mass 1.98M⊙.The abbreviation of the neutrino reactions are as follows: the electroncapture on nucleon (ecp), the positron capture (aecp), the electron-positron process (pap) and the nucleon-nucleon bremsstrahlung (nbr).

Figure 6 .
Figure 6.Same as figure 5, but the horizontal axis is the density.

Figure 8 .
Figure 8. Same as 7, but for the mean energy of emitted neutrinos.

Figure 9 .
Figure9.The event rate of neutrinos for the Super-K (top) and DUNE (bottom), assuming the distance of 10 kpc.The horizontal lines denote the event rates assuming thermal emission.

Figure 10 .Figure 11 .
Figure10.The event rate per 1 MeV energy bin for MPNS = 1.98M⊙ model, assuming the distance of 10 kpc.The top and bottom panel corresponds to normal and inverted mass hierarchy, respectively.The left panels are for Super-K and right panels are for DUNE.The gray plots corresponds to the numerical results, and the purple and blue plots corresponds to thermal values.

Figure 12 .
Figure 12.The event rate per 1 MeV energy bin for different PNS temperatures; 2 MeV (blue), 3 MeV (purple) and 4 MeV (red), assuming the distance of 10 kpc.Top panel is for Super-K and bottom panel for DUNE.Mass hierarchy is inverted (top) and normal (bottom).PNS mass is MPNS = 1.98M⊙ and the accretion rate is Ṁ = 10 −3 M⊙ • s −1 .