Neutrinos from the Galactic Center Hosting a Hypernova Remnant

Similar to star-forming galaxies or starburst galaxies, star-forming regions in our Galaxy can host cosmic-ray (CR) accelerators and rich gas as targets of hadronuclear interaction. By our estimations, the IceCube neutrino observatory might detect muon neutrinos from a CR accelerator associated with a molecular cloud complex in our Galaxy. The associated high-energy gamma-ray emission might be observed by the Cherenkov Telescope Array (CTA), High-Altitude Water Cherenkov Gamma-Ray Observatory (HAWC), and Large High Altitude Air Shower Observatory (LHAASO). Furthermore, taking the Galactic Center (GC) region as an example, we assume that a hypernova exploded in the past in the GC. We simulate the acceleration of CRs in the hypernova remnant (HNR) as well as their confinement and escape. The high-energy protons escape from the HNR, diffuse around the GC, interact with molecular clouds, and then produce gamma-rays and neutrinos. In the optimal cases, the GC would be a promising 100 TeV gamma-ray source for LHAASO’s one-month observation. We propose that neutrino-induced searching for starting track-like and high-energy starting events (HESEs) observed by IceCube, from the GC region with a radius of 1.°8, would help us discover the particle accelerator in the GC or constrain our models. Under the constraint from high-energy gamma-ray observations by the H.E.S.S. telescope, we estimate the exposure time needed to make a significant discovery for the optimal cases. The analysis combining observations of IceCube and ANTARES, starting track-like events and HESEs, future observations by neutrino detectors IceCube-Gen2 and KM3net, and gamma-ray telescopes CTA, HAWC, and LHAASO would help to constrain our models.


