Does or Did the Supernova Remnant Cassiopeia A Operate as a PeVatron?

For decades, supernova remnants (SNRs) have been considered the prime sources of Galactic cosmic rays (CRs). But whether SNRs can accelerate CR protons to PeV energies and thus dominate CR flux up to the knee is currently under intensive theoretical and phenomenological debate. The direct test of the ability of SNRs to operate as CR PeVatrons can be provided by ultrahigh-energy (UHE; E γ ≥ 100 TeV) γ-rays. In this context, the historical SNR Cassiopeia A (Cas A) is considered one of the most promising targets for UHE observations. This paper presents the observation of Cas A and its vicinity by the LHAASO KM2A detector. The exceptional sensitivity of LHAASO KM2A in the UHE band, combined with the young age of Cas A, enabled us to derive stringent model-independent limits on the energy budget of UHE protons and nuclei accelerated by Cas A at any epoch after the explosion. The results challenge the prevailing paradigm that Cas A–type SNRs are major suppliers of PeV CRs in the Milky Way.


INTRODUCTION
Within the current Galactic CR paradigm, SNRs are considered the major contributors to the observed CR flux.
For decades, this conviction has been based on phenomenological arguments and supported by the Diffusive Shock Acceleration (DSA) theory (for a review, see Malkov & Drury (2001)).The dominance of SNRs' contribution up to the "knee", a distinct break in the CR spectrum around 3-4 PeV, would imply that at certain epochs of their evolution, SNRs should operate as PeVatrons.However, as has been recognised long ago (Lagage & Cesarsky 1983), within the framework of the standard DSA applied to young SNRs, the energy of accelerated protons could hardly exceed 100 TeV.A possible solution was proposed by Bell (2004) by amplifying the magnetic field upstream of the shock through the instabilities driven by CRs.Combined with two other critical conditions, (1) shock speed of thousands of km/s, and (2) particle diffusion in the extreme (Bohm) regime, enhancing the magnetic field to B ≥ 100µG allows acceleration of protons and nuclei to 1 PeV/nuc.The interactions of the shock accelerated protons and nuclei with the ambient gas make young SNRs potentially detectable sources of γ-rays and neutrinos.Therefore, the detection of TeV γ-rays from more than a dozen young SNRs1 is a direct proof of effective acceleration of protons and/or electrons to very high energies.In general, the available data do not allow us to distinguish between the hadronic ("π 0 -decay") and leptonic ("inverse Compton") origin of TeV γ-rays.Even ignoring this ambiguity, i.e., assuming that γ-rays are produced by accelerated protons, we face a problem with the explanation of TeV γ-ray measurements.The reported steep γ-ray spectra with a differential power-law index Γ ∼ 2.4 − 2.7 (e.g., as summerised in the supplementary of Aharonian et al. 2019) do not readily agree with predictions of the standard DSA theory for strong shocks of young SNRs.This result can be interpreted in three different ways: (1) A steep power-law proton spectrum is mimicked by the result of the combination of a hard (E −2 type) power-law spectrum and an "early" exponential cutoff, E 0 ≪ 100 TeV.This simple explanation supports the predictions of the standard DSA and the conclusion of Lagage & Cesarsky (1983) for the case of a non-amplified magnetic field.This mathematical trick works in a relatively narrow (one decade or so) energy interval around and below the cutoff energy.If this explanation is correct, at higher energies the measurements should reveal a gradual spectral steepening in the deep exponential cutoff region.
(2) The claim that strong shock acceleration results in a hard acceleration spectrum is not always correct.Indeed, steep power-law proton spectra can be formed over a wide energy interval in certain realistic environments and scenarios (see, e.g., Bell et al. 2011Bell et al. , 2019;;Malkov & Aharonian 2019;Hanusch et al. 2019;Caprioli et al. 2020;Xu & Lazarian 2022).The steep power-law spectra can extend over a wide energy interval, formally up to 1 PeV, assuming that such energetic protons are effectively confined in the shell.Correspondingly, the resulting steep power-law γ-ray spectra could continue up to 100 TeV.The detection of steep multi-TeV γ-ray spectra in both above scenarios is a rather difficult but feasible task for detectors like CTA and LHAASO.
(3) Acceleration of protons to PeV energies at the very early (≪ 100 yr) epoch of the SNR evolution is preferable (Bell et al. 2013;Zirakashvili et al. 2014), in particular, because of the shock speed exceeding 10,000 km/s and the high density of the progenitor wind.However, the CRs with highest energy have already escaped from the shock upstream and thus the remnant is currently emptied of multi-TeV CRs.This can naturally explain the steep spectra of γ-rays from SNRs, including from very young ones, in particular from SNR Cassiopeia A (Aharonian et al. 2001;Ahnen et al. 2017;Acciari et al. 2010) and Tycho (Acciari et al. 2011).
The γ-ray spectrum in scenario (2) is predicted to have a steeper power law shape without cutoff, which is significantly different from the other two scenarios.The γ-ray spectra from the remnants in scenarios (1) and (3) could be similar, although for different reasons.On the other hand, scenario (3) differs principally from scenario (1) by the γ-radiation outside the remnant.While in scenario (1), PeV particles are not produced at all, scenario (3) allows proton acceleration up to 1 PeV, but assumes that presently the highest energy particles already left the remnant and presently occupy regions beyond the shell.Any outcome of γ-ray measurements at tens to hundreds TeV performed with detectors of adequate sensitivity -either detection of positive signals or upper limits -would provide definite answers to whether specific SNRs could accelerate protons 1 PeV at any epoch of their evolution.
The detection or upper limits from SNRs of different types and ages could finally establish whether SNRs, as a source population, can be responsible for the CR fluxes up to the "knee".The results would be relevant to both scenarios (2) and (3).
The historical SNR Cas A represents a key target for such studies.It has a well-defined age and distance.It also belongs to a class of core-collapse supernova that explodes in a dense red supergiant wind.Such systems are believed to be able to accelerate ions to the highest energy due to the efficient amplification of magnetic fields by CR streaming in dense environments.
In this paper, using LHAASO KM2A data, we derive γ-ray flux upper limits towards Cas A and test its ability to act as a PeVatron.In Sec.2, we discuss the status of previous gamma-ray observations of Cas A and the theoretical challenges of acceleration of protons to ultrahigh energies at different epochs of its evolution.In Sec.3, we review the physical properties of Cas A, including the gas distribution in proximity to the remnant.In Sec.4,we describe the procedure of extraction of flux upper limits based on the LHAASO KM2A data.In Sec.5, we discuss the astrophysical implications of the obtained results.