Introduction
The origin of high-energy neutrinos observed by the IceCube neutrino observatory is still unresolved. Generally, high-energy astrophysical neutrinos can be produced by the accelerated cosmic rays interacting with matter via hadronuclear interactions and with radiation via photohadronic interactions.
Similarly, in our galaxy, star-forming regions, for instance, the Cygnus-X region (Yoast-Hull et al. 2017), and compact young star clusters, for instance, Westerlund 1 (Bykov et al. 2015;Aharonian et al. 2019), can also host massive star explosions as CR accelerators and provide the dense ISM or molecular clouds (MCs) that interact with the accelerated CRs, leading to the production of gamma-rays and neutrinos. Because they are much closer to Earth than extragalactic sources, their contributions to neutrinos might be resolved more easily than those distant SFGs/SBGs, depending on the total energy of the accelerated CRs and the mass of MCs. Furthermore, past HNe and other unusual SNe in Milky Way have been considered in connection with ultrahigh-energy nucleus acceleration (Calvez et al. 2010). Young TeV supernova remnants (SNRs) associated with MCs in the Galaxy, such as W51C, W44, IC 443, and W49B in the northern sky, and Sgr A East, Kes 75, 3C391, RX J1713.7-3946, CTB 37A, and 1FGL J1717.9-3729 in the southern sky are selected as possible neutrino source candidates in Aartsen et al. (2013Aartsen et al. ( , 2014. It is reasonable to look for contribution of past galactic HNe to the observed high-energy neutrino flux. In this paper, we first estimate a general contribution to neutrinos and gamma-rays from a past CR accelerator associated with an MC complex in the star-forming regions in our galaxy in Section 2. In Section 3, we apply the similar estimation to the GC region. Then, we assume that an HN in the GC region exploded in the past and injected protons to the central molecular zone (CMZ) and the ISM around the GC region. Under constraints from current gamma-ray observations on the GC region by the High-Energy Stereoscopic System (H. E.S.S.) Collaboration, we derive optimal spectra and templates of muon neutrinos from the GC. We predict the expected counts of neutrino-induced events and estimate the exposure time needed to make a significant discovery for the optimal cases.

Neutrinos from a Past CR Accelerator Associated with an MC Complex
We assume that, in a star-forming region of our galaxy, at a time T ago, a petaelectronvolt (PeV) CR accelerator started to accelerate protons to energies as high as  p,max . According to the standard theory of diffusive shock acceleration (DSA), the energy distribution of accelerated particles follows a power-law spectrum s 2 1 , with γ as the compression ratio (Fermi 1949;Drury 1983;Malkov & Drury 2001). Adopting γ=4 for strong shocks, the powerlaw index is s=−2; therefore, we assume a flat spectrum for the accelerated particles here. When protons are accelerated to higher than a certain energy ò p,esc , their diffusive length is larger than the scale of the acceleration site, then the chance that they escape the acceleration site significantly increases, while the lower-energy protons are still confined in the acceleration site as their diffusive length is shorter. Therefore, the spectrum of escaped protons would have a lower-energy cutoff at an energy of  p,esc (Ellison & Bykov 2011). Assuming the PeV accelerator is a young HN remnant, the energy range of the escaped protons depends on the age of the remnant, the energy of the SN explosion, the ejecta mass, the density of the medium, and so on (Ellison & Bykov 2011). According to our simulation in Section 3, it is reasonable to assume that   100 TeV p,esc and   1 PeV p,max for a young HN. For the sake of simplicity, we assume that the escaped protons follow a flat spectrum, the same as accelerated protons, while a more realistic spectrum does not follow a power-law form (Ellison & Bykov 2011). We denote the total energy of protons which escaped from the accelerator and are injected into the MCs by E inj , and the ratio between E inj and the total energy of accelerated protons E acc as as the diffusion coefficient for protons with energy of ò p =100 TeV. The diffusion coefficient D o and the diffusion index δ are uncertain parameters; their values are usually set to be ( )  -´-D 1 8 10 cm s o 29 2 1 and δ;0.3-0.6 for different empirical diffusion models of CR propagation (Strong et al. 2007). In Section 3, we will discuss cases with different values of δ. In this section, we will fix δ=0.6 and calculate the diffusion and interactions of petaelectronvolt protons. For the diffusion time of T=10 4 T 4 yr, the diffusion radius of protons with energy ò p are approximated to be is the distance of the accelerator from Earth. In the star-forming region, the accelerator is plausibly associated with the MC complex, and the accelerated protons will interact with matters in the MC complex via inelastic pp collision, leading to the production of neutrinos and gammarays. If the total mass of MCs contained within the diffusion region with a radius of R diff,PeV is M, then the average density of the matter within the diffuse region is approximated to be  (Gaisser 1990;Loeb & Waxman 2006), we derive the timescale of the proton energy loss dominated by the inelastic pp collision to be as the inelasticity (Gaisser 1990;Loeb & Waxman 2006). Here we assume that the energy loss timescale t loss is much longer than the diffusion time T of those injected protons; thus, the energy loss of protons during the propagation can be ignored.
A third of the proton energy is converted into neutral pions via inelastic pp collision (Gaisser 1990;Loeb & Waxman 2006), thus the neutral pions decay into gamma-ray photons, i.e., p gg  0 . According to Kelner et al. (2006), the energy distribution of secondary gamma-rays peaks at about  g   0.1 p . The luminosity of gamma-rays at the energy ofg  -10 100 TeV is approximated to bē ( ) The gamma-ray flux at the energy of g  10 100 TeV is about where the factor ln(10) denotes that the energy of gamma-rays are distributed over one energy decade from 10 to 100 TeV following a flat spectrum distribution. From Equation (5), one can see that the observed gamma-ray flux depends on the age of the CR accelerator T, the total energy of injected protons E inj , the total mass of MCs in the diffuse region M, the diffusion coefficient D o , and the distance of the source d s .
The Cherenkov Telescope Array (CTA) South is most sensitive at energies around 3-10 TeV, with the 50 hr sensitivity of´-  (5), as long as the source is located within their field of view. Because LHAASO is background free above 100 TeV, simply, the exposure time needed to detect the gamma-ray emission from the source can be estimated to be 1 for 10-100 TeV photons (Amenomori et al. 2019), which is comparable with our prediction in Equation (5) for benchmark parameters, indicating that the Tibet AS array with the underground water-Cherenkov-type MD array also has the capability to detect high gamma-ray emission from a CR accelerator associated with an MC complex in the Galaxy.
About two-thirds of the proton energy is carried by charged pions. The charged pion will then decay to produce four leptons, via processes (¯) p n n m   , which share the pion energy equally, in approximation. Then, the fraction of the total proton energy deposited in neutrinos is ´= 2 3 3 4 1 2 , with the energy of each neutrino  n   0.05 p (Loeb & Waxman 2006). Then, the luminosity of neutrinos at the energy of ò ν ;5-50 TeV is The flux of neutrinos at the energy of ò ν ;5-50 TeV is approximated to be Considering equipartition among the three neutrino flavors after their oscillations during propagation, the luminosity of muon neutrinos is n L 3. The counts of muon neutrinos at  n m  50 TeV detected by IceCube is 2 is the effective area of the detector for muon neutrinos with energyñ m  50 TeV, and Δt= 10Δt 1 yr is the exposure time. We assume that the source distance is 10 kpc, and the injected proton energy is = E inj -5 10 erg s 50 1 . As the effective area of IceCube is a function of the source declination (Yacobi et al. 2014), about one muon neutrino might be detected from a source in the northern hemisphere during 10 yr operation of IceCube, adopting the benchmark parameter set, while more operation time of IceCube is needed to discover sources in the southern hemisphere as through-going muon neutrinos from the southern hemisphere are contaminated by muon backgrounds.
HNe, with a typical kinetic energy of E k ∼(0.5-5)× 10 52 erg and a typical ejecta velocity of~-V 10 cm s  SN 1997ef, SN 1997dq, SN 1998bw, and SN 2002ap (Iwamoto et al. 1998Mazzali et al. 2002Mazzali et al. , 2004, are believed to be able to accelerate protons to PeV. The energy of injected protons of -100 TeV 1 PeV is 10 erg , where = -f f 0.1 conv conv, 1 is the fraction of kinetic energy that is converted into accelerated protons, and = f 0.17 inj is the ratio of the total energy of injected protons to that of accelerated protons, as defined in the beginning of this section. In our Galaxy, the estimated supernova rate is about one to three SNe per century (Adams et al. 2013). The ratio of the HN rate to SN Ibc rate is 7% in the local universe (Guetta & Della Valle 2007), and SNe Ibc contribute to 11% of the total SN rate (Cappellaro et al. 1999); therefore, the HN rate is about 1% of the SN rate. As a result, in our Galaxy, the total HN rate can be estimated to be about --10 yr 4 1 . So, in 10 4 yr, there might be one HNe that explodes in the Galaxy. If it is associated with a high-mass MC complex, IceCube, CTA, HAWC, and LHAASO have a chance of observing the association between neutrinos and gamma-rays.
Therefore, we propose to analyze the gamma-ray data of LHAASO, HAWC, and CTA in the direction of IceCubeobserved muon neutrinos associated with MC complexes or star-forming regions around the galactic plane. The future detection of >10-100 TeV photons associated with >5-50 TeV muon neutrinos will indicate the existence of a petaelectronvolt accelerators (PeVatrons), for example, an HN, embedded in the MC complex.

The Implication for the Galactic Center
The GC region, with distance to Earth d GC ∼8.15 kpc (Reid et al. 2019), contains a significant group of massive supergiants and Wolf-Rayet stars (Kauffmann 2017;Lu 2018), indicating a high rate of SNe/HNe/GRBs. The existence of the Fermi Bubble might be evidence of past star formation activity or central supermassive BH activity (Su et al. 2010;Fujita et al. 2015). The study by Yusef-Zadeh et al. (2009) suggests that the star formation rate in the GC region peaked around 10 5 yr ago.
These facts indicate that there might be possible cosmic-ray accelerators in the GC region, for example, an HN remnant, which will be discussed in Section 3.2; the blast wave formed in the tidal disruption event (TDE) caused by the supermassive BH Sagittarius A * (Sgr A * ; Liu et al. 2016); and a radiatively inefficient accretion flow (RIAF) of the low-luminosity AGN (LLAGN) Sgr A * (Fujita et al. 2015). The accelerated CRs are ejected into the surrounding environment and then interact with the dense matter via inelastic pp collision when propagating around the GC region, leading to the production of high-energy gamma-rays and neutrinos simultaneously. This assumption makes the GC region one of the important galactic source candidates of high-energy gamma-rays and neutrinos.
The H.E.S.S. Collaboration observed gamma-ray emission with energy up to ∼10 TeV, associated with the MCs around the GC, implying that a possible PeVatron exists in the GC region (HESS Collaboration et al. 2016). Adopting the constraint from the very high-energy gamma-ray observations by H.E.S.S., we will calculate the optimal contribution to neutrinos from a PeVatron in the GC.

The Estimation
As protons propagate around the GC, they interact with the interstellar gas via inelastic pp collision along their path, producing charged and neutral pions. We assume that protons, with energy of ò p and amount of N inj (ò p ), are injected to the GC region instantly at a time T ago. For this impulsive injection, the density of protons with energy of ò p at radius r is approximated to be (Aharonian 2004) where R diff is the diffusion radius of protons with energy ò p , as defined in Equation (1). For r<R diff (ò p ), the proton density is approximated to be a constant, i.e., ( ) The luminosity of gamma-rays and neutrinos is approximated to be p inj p is the total energy of injected protons with energy ò p , n MZ is the particle density of the molecular zone averaged over its volume V MZ within the CR diffuse region, and m p is the proton mass. Therefore, the luminosity of gamma-rays and neutrinos is proportional to the mass of MCs = M m V n MZ p MZ MZ within the CR diffuse region along our line of sight.
Based on observations, Nakanishi & Sofue (2006 derived the three-dimensional distribution map of the ISM, i.e., H I and H 2 maps, throughout the Galaxy. As the resolution of the ISM map is not high enough to resolve the structure of the CMZ in the GC region, we adopt the 3D map from Nakanishi & Sofue (2006 only for the ISM outside the CMZ. According to HESS Collaboration et al. (2016), the total mass in the inner 150 pc of the CMZ is approximated to be ( ) , corresponding to a uniform density of n 150 cm CMZ 3 . We adopt a disk-like CMZ model with radius of R CMZ =250 pc, height of h CMZ =70 pc, and uniform density of n CMZ =150 cm −3 , leading to the total mass of the CMZ of about M CMZ =9×10 7 M e . Hereafter, we mark this circular area with radius of R OZ =1°.8 as Area A 1°. 8 .
In a larger scale outside the CMZ, there is also a disk-like zone of dense ISM with the height of h OZ ∼0.25 kpc (Nakanishi & Sofue 2006, and the total mass of ISM in the outer zone (OZ) with the radius of R OZ =1 kpc (R A =7°.0 in angular scale) and the height of h OZ =0.25 kpc is about M OZ =5×10 8 M e , corresponding to an average density of = n 28 cm OZ 3 . Hereafter, we mark this larger circular area containing both the CMZ and OZ as Area A 7°. 0 .
If protons diffuse to 1 kpc, the proton density within the radius of 1 kpc is approximated to be a constant. Aside from the contribution from the CMZ, the ISM in the OZ along the line of sight also contribute to the gamma-ray and neutrino emission. According to Equations (11) and (12), the ratio of the contributions from two different regions is the same as the ratio of the MC mass in the two regions along the line of sight, i.e., M OZ /M CMZ ∼2.6 and M OZ /M CMZ ∼0.9, for Area A 1°. 8 and Area A 7°. 0 , respectively. Therefore, the contribution from the OZ cannot be ignored if protons diffuse to 1 kpc.
Assuming that protons diffuse to the radius of 1 kpc, within the diffusion region, along our line of sight, Area A 1°. 8 contains MCs with a mass of . Therefore, besides the gamma-ray and neutrino emission from the central region with radius of 1°. 8, future observatories, such as CTA or IceCube Gen II, might also be able to observe the diffusive gamma-ray and neutrino emission with a comparable luminosity outside the central region, with the image extended to more than a few degrees.
The H.E.S.S. Collaboration reported the gamma-ray spectrum from the annulus centered at Sgr A * with inner radius of 0°. 1 and outer radius of 0°. Then we can constrain the diffusion radius as a function of the MC mass M annu and the energy of the injected 100 TeV protons E inj (100 TeV), by setting ( Adopting the total MC mass in the direction of the annulus along the line of sight M annu =3×10 7 M e , we derive the constraint as Furthermore, according to Equation (1), the product of the diffusion coefficient D o and the propagation time T is constrained to be The above estimation does not depend on the detailed model of the CR accelerator. In Section 3.2, we will discuss the situation with an HNR as the CR accelerator.

Simulations
Yusef-Zadeh et al. (2009) (Fukugita & Kawasaki 2003), one can derive the SN rate in the GC aś --1.7 10 yr 3 1 . As mentioned in Section 2, the ratio of the HN rate to the normal SN Ibc rate is 7% in the local universe (Guetta & Della Valle 2007), and the normal SN Ibc rate is about 11% of the total SN rate (Cappellaro et al. 1999); therefore, the HN rate in the GC is about´--1.3 10 yr 5 1 around 10 5 yr ago. If an HN exploded a certain time T ago, the ejecta from the HN explosion swept up the ambient medium, then formed an HNR. The shocks in the HNR are able to accelerate CRs. Meanwhile, CRs at the high-energy end are difficult to confine, and thus are able to escape from the acceleration site, leading to a continuous injection of CRs with a duration of -10 10 yr 4 5 . By performing hydrodynamic simulations coupled with nonlinear DSA (Lee et al. 2012) for the evolving HNR, we derive the spectra of protons that are accelerated in the remnant, and those that have escaped from the remnant and been injected into the surrounding environment. The simulation code is composed of a spherically symmetric Lagrangian hydrodynamics and a semianalytic treatment of nonlinear DSA with magnetic field amplification by streaming instability of CRs in the shock precursor. The code offers end products such as the space distribution and time evolution of CR and electromagnetic spectra.
The acceleration of protons is described by the diffusionadvection equation, which is written in the shock rest frame as follows: A A Here in our simulations, Bohm diffusion is considered, i.e., the diffusion coefficient ( 3 is a function of its position x and the momentum p, with v(p) as the speed of the proton, ( ) dB x as the strength of the turbulent magnetic field at position x, and e is the charge of a proton. u(x) is the speed of the upstream flow and ( ) v x A is the speed of the Alfvén waves in the shock rest frame. Q(x, p) is the injection rate from the shock-heated plasma. f (x, p) is the phase-space distribution of the DSA-accelerated protons, which is rewritten as with β cut =1 and α cut =0.5 to account for the cutoff at the maximum momentum p max . The maximum momenta of the accelerated protons p max can be calculated by equating either the acceleration timescale with the remnant age (age limited) or the upstream diffusion length with the length of a free escape boundary L feb ≡f feb R HNR (t) (escape-limited), where f feb =0.2 is fixed in this simulation and R HNR (t) is the radius of the HNR. A uniform ambient medium around the HNR is adopted for simplicity. Because it is uncertain whether the ambient medium surrounding the HN is similar to the ISM, or is very dense due to the early mass loss of the star before the explosion, we study two cases with different densities of the ambient medium for the HN explosion, i.e., = n 1 cm 3 and = n 100 cm 3 . For the later case, the ejecta interacts with a denser environment, and its kinetic energy will dissipate into the accelerated protons faster, which is marked as the "fast dissipation case"; the former case is marked as the "slow dissipation case." The total kinetic energy of the ejecta we adopt in the simulations here is E HN =3×10 52 erg and the total ejecta mass is M ej =14M e . The other parameters are the same as Model II in Lee et al. (2012; see their Table 1), which is calibrated to the multiwavelength observational data of the TeV-bright young Galactic SNR RX J1713.7-3946, except that a stellar wind is replaced by a uniform ambient medium for simplicity and that the ejecta is swapped with an HN progenitor model as mentioned above. In particular, an injection parameter, defined as the ratio between the momentum of the proton that is injected into the shock p inj to that of the thermal proton p th , i.e., c º = 3.75 p p inj inj th for DSA is used, which is typical for young SNRs with strong shocks (e.g., Lee et al. 2013;Slane et al. 2014). For more details of the treatment of hydrodynamics, DSA, magnetic field amplification, and other physics included in our model, please refer to Lee et al. (2012) and references therein.
We approximate the continuous injection of protons from an evolving HNR into a series of impulsive injections. We mark the instant of the explosion as T 0 =0 and the instant of the i th impulsive injection as t i . The energy distribution of the injected protons ( )   N p p p for each impulsive injection is shown in Figure 1; the color code represents the value of t i . The peak energy of those injected protons increases first and then decreases due to the expansion of the HNR. Our simulation shows that the HNR can inject protons withPeV energy for a duration of a few thousand years for the slow dissipation case and a few hundred years for the fast dissipation case, and the total injected energy of protons with energy above 1 PeV is about ( ) E PeV 10 erg inj 51 . For the fast dissipation case, protons are injected with a higher peak energy and a shorter injection duration. The different features of the proton injection for the two cases will lead to different features of gamma-ray and neutrino productions.
The diffusion time of protons from the i th impulsive injection is , then the corresponding diffusion radius for protons with energy ò p is . As a result, at the current time instant T, the number density of protons at the radius of r from the ith impulsive injection is approximated to be (Aharonian 2004) p is the count of protons with energy ò p from the ith impulsive injection.
One can calculate the contribution to photons or muon neutrinos at the coordinate of (x, y, z) from the ith impulsive injection via where the coordinate of the GC is (0, 0, 0), ) n x y z , , H is the matter density at the coordinate (x, y, z), and s denotes gamma-ray photons γ or muon neutrinos ν μ . The function ( ) Kelner et al. (2006), and F ν ( f x , ò ν /f x ) is defined by the sum of Equations (62) and (66) in Kelner et al. (2006). The total gamma-ray or neutrino flux observed on Earth is the sum of the contribution from each impulsive injection, i.e., where V is the volume of the emitting region along the line of sight with dV=dxdydz. The HN age, T, the diffusion coefficient of protons at energy of 100 TeV, D o , and the diffusion index, δ, are free parameters.

The Resulting Spectra
Adopting the IceCube effective area for muon neutrinos from the direction of the GC, as adopted in Yacobi et al. (2014), we can derive the expected counts of the astrophysical muon neutrinos. We calculate the counts of muon neutrinos N νμ from Area A 1°. 8 observed by IceCube over 10 yr of operation, for different HNR age T and diffusion coefficient D o , with a fixed diffusion index δ=0.6. The contours of muon neutrino counts as a function of T and D o are shown in Figure 2. The black solid lines are constraints from the H.E.S.S. observation of gamma-rays from the central annulus at 10 TeV, i.e., 10 TeV 10 TeV 1 10 GeV cm s obs 9 2 1 . The parameter space above the black solid line is allowed. We  can tell that the constraint on the product of T and D o is roughly a constant, which is consistent with our estimation in Equation (15). Because the proton injection duration of our simulation is as long as 6×10 4 yr in the slow dissipation case, if the HN age <T 6 10 yr 4 , the total injected proton energy is not a constant. Thus, the product of T and D o is not a constant for small T, as seen in the left panel of Figure 2.
Under the constraint of the H.E.S.S. observation on gammarays from the central annulus, we plot the gamma-ray spectra (the solid black lines) from the same region for four optimal cases in Figures 3 and 4. The value of the diffusion index δ is uncertain. For a fixed diffusion coefficient for 100 TeV protons, a smaller δ indicates that PeV protons will diffuse slower than in the case with a larger δ. As a result, the density of PeV protons in the case with a smaller δ is higher than that in the case with a larger δ, under the same constraint on 10 TeV gamma-ray photons, leading to higher gamma-ray and neutrino fluxes at the high-energy end. In Figure 4, the solid black lines and gray lines represent the cases with δ=0.3 and δ=0.6 respectively. As we can see, with the similar gamma-ray flux at around 10 TeV, the gamma-ray flux at ò γ 100 TeV in the case with δ=0.3 is higher than that in the case with δ=0.6. Figures 3 and 4, gamma-rays from the GC would be detected by the CTA south with a 50 hr exposure. The observation of the GC by the CTA in the future would further constrain our model. The sensitivity of LHAASO to the GC is also reduced because the GC is located outside the sensitive region of LHAASO. A detection of 10 gamma-ray events from the GC by LHAASO will lead to a 5σ discovery, due to the zero background at 100 TeV (Bai et al. 2019). To accumulate more 100 TeV gamma-ray photons, it is better to count gamma-rays from a more extended area instead of from the small annulus region. The gamma-ray spectra from Area A 1°. 8 for four optimal cases, and the sensitivities of LHAASO detecting 10 events with 1 yr observation and 1 month observation, are plotted in Figure 5. The sensitivity of LHAASO for the GC is adopted from the LHAASO's white paper (Bai et al. 2019). According to the white paper, simulations on the Tibet AS Gamma experiment detecting high-energy gamma-ray photons from the GC are done to derive its effective area. Since the geometry area of LHAASO is about 46 times larger than that of the Tibet AS Gamma experiment, the sensitivity of LHAASO for the GC is approximately 46 times better than the simulated sensitivity of the Tibet AS Gamma. Note the actual sensitivity of LHAASO for the GC would be  1 , T=4×10 5 yr, and δ=0.3. The black solid line denotes the overall gamma-ray spectrum. The gray line denotes the overall gamma-ray spectrum adopted from Figure 3 for comparison. different for the real data analysis, which might be updated in the future. From this figure, we can tell that LHAASO would be a promising telescope to observe 100 TeV photons from the GC, and the future observation on the GC by LHAASO will constrain our models. Though the GC is located at the edge of HAWC's sensitive region, HAWC is still sensitive to the highest energy gamma-rays from the GC (Abeysekara et al. 2017); thus, HAWC's future observation of gamma-ray photons with energy of 100 TeV from the GC can also help to constrain our model.

As plotted in
In Figure 6, we plot the distribution of 100 TeV muon neutrinos around the GC. One can see that the neutrino flux per solid angle is higher in the central region, due to the denser MC in the CMZ. The distribution of neutrinos extends to the region with radius larger than 1°.8, which might be observed as a disklike extended neutrino source around the GC. As we can see, the 100 TeV muon neutrino image for the case with δ=0.3 is brighter than that for the case of δ=0.6.
Adopting the same parameter sets as in Figures 3 and 4, we calculate the neutrino spectra from Area A 1°. 8 and Area A 7°. 0 , respectively, shown in Figure 7. We plot the conventional atmospheric muon neutrinos (blue dotted line), which are 1.07 times those predicted by Honda (2006), averaged over the zenith angle, and the upper limit of the prompt atmospheric muon neutrinos (green dotted line; IceCube Collaboration et al. 2017). The black points are the astrophysical muon neutrino flux derived from the analysis of high-energy starting events (HESEs), observed by IceCube over 6 yr, averaged over the full sky (IceCube Collaboration et al. 2017). Current observations on the whole sky show no significant anisotropy or clustering (IceCube Collaboration et al. 2017), indicating a dominant contribution from extragalactic origins. Therefore, from the direction of the GC, aside from neutrinos produced by our model, there is an astrophysical neutrino background contributed by extragalactic sources, with flux the same as the observed astrophysical muon neutrino flux averaged over the full sky. We plot the best fit of the isotropic astrophysical muon neutrino spectrum in the black dotted line, following the power-law spectrum f µ Collaboration et al. 2017). Here we assume that this spectrum can be extrapolated to the energy as low as 30 TeV.
The predicted muon neutrino spectra for the fast dissipation case (red solid lines) peak at higher energy than those for the slow dissipation case (purple dashed lines), because the protons are accelerated to a higher energy in the fast dissipation case. The neutrino flux per solid angle averaged over Area A 1°. 8 is higher than that averaged over Area A 7°. 0 , which is consistent with the neutrino distribution shown in Figure 6. The former dominates over the atmospheric and isotropic astrophysical background at the energy band of < <ñ m  30 TeV 1 PeV. The fast dissipation case with δ=0.3 contributes more PeV neutrinos than the other cases.

The Detectability and Prospect
The muon-neutrino-induced track-like events observed by IceCube can be reconstructed with subdegree angular resolution. From Figure 2, we can tell under the constraint of H.E.S.S. observation that we expect IceCube to detect more than one muon neutrino from Area A 1°. 8 around the GC for the optimal cases, including through-going track-like events and starting track-like events. The through-going track-like events can be reconstructed with smaller angular uncertainty. However, from the direction of the GC, the IceCube-detected through-going track-like events, which start outside the detector, are dominated by atmospheric muon backgrounds, resulting in a low discovery potential (Aartsen et al. 2014. If there is a way to remove background muons, searching for through-going muon neutrinos from the central region with radius of 1°. 8 is one promising method to test our model.   The starting track-like events with vertex contained inside the detector would reduce the atmospheric muon backgrounds significantly, though the effective detector volume is decreased consequently (Aartsen et al. 2016. Aside from starting track-like events, selecting cascade-like events with the neutrino interaction vertex occurring inside the detector also rejects down-going atmospheric muons efficiently. What's more, the self-veto method can veto an atmospheric neutrino that is accompanied by a detectable muon from the same cosmic-ray air shower, then the atmospheric neutrino background is significantly reduced (Argüelles et al. 2018). The effective area for detecting both starting track-like and cascade-like events for neutrinos of all three flavors from the GC (IceCube Collaboration 2013) is higher than that for detecting the starting track events , though the angular uncertainty of cascade-like HESEs is relatively large.
The IceCube Collaboration has released a sample of 82 HESEs with energy larger than 60 TeV using 6 yr of data, including both tracks and cascades (IceCube Collaboration et al. 2017), with expected atmospheric muon backgrounds of 25.2±7.3 events and expected atmospheric neutrino backgrounds of -+ 15.6 3.9 11.4 events. A cascade event with energy of about 1 PeV is reconstructed to point toward 1°. 2 from the GC with a median angular uncertainty of 13°.2 (IceCube Collaboration 2013).
In Table 1, we list expected counts of the neutrino background for starting track events and HESEs with energy larger than 30 TeV from Area A 1°. 8 and Area A 7°. 0 for IceCube 10 yr operation, including the atmospheric neutrinos and isotropic astrophysical neutrinos derived from the IceCube observations. A self-veto method is adopted to suppress the atmospheric neutrino background. The significance of discovery if one event is observed during 10, 20, and 30 IceCube 86 strings (IC86) equivalent year exposures is also listed. One can see that, for the IceCube 10 yr operation, the significance of the discovery is low if counting neutrinos from Area A 7°. 0 , as more neutrino backgrounds are included from the larger area. The counts of neutrino background are suppressed to be much lower than 1 for Area A 1°. 8 . For the same operation time, the probability of discovery if detecting one starting track neutrino from Area A 1°. 8 is higher than detecting one HESE, because there is less background detected for starting track events than HESEs. If one starting track event is detected for the IceCube 10 yr operation from Area 1°. 8, the significance of discovering a neutrino source is about 98%, while if one starting track is detected with a 40 IC86 equivalent year exposure, the significance is decreased to about 92%. If one HESE from Area 1°.8 is detected in IceCube 10 yr operation, the significance of discovery is about 91%, while if one HESE is detected with a 40 IC86 equivalent year exposure, the significance decreases to about 71%. This conclusion only depends on the background models, but not dependent on the neutrino source models in the GC.
In Table 2, we list the expected neutrino counts for our optimal models of the neutrino source in the GC, comparing with the expected background neutrinos. Our optimal models predict 0.5-1.4 HESEs with energy larger than 30 TeV from Area A 1°. 8 for IceCube 10 yr operation, which are not violated by the current observation. From Table 2, we can see that the lower diffusion index case, i.e., δ=0.3, contributes more highenergy neutrinos than the higher diffusion index case. This is Note. Self-veto on the atmospheric background neutrinos for HESEs is adopted. P 10 , P 20 , P 30 , and P 40 are the significance of discovery if one event is observed during 10, 20, 30, and 40 IC86 equivalent year exposures.  because high-energy protons diffuse slower for a lower diffusion index. If we count neutrinos from a larger region around the GC, for instance, Area A 7°. 0 , the ratio of the expected signals to the background is smaller, making it more difficult to discover the source. For starting track-like events with energy above 30 TeV from Area A 1°. 8 observed by IceCube 10 yr observations, the background is suppressed to be much less than 1, and the signal-to-noise ratio is larger than 10. The expected neutrino counts of starting tracks and HESEs from Area 1°. 8 as a linear function of equivalent exposure time are plotted in Figure 8. A 16 to 37 IceCube equivalent year exposure is needed to detect one starting track event for the four optimal cases listed in Table 2. As listed in Table 1, if one starting track event is detected within 40 IceCube equivalent year exposure, the significance of discovery is larger than 90%. The effective area of the future detector IceCube-Gen2 on muon neutrinos is expected to be a factor of 2 times larger than that of IceCube for events of energy 10 TeV, and a factor of more than 3 times larger for events with energy above 100 TeV to the GC direction (van Santen 2018). Therefore, with the current 10 yr operation of IceCube plus a few more years of future operation of IceCube and IceCube-Gen2, one starting track event is expected for the optimal cases. The observation then would make a discovery of significance larger than 90% if one starting track is observed, or constrain our models if no starting track is observed.
About 7-20 IC86 equivalent year exposure is needed to detect one HESE, for the four optimal cases listed in Table 2. As mentioned earlier in this section, the IceCube reported one HESE pointing toward Area A 1°. 8 with 6 yr operation; however, the angular uncertainty of this cascade event is very large, thus, from current observations of HESEs, one cannot announce any discovery. A larger exposure accumulating more HESEs from the GC, or an improved angular resolution of HESEs, is needed to make a significant discovery. IceCube-Gen2 phase I (The IceCube-PINGU Collaboration 2014) will improve the angular resolution for track-like and cascade-like events in the near future, and thus would help to select HESEs from the GC.
The track-like HESE sample is overlapped with the starting track sample, according to Aartsen et al. (2016). Future analysis combining the starting track events and cascade-like HESEs will enhance the discovery potential.
ANTARES can detect up-going track events from the GC with atmospheric muons absorbed by Earth, with effective area smaller than that of IceCube for through-going events by a factor of 10 (Adrián- Martínez et al. 2012). The combined search with the ANTARES and IceCube neutrino telescopes for neutrino point sources in the southern sky would improve the discovery potential (Adrián-Martínez et al. 2016). The nextgeneration neutrino telescope KM3Net can also detect upgoing track events with no atmospheric muon backgrounds, and the effective area of KM3NeT-ARCA will be about 10 times larger than IceCube (Celli et al. 2017). Therefore, the future observation of the GC by combining KM3NeT-ARCA and IceCube is important to constrain our model.
In Pinat & Sánchez (2018), adopting both the spatial and energy distributions of backgrounds and signals, the IceCube Collaboration searched for extended sources using 7 yr of though-going track events, where the spatial distribution of signals is assumed to be a two-dimensional Gaussian distribution, and the energy distribution of signals is assumed to be following a flat spectrum. However, this analysis is sensitive to PeV-EeV neutrinos in the southern sky, due to the need to reject atmospheric muon backgrounds from through-going track events. The sensitivity of this analysis for a 2°source in the direction of the GC is about ---10 TeV cm s 11 2 1 , while our predicted neutrino flux for a 1°. 8 source at the GC is about ---8 10 TeV cm s 13 2 1 at 1 PeV for the optimal case. Therefore, it is not surprising that the result of extendedsource searches is consistent with the background-only hypothesis. The extended-source searches combining the through-going track event data with the starting track event data and cascade-like HESEs would lower the threshold energy and improve the sensitivity of the search . Our work presents the spatial and energy contribution to the probability distribution function of the signal events, as shown in Figures 6 and 7, which can be used in future searches for extended sources at the GC.

Conclusions and Discussions
Similar to the SFGs/SBGs, star-forming regions in the Galaxy might be the host of CR accelerators and contain rich gas as the target of inelastic pp collision. Their contribution to Figure 8. The expected counts of starting tracks and HESEs with energy larger than 30 TeV, from Area A 1°. 8 , as a function of exposure time in IC86 equivalent years, for optimal cases. Self-veto on the atmospheric background neutrinos for HESEs is adopted. The other parameters are the same as in Table 2. neutrinos might be resolvable, as they are much closer to Earth than distant SFGs/SBGs. We estimate the contribution of an HN associated with an MC complex in the star-forming region in our galaxy to neutrinos, and predict that IceCube might detect neutrinos with energy around 100 TeV from a typical galactic star-forming region in 10 yr of operation, if the source is located in the northern sky, and injected sub-PeV to PeV protons with total energy of 5×10 50 erg within a distance of 10 kpc. The corresponding 100 TeV gamma-ray emission might be detected by CTA, HAWC, LHAASO, and the Tibet AS array, if the source is located in the field of view of those telescopes. To be noted here, the total energy of protons injected into MCs, the injected time, the diffusion coefficient of the injected protons, the source distance, and the mass of MCs at the star-forming regions in the Galaxy are uncertain. The future observation by CTA, HAWC, and LHAASO would be able to constrain those uncertainties.
Furthermore, we assume that an HN exploded in the past at the GC, and simulate proton acceleration and confinement in the HNR. Protons that escape from the HNR will interact with MCs when diffusing around the GC region, then produce gamma-rays and neutrinos. Under the constraint from the H.E.S.S. Collaboration's observation on the high-energy gamma-rays from the GC region, we propose that LHAASO can discover the GC source at the energy of 100 TeV via a onemonth observation for optimal cases in the future. CTA South would help constrain our model.
We predict that more than one muon neutrino with energy higher than 30 TeV can be observed from the central region with a radius of 1°.8, by IceCube 10 yr operation, for optimal situations. However, atmospheric muons from the direction of the GC would dominate over the neutrino-induced events. If there is a way to remove background muons, searching for through-going muon neutrinos from the central region with a radius of 1°. 8 is one promising method to test our model. The through-going muon neutrino events detected by ANTARES and the future KM3NeT-ARCA from the direction of the GC contain no background muon contaminations, and thus will help to constrain our models.
If selecting starting track events from Area A 1°. 8 , atmospheric muons can be reduced significantly, and the background muon neutrinos are reduced to be lower than 1. However, the method of selecting starting track events reduces the volume of the detector. About 16 to 37 IceCube equivalent year exposure is needed to detect one starting muon neutrino with energy higher than 30 TeV from Area A 1°. 8 , for the optimal situations.
The third way is selecting both cascade-like and track-like HESEs for neutrinos of all flavors, with the neutrino interaction vertex contained inside the detector, adopting the self-veto method. In this search, the atmospheric backgrounds can be reduced significantly, but the angular uncertainty for the cascade-like HESEs is much larger than the track-like events. So far, a cascade HESE with energy of about 1 PeV observed by IceCube is reconstructed to point toward 1°.2 from the GC with a large angular uncertainty. No discovery can be announced. With the help of IceCube-Gen2 Phase I, the angular resolution will be improved. Together with a larger exposure provided by IceCube and IceCube-Gen2, a tighter constraint on our models via the observation of HESEs would be provided in the future.
A further sophisticated method to test or constrain our models is searching for an extended source (IceCube Collaboration et al. 2017) around the GC, combining through-going track events, starting track events, and cascade-like HESEs by adopting our predicted spatial templates and spectra of neutrinos.
In Section 3.1, adopting the H.E.S.S. observations on 10 TeV gamma-rays, we constrain the diffusion radius of 100 TeV protons to be larger than . Under this constraint, a larger diffusion coefficient is needed for a younger HN, according to Equation (1). On the other hand, if the total energy of the injected protons is lower, the constrained minimal diffusion radius of protons would be smaller. However, an HN or SN with lower injected energy might not be able to accelerate protons to energies higher than PeV, according to the simulation in Section 3.2. Therefore, according to our model, if gamma-rays and neutrinos with energy 100 TeV exist, the scale of the gamma-ray and neutrino-emitting region should be larger than a few degrees, to avoid overshooting the H.E.S.S. observation.
In this paper, we assumed an isotropic diffusion with the diffusion coefficient of protons depending on the proton energy following where the value of the diffusion index δ is still highly uncertain. Different δ would lead to different contributions to high-energy photons and neutrinos. Adopting the same D o , a smaller index implies that protons with energy ò p >100 TeV will diffuse more slowly, compared to the case for a larger diffusion index, thus there are more high-energy protons contained in the central region. As a result, more high-energy gamma-rays and neutrinos will be produced from the central region if the diffusion index is smaller.
What's more, as we cannot resolve the MCs in the CMZ, we assume a disk-like molecular zone in the central region. Future detections on MCs around the GC with a better angular resolution can help us to get a more precise distribution of neutrinos and gamma-rays.
Another uncertainty in our simulations is the ambient medium surrounding the exploding HN. We assume a homogeneous ambient medium surrounding the exploding HN. If the ambient medium is very dense due to the early mass loss of the star before the explosion, the kinetic energy of the ejecta will dissipate into the accelerated protons faster. If the ambient medium is similar to the ISM, the dissipation of the kinetic energy is slower. We discuss both of these situations and mark them as "fast dissipation case" and "slow dissipation case," respectively. In the fast dissipation case, the HNR is able to accelerate more protons to energies larger than PeV, leading to more high-energy neutrinos and gamma-ray photons.
The future observation by CTA South and LHAASO of gamma-ray emissions with energy larger than 10 TeV from the GC region, and the future detection of neutrino emission from the GC region by IceCube will further constrain our models and the dependence of the diffusion coefficient on the proton energy.
The H.E.S.S. observation of >10 TeV gamma-ray photons indicate a possible PeVatron in the GC. In this paper, we assumed the PeVatron to be an HNR formed in the past. As plotted in Figures 3 and 4, our simulation can only explain the H.E.S.S.-observed spectrum of the annulus region above 1 TeV, as it is difficult for protons with energy lower than 10 TeV to escape from the HNR, as shown in Figure 1. The gamma-ray photons with energy lower than 1 TeV might be explained by the interaction between protons with energy lower than 10 TeV and MCs. Protons with energy lower than 10 TeV around the GC region can be produced by other sources such as SNRs (HESS Collaboration et al. 2016) or the stellar winds of the compact stellar clusters (Crocker et al. 2011). Because the acceleration of PeV protons by individual SNRs or stellar clusters cannot last much longer than 100 yr (HESS Collaboration et al. 2016), a continuous PeVatron is needed to explain high-energy gamma-ray photons, for example, the HNR in our simulation, which is able to accelerate PeV protons with a duration longer than 10 3 yr. Other possible accelerators, such as the TDE of Sgr A * (Liu et al. 2016) and the RIAF of LLAGNs (Fujita et al. 2015), also suggest past activity in the GC, with uncertain maximum accelerated energy and total injected energy of protons, which can be tuned to satisfy the constraint by the current observations of H.E.S.S. The TDE of Sgr A * provides a fading proton injection (Liu et al. 2016), the RIAF of LLAGN provides an instantaneous proton injection (Fujita et al. 2015), while the proton injection from the HNR according to our simulation is plotted in Figure 1. The different evolution of the proton injections for the three accelerators indicates different distributions of protons, gammarays, and neutrinos around the GC region. Future observations, which can derive the gamma-ray and neutrino distribution around the GC with higher angular resolution, can derive the CR distribution around the GC, and can further distinguish among those accelerator candidates.