BASIC FACTS ABOUT CAS A
Cassiopeia A (Cas A, SNR G111.7-02.1) is one of the youngest SNRs in our Galaxy.This ≈ 340 year old (Reed et al. 1995) structure is the remnant of a type IIb supernova (Krause et al. 2008) located at a distance of 3.4 kpc.It is one of the brightest galactic radio sources with a shell structure of angular radius 2.5' corresponding to the physical size of 2.5 pc (Kassim et al. 1995).The synchrotron radiation extends from radio (Tuffs et al. 1997) to hard X-rays of about 100 keV (Grefenstette et al. 2015).Fermi LAT observations reveal a hint of a hadronic origin of the γ-ray emissions (Yuan et al. 2013;Zirakashvili et al. 2014).TeV γ-ray emission from Cas A has been discovered by the HEGRA IACT array (Aharonian et al. 2001), and later confirmed by the MAGIC (Ahnen et al. 2017) and VERITAS (Acciari et al. 2010) collaborations.
The hadronic models of GeV-to-TeV radiation require a very high, but still acceptable, acceleration efficiency of about 25% (Zirakashvili et al. 2014), while in the two-zone models the acceleration efficiency required can be lower (Liu et al. 2022).The significant steepening or a cutoff in the spectrum reported around a few TeV by Ahnen et al. (2017) and Abeysekara et al. ( 2020) implies a corresponding steepening in the spectrum of parent protons at energies of tens of TeV.This constrains but does not exclude the presence of PeV protons in the nebula or outside the shell.
Due to the limited age of Cas A, even if PeV protons have been accelerated at the early (≪ 100 years) epoch of the SNR, they could not propagate too far from the SNR.Thus, for robust conclusions, one should probe both the SNR itself and the ≈ 100 pc environment surrounding Cas A, for γ-rays of energy exceeding 100 TeV.

GAS AND CRS IN THE VICINITY OF CAS A
Since we are searching for γ-rays produced by PeV CR protons through interactions with the ambient gas, the information about the gas is crucial for our studies.Cas A is located within a dense gas environment.Zhou et al. (2018) have detected molecular clouds of total mass 200 M ⊙ in front of Cas A. Within 200 pc around Cas A, Kim et al. (2008) has derived the average gas density of about 3 cm −3 by investigating the infrared emission.In this region, compact molecular clumps with a density as high as 10 5 cm −3 have also been detected (Reynoso & Goss 2002).Ma et al. (2019) have performed a detailed J = 1-0 survey of CO in a large field near Cas A using the Purple Mountain Observatory (PMO) 13.7 m millimetre telescope.The total mass of molecular gas was ≈ 9.5 × 10 5 M ⊙ , which corresponds to an average density of 10 cm −3 within 100 pc proximity of Cas A. The line profiles have asymmetric or broadened shapes, which indicates a possible interaction between the SNR and the molecular clouds.In this work, we used the J = 1-0 12 CO data as described in Ma et al. (2019).The derived molecular gas column densities are shown in Figure .1 by assuming 2.0 × 10 20 (K km s −1 ) −1 for the CO − to − H 2 conversion factor (Bolatto et al. 2013).
As mentioned above, the average gas density is estimated as 10 cm −3 in the vicinity, so to estimate the gas mass we need to calculate the volume that is occupied by CRs.The latter depends on the CR propagation in the vicinity of Cas A. Since we deal with relatively small spatial scales and ultra-high energies, we cannot a priori assume that the propagation proceeds in the diffusive regime.One cannot exclude that close to the accelerator, highest energy CRs propagate ballistically near the source and enter the diffusive regime when r > D/c, where r is the distance from the source, D is the diffusion coefficient and c is the speed of light.The diffusion coefficient in the Interstellar Medium (ISM) is estimated as D(E) ∼ 2 × 10 28 ( E 1 GeV ) 0.5 cm 2 /s (Yuan et al. 2017).Its extrapolation to E = 1 PeV gives D ∼ 2 × 10 31 cm 2 /s.This implies r D/c ≃ 200 pc.For comparison, rectilinear propagation length of protons over 340 years is ≃ 110 pc.Thus, even in the case of acceleration of PeV protons at the early epochs of the SNR, they would still continue to move ballistically.Due to the beaming effects, only the protons propagating towards us can produce visible γ-rays.Then the γ-ray source should have a similar angular size as the accelerator, i.e., the SNR shell.This is about 2.5 ′ , so it can be detected by LHAASO-KM2A as a point-like source.
The diffusion of CRs in the vicinity of the accelerators could be strongly suppressed due to enhanced turbulence introduced, e.g., by the streaming instability of injected CRs (see e.g.(Malkov et al. 2013)).The interaction of the cosmic ray precursor and dense clumps can also drive turbulence.By either following the turbulent tangling field lines or mirror diffusion in compressible magnetic fluctuations (Lazarian & Xu 2021), high-energy CR protons might undergo very slow diffusion.In this case, PeV CRs could propagate in the diffusive regime and the region occupied by VHE CRs would be of the order r ∼ √ 2DT ∼ 20 pc( χ 10 −2 ) 0.5 , where χ is the suppression factor of the diffusion coefficient near the source compared with the Galactic diffusion coefficient, D source (E) = χD(E).The angular size of the γ-ray emission is estimated as r/d ∼ 0.33 • ( χ 10 −2 ) 0.5 .The maximum diffusion length, in any case, should be smaller than 110 pc, the ballistic propagation length over 340 years.And the maximum angular extension of the γ-ray emission when CRs are in the diffusive regime should be smaller than 1.8 • .To estimate the total CR energy budget from the γ-ray luminosity, only the average density matters.As estimated in Ma et al. (2019), the average density in the nearby 100 pc of Cas A is 10 cm −3 .In this case, due to the lower average density and larger integration area in deriving the γ-ray flux, the derived CR energy budget would be more conservative than the ballistic propagation case we considered above.To derive the robust and model independent estimations we consider only the diffusive regime in the following discussion.

KM2A DATA ANALYSIS
The data used in this analysis were collected from December 2019 to September 2022 by the half KM2A, three quarters KM2A and full KM2A.Quality selections have been applied to ensure a reliable data quality.We select the periods when at least 95% of working sub-detectors were in good condition.The total effective observation time is 932.74 days after data quality selection.To achieve a better performance of the array, several event cuts are also implemented.The cuts on data and event selections are the same as listed in the KM2A performance paper (Aharonian et al. 2021).
Considering the energy resolution and statistics, one decade of energy is divided into 5 bins with a bin width of log 10 E = 0.2.The sky in celestial coordinates (right ascension and declination) is divided into cells of size 0.1 • × 0.1 • and filled with detected events according to their reconstructed arrival directions for each energy bin.To extract the excess of γ-rays from each cell, the "direct integration method" (Fleysher et al. 2004) is adopted to estimate the number of cosmic ray background events.The test statistic (TS) used to evaluate the significance of the source in this work is TS = 2 log(λ), where λ = L 1 /L 0 , L 1 is the maximum likelihood value for the alternative hypothesis, and L 0 is the maximum likelihood value of the null hypothesis.The TS map of the ROI is shown in the left panel of Fig. 1.We found that there is no significant excess in the vicinity of Cas A.
As mentioned in the last section, we consider both the ballistic propagation and diffusive regime.In the former case, the expected γ-ray emission produced by CRs escaping from Cas A is point-like.In KM2A analysis, taking into account the point spread function (PSF), we integrate all photons in the 90% contamination radius from a point source.
For diffusive regime, the intrinsic size of the source varies with the diffusion coefficient; however, the maximum angular radius is 1.8 • as estimated in the previous section.Taking into account of the angular resolution of the detector, we use an integration radius of 2.2 • for the lower energy bin and 1.85 • for energies higher than 100 TeV, to estimate the upper limits.No significant excess has been detected.We derived an upper limit for each energy bin at the 95% confidence level using the method proposed by Helene (1983).For a more conservative estimate of the upper limit, we use the CERN ROOT tools following Feldman & Cousins (1998) to calculate the 95% upper limit of the signal when the background count is less than 20.The results are shown in Fig. 2 and Fig. 3.The main uncertainties come from the uncertainty in the differential spectral index of the γ-ray emission, which we assume to derive the upper limit.In order to estimate such systematic uncertainties, the analysis is separately performed for the γ-ray spectral index of -2 or -3.The uncertainties found from this test are small: about 10% and 3% for the first two energy bins, and less than 1% at higher energies.
The systematic uncertainties affecting the flux upper limits are similar to those studied in Aharonian et al. (2021).
The number of operating units varied with time due to debugging a few percent of detector units.The variable layout of the array, affecting the γ-ray/background separation, has a slight effect on the flux.The systematic uncertainty also comes from the atmospheric model used in the simulation, which affects the detection efficiency.The overall systematic uncertainty affecting the flux is estimated to be ∼7%.
In principle, the diffuse galactic γ-ray emission (DGE) should be subtracted before we derive the upper limit for sources.But the inclusion of such a component in deriving the upper limit only has a marginal (10%) impact on the   The UHE γ-ray flux upper limits within 1.8 degree from Cas A. The black, red and blue curves correspond to the predicted pion-decay γ-ray fluxes for the proton index of 2.0, 2.4 and 2.7, respectively.The total energies of CR protons within 100 − 1000 TeV are 5 × 10 49 erg, 3 × 10 48 erg and 1 × 10 47 erg, respectively, which are the CR energy budgets expected if Cas A type SNR produce all the PeV CRs in our Galaxy.
final results.Thus to avoid further systematic errors, we do not subtract the DGE in this work and use the most conservative upper limits.
the possibility that all CRs above 100 TeV in our Galaxy are injected from Cas A type objects.As an alternative way to express this result, we plotted the predicted γ-ray flux in Figure 3 assuming the CRs injected by Cas A type SNRs dominated the PeV CRs in our Galaxy.We again see that the limits from LHAASO KM2A rule out that Cas A type SNRs produce all PeV CRs in the Galaxy for γ = 2.0 or 2.4,but do not rule out if γ = 2.7.Indeed, as shown in the last paragraph, the upper limit of the rate of CR above 100 TeV injected by Cas A only slightly depends on the injection spectrum and is estimated as ∼ 1.3 × 10 38 (10 cm −3 /n) erg/s.On the other hand, the total energy budget required for CRs above 100 TeV depends only on the injection spectrum.If we assume a pure power law, the injection energy budget above 100 TeV is W (E > 100 TeV) = (1 − 3) × 10 (41−5(Γ−2)) erg/s for Γ > 2.0, where Γ is the injection spectral index.By comparing with the limit set by KM2A observations on Cas A, we found that for n ≃ 10 cm −3 as estimated from CO observations, our results can rule out the hypothesis that Cas A type SNRs are the only PeVatrons contributing to CRs above 100 TeV for a CR injection spectrum index less than 2.5.A smaller gas density would relax the constraints.But note that an average gas density well below 10 cm 3 seems unrealistic in this region.In the most conservative case, assuming a gas density of 1 cm −3 , the lower limit for the injection spectral index is 2.4.
In scenario (3), the highest energy CRs have already escaped from the shock upstream.Cas A still lies in its progenitor's wind bubble (Weil et al. 2020), which has a size of about 10 pc and a much lower density compared with the ISM.In the most extreme case, it is possible that CRs escaped from the shell and confined in the bubble.The angular size of the bubble is about 5 arcminutes, which is much smaller than the PSF of KM2A.As a result, the bubble is also a point source in KM2A, and the upper limit derived from a point source is about 1/5 of the extended source scenario considered above.If we assume an average density of 0.1 cm −3 in the wind bubble, the constraints on the total CR energy budget can be 20 times larger than our estimation above.
In the derivation of the γ-ray upper limit, the Galactic diffuse γ-ray emission is not considered, which would make the upper limit even lower and the constraint more stringent.Thus we argue that our calculation implies that the Cas A-type SNRs cannot be the only sources for CRs up to PeV in the regime of standard DSA (hard injection spectrum).
If we assume a softer injection spectrum, Cas A type SNRs can still be major PeVatrons in our Galaxy.But such a softer injection spectrum would be in tension with most CR propagation models, in which a strong energy dependence of the CR diffusion coefficient is used, and the CR injection spectrum should be hard to fit the locally observed CR spectrum.
We also note that a more complex spectral shape, such as a broken power law, is possible for the total CR flux injected by Cas A. But in this work our constraints mainly come from a rather narrow energy band, and the constraints on the exact spectral shape are limited.To avoid uncertainties and degeneracy, we used a single power law spectra in this study.
In conclusion, based on the LHAASO KM2A measurements, we set stringent upper limits on the total energy budget of UHE protons accelerated by Cas A. Our constraints here are nearly model independent, unless in the hypothetical, and in our view, not very realistic situation, when cosmic rays are confined in a compact low density wind bubble after escaping from the shell.Although we cannot formally rule out the possibility that the Cas A type SNRs are the Galactic PeVatrons, the stringent limits on the accelerated (injection) spectrum have important implications.Further UHE observations will investigate more parameter spaces.We note that the most stringent limits come from the γ-ray fluxes above 100 TeV, for which the sensitivity of LHAASO KM2A will remain unrivalled in the foreseeable future.
We anticipate that the accumulation of LHAASO KM2A data for another 5 to 10 years will improve the upper limit substantially and would finally rule out or confirm Cas A type SNRs as a PeVatron population.

Figure 1 :
Figure 1: left panel: The LHAASO-KM2A TS map above 25 TeV in the 4 • × 4 • region around Cas A, the colorbar shows the TS value.right panel:The gas column density near Cas A derived from 12 CO line observations (Ma et al. 2019).The green circle labels the position of Cas A.

Figure 2 :
Figure2: The UHE γ-ray flux upper limits assuming Cas A is a point source.The data points are from MAGIC and VERITAS obervations.The grey curve is the power law function with the index of 2.8.A harder spectrum would violate the LHAASO upper limit.

Figure 3 :
Figure3: The UHE γ-ray flux upper limits within 1.8 degree from Cas A. The black, red and blue curves correspond to the predicted pion-decay γ-ray fluxes for the proton index of 2.0, 2.4 and 2.7, respectively.The total energies of CR protons within 100 − 1000 TeV are 5 × 10 49 erg, 3 × 10 48 erg and 1 × 10 47 erg, respectively, which are the CR energy budgets expected if Cas A type SNR produce all the PeV CRs in our Galaxy.