Electron Acceleration in Middle-age Shell-type γ-Ray Supernova Remnants (2024)

Since it was first proposed that supernovae may be responsible for the production of cosmic rays (CRs) in the 1930s (Baade & Zwicky 1934), it has been widely accepted that CRs up to the spectral knee energy of ∼3×10¹⁵ eV mainly originate in the Milky Way, presumably being accelerated by shocks of supernova remnants (SNRs). Now we know that energetic electrons and protons can indeed be accelerated in the shock waves of SNRs. Evidence for acceleration of electrons first comes from the SNRs' radio continuum spectra, which are produced by relativistic electrons through the synchrotron process. Moreover, the detection of nonthermal X-ray emission in the SNR SN1006 shows that SNR shocks can boost electrons up to TeV energies (Koyama et al. 1995). With regard to protons, evidence of their acceleration in SNRs was confirmed recently by the discovery of the pion-decay γ-ray signal in SNR IC443 (Ackermann et al. 2013), W44 (Giuliani et al. 2011; Ackermann et al. 2013; Cardillo et al. 2014) and W51C (Jogler & Funk 2016).

The presence of both energetic electrons and protons in SNRs makes the origin of their observed γ-ray emission ambiguous, because both leptonic and hadronic processes can generate γ-ray emission efficiently. In the leptonic scenario, the γ-rays are produced via inverse Compton (IC) scattering on the background photon field by relativistic electrons and/or via the bremsstrahlung process. In hadronic models, the γ-rays come from the decay of π⁰ produced mainly through inelastic collisions between relativistic protons and ambient nuclei. Based on these two scenarios, various models have been proposed to explain the observed γ-ray emissions from individual sources (e.g., Li & Chen 2010; Uchiyama et al. 2010; Ellison & Bykov 2011; Ohira et al. 2011; Tang et al. 2011; Atoyan & Dermer 2012; Morlino & Caprioli 2012; Berezhko et al. 2013; Zhang et al. 2013; Gabici & Aharonian 2014; Tang & Chevalier 2014; Zhang & Chen 2016; Ohira & Yamazaki 2017).

With the development of ground-based γ-ray telescopes like H.E.S.S., VERITAS, MAGIC, etc., about two dozen TeV-bright SNRs⁴ have been detected. With the accumulation of data from the Fermi satellite, the first Fermi-LAT SNR catalog with identified SNRs and candidates has been built. These observations reveal diversified γ-ray spectra, which Yuan et al. (2012) attributed to variations in the density of the emission region. The increasing number of γ-ray-bright SNRs also makes it possible to perform a population study (Mandelartz & Becker Tjus 2015; Zeng et al. 2019).

Traditionally, the exploration of the SNR origin of galactic CRs has focused on the maximum energy to which shocks of SNRs can accelerate particles (Lagage & Cesarsky 1983). Although it has been well recognized that strong shocks of young SNRs likely dominate the acceleration of the highest energy CRs (Bell 2004), the nondetection of SNRs above 100 TeV challenges SNRs as PeVatrons. Moreover, recent deep TeV observations of the brightest TeV SNR RX J1713.7−3946 reveal a convex γ-ray spectrum, which can be attributed to a broken power-law distribution of energetic particles not anticipated by most particle acceleration models. Ohira & Yamazaki (2017) first proposed that such a broken power-law particle distribution may be attributed to the time integration of accelerated particles by a gradually weakening shock that accelerates more less-energetic particles as the shock slows down. Such a scenario is in line with a time-dependent particle acceleration model proposed by Zhang et al. (2017) for the anomalous energy distribution of CR ions.

Using a simple one-zone emission model, Zeng et al. (2019) fit the spectral energy distribution of nonthermal emission from a sample of 35 middle-age or old SNRs. They found that in general, the particle distribution in these SNRs can be described by a broken power-law function with a high-energy cutoff and as the SNR ages, the break energy decreases and the low-energy spectrum becomes harder, in agreement with the model proposed by Ohira & Yamazaki (2017) and implying very efficient acceleration of low-energy CRs in old SNRs (Zhang & Liu 2018). In this paper, we focus on the spectral modeling of seven shell-type bright γ-ray SNRs, namely SN1006, RXJ1713.7−3946, RCW86, RX J0852.0−4622, HESSJ1731−347, G150.3+4.5, and HESSJ1534−571. The basic properties of these SNRs are listed in Table 1. All of them have relatively harder GeV spectra and have spatially well-correlated nonthermal emissions from radio to γ-rays so that a simple one-zone leptonic emission model is applicable, and in general the high-energy electron distribution is assumed to follow a broken power law with a high-energy cutoff. The model is described in Section 2. In Section 3, we discuss spectral fits to each SNR. A unified model is proposed in Section 4 for the evolution of high-energy electron distribution in these SNRs. We draw conclusions in Section 5.

To study the overall spectral characteristics, we adopt a one-zone emission model. In this model, electrons accelerated at the shock front obey a single power-law energy distribution with a maximum energy E_max, and are injected into the main emission region. The maximum energy is determined by the balance of acceleration and energy loss in early stages of the SNR shock evolution. The acceleration timescale is given by , where η_acc≈10 is a numerical factor, u_s is the shock speed, and D(E, t) is the diffusion coefficient (Drury 1983). Following Ohira et al. (2012), we assume that the diffusion coefficient has a Bohm-like form , where c, e, E, and B(t) are the speed of light, the elementary charge, the particle energy, and the magnetic field, respectively. The gyrofactor η_g in the diffusion coefficient is still not well understood, and maybe changing with the SNR age (Aharonian et al. 2017). In this paper, we assume it has a power-law behavior η_g∝t^λ (e.g., Ohira et al. 2012 and references therein). Due to the dominance of synchrotron loss, the cooling time is given by τ_syn=1.25×10⁴(E/1 TeV)⁻¹(B_SNR/1 μG)⁻² kyr. By setting τ_acc=τ_syn, we can get the cooling-limited maximum energy . In very early stage of SNRs, the energy loss time is longer than the SNR age t and the maximum energy is constrained by t, which gives an age-limited maximum energy . Therefore, the maximum energy up to which relativistic electrons can be accelerated at the shock front is given by , where and are given, respectively, by

Electrons injected into the emission region are subject to radiative energy loss. The accumulated electron distribution may be described with a broken power-law function. In the Sedov phase of SNRs, shock speed decreases, and so does the . In combination with an efficient injection, the time-integrated electron distribution may also be described by a broken power-law function (Ohira & Yamazaki 2017).

Without loss of generality, we assume that the time-integrated electron spectrum in the emission zone has a broken power-law form with a high-energy cutoff:

where E_bre is the break energy, α is the index below the break, β characterizes the cutoff shape, σ describes the smoothness of the transition at the break energy with σ=5 adopted here, and Δα is the index change. The radio and X-ray emissions are produced by energetic electrons via synchrotron process while the γ-rays are produced via IC. We adopt the photon emissivity given in Blumenthal & Gould (1970) for synchrotron and IC processes. Both the hadronic and bremsstrahlung processes are ignored because the density of the emission region for these seven sources is very low.

To calculate the IC emission, besides the microwave background radiation, we include an IR photon background with a temperature of 35 K and an energy density of 0.7 eV cm⁻³ except specified otherwise below. In general, there are seven model parameters that need to be determined by the spectral fit: α, Δα, E_bre, E_cut, β, the total energy of electrons above 1 GeV W_e, which determines the normalization of the electron distribution, and B_SNR. To facilitate the comparison of these SNRs, we assume β=2 and Δα=0.8 except for SN1006, which is the youngest SNR in our sample. Justification for these assumptions will be given in Section 4. There are therefore five parameters.

3.1.SN1006 (G327.6+14.6)

SN 1006 is the remnant of a historical supernova that exploded on AD 1006 May 1 (Stephenson & Green 2002). It is a typical bilateral SNR with most nonthermal emissions originating from two large segments on opposite sides of the remnant. In the GeV band, only the northeast (NE) shell was detected by Fermi-LAT (Xing et al. 2016; Condon et al. 2017), showing a hard GeV spectrum. Due to lack of GeV data for the southwest shell, we only fit γ-ray data of the NE shell and the synchrotron spectrum is multiplied by a factor of 2 to match the radio and X-ray fluxes of the whole remnant. We find that a single power-law distribution with a high-energy cutoff can fit the nonthermal emission spectrum and for β=1, the energy loss time at the cutoff energy due to synchrotron radiation is comparable to the age of the remnant. The corresponding model parameters are given in Table 2. The spectral fit is shown in Figure 1.

Table 2.Fitting Parameters

Source	β	Δα	α	Ebre	Ecut	We	BSNR	τsyn(Ecut)
				(TeV)	(TeV)	(1047 erg)	(μG)	(kyr)
SN 1006	1.0	⋯	2.2	⋯	16	1.7	24	1.4	30
RX J1713.7−3946	2.0	0.8	2.0	3.5	80	4.0	16	0.6	86
RCW 86	2.0	0.8	2.2	2.5	50	9.0	26	0.4	113
RX J0852.0−4622	2.0	0.8	1.9	1.0	50	8.0	9	3.1	192
HESS J1731−347	2.0	0.8	1.9	3.0	30	2.8	26	0.6	99
G150.3+4.5	2.0	0.8	1.7	0.2	60	0.4	3	23.2	96
HESS J1534−571	2.0	0.8	1.3	1.0E−2	15	6.0	5	33.5	4731

Note.

^aCalculated from .

Download table as: ASCIITypeset image

3.2.RXJ1713.7−3946 (G347.3–0.5)

RXJ1713.7–3946 is a Galactic SNR discovered in X-rays by ROSAT (Pfeffermann & Aschenbach 1996) and has an accurate age measurement due to its association with the historical SN AD393 (Wang et al. 1997). It has extremely faint radio emission and prominent nonthermal X-ray emission, both showing shell morphology with enhancement in the northwest portion. Observations of the H.E.S.S. experiment first resolved TeV γ-rays from the shell, and its ringlike shape closely matches the nonthermal X-ray emission (Aharonian et al. 2007a). The GeV γ-ray spectrum has a photon index of 1.5 (Abdo et al. 2011) and can be smoothly connected with TeV data. The GeV–TeV combination analyses performed by the H.E.S.S. Collaboration et al. (2018a) with deeper observations show that there is a spectral break around 100 GeV. For Δα=0.8, β=2.0, and infrared photons with a temperature of 26.5 K and an energy density of 0.42 eVcm⁻³ (Finke & Dermer 2012), we obtain the following best-fit model parameters: α=2.0 W_e=4.0×10⁴⁷ erg, B_SNR=16 μG, E_bre=3.5 TeV, and E_cut=80 TeV.

3.3.RCW86 (G315.4−2.3)

RCW86 is a radio shell-type Galactic SNR with a nearly circular shape and the diameter is about 42' and the radio spectral index is about 0.6 (Caswell et al. 1975; Kesteven & Caswell 1987). Although the inferred young age of RCW86 seems to be inconsistent with the large physical size at a distance of ∼2.5 kpc (Rosado et al. 1996; Sollerman et al. 2003; Helder et al. 2013), it is likely the remnant of the first historical SN seen by the Chinese in AD 185 (e.g., Vink et al. 2006). At the GeV band, only upper limits were derived toward this source by Lemoine-Goumard et al. (2012), and γ-ray emission from a point-like source coincident with the position of RCW86 was subsequently reported by Yuan et al. (2014). With 6.5 yr data, a study performed by Ajello et al. (2016) showed that the γ-ray emission of RCW86 probably originates from a spatially extended ring region with a hard GeV spectrum with a photon index of Γ∼1.42. The spectrum extends into the TeV range with a high-energy cutoff as revealed by H.E.S.S. observations (Aharonian et al. 2009; H.E.S.S. Collaboration et al. 2018e). For Δα=0.8 and β=2.0, we obtain the following best-fit model parameters: α=2.2, W_e=9.0×10⁴⁷ erg, B_SNR=26 μG, E_bre=2.5 TeV, and E_cut=50 TeV.

3.4.RXJ0852.0−4622 (G266.2−1.2)

SNR RXJ0852.0−4622, also known as Vela Jr., overlapping the southeast corner of the Vela SNR, was discovered in X-rays by ROSAT satellite (Aschenbach 1998). It has faint radio emission with a spectral index of ∼0.3 and a rim-brightened morphology similar to that in the X-rays (Combi et al. 1999). Like radio and X-ray images, the TeV γ-ray emission also has an almost circular shape with rim enhancement, namely, showing the shell structure (Aharonian et al. 2007b). A detailed analysis of the H.E.S.S. data (H.E.S.S. Collaboration et al. 2018b) shows that the TeV γ-ray spectrum can smoothly connect to the hard spectrum in the GeV band detected by Fermi-LAT (Tanaka et al. 2011). In combination with the lack of thermal X-ray emission (e.g., Slane et al. 2001), these observations strongly favor a leptonic origin for the γ-ray emission. For Δα=0.8 and β=2.0, we obtain the following best-fit model parameters: W_e=6.5×10⁴⁷ erg, α=2, B_SNR=9 μG, E_bre=3.0 TeV, and E_cut=53 TeV.

3.5.HESSJ1731−347 (G353.6−0.7)

3.6.G150.3+4.5

The eastern segment of G150.3+4.5 was first reported by Gerbrandt et al. (2014, called G150.8+3.8) and was considered a strong SNR candidate due to its semicircular shape, nonthermal radio spectrum, and red optical filaments. Using Urumqi λ6 cm survey data, Gao & Han (2014) found a western and a northern shell beside the eastern shell, forming an enclosed oval shell structure with size of 25×30 and, combining with Effelsberg 11 and 21 cm data and CGPS 1420 and 408 MHz observations, they confirmed it to be an SNR. Based on the H i observations and the Galactic rotation curve, Cohen (2016) obtained three possible distances and suggested the nearest distance of 0.4kpc, due to its large angular size. Using this distance and the Sedov evolution model, the age is found to be between 0.5 and 5 kyr (Cohen 2016). Here we adopt 3 kyr as its age and derive a velocity of ∼1500 km s⁻¹. At high-energy bands, it was spatially coincident with the extended Fermi-LAT source 2FHL J0431.2+5553e (Ackermann et al. 2016). Based on seven years of Pass 8 data and theROSAT data, Cohen (2016) performed a detailed analysis and obtained a hard spectrum with a photon index ofΓ≈1.8 in the 1 GeV to 1 TeV band and an X-ray flux upper limit of ∼5×10⁻¹¹ erg cm⁻² s⁻¹ in the 0.5–2 keV band. For Δα=0.8 and β=2.0, we obtain the following parameters: α=1.7 W_e=0.3×10⁴⁷ erg, B_SNR=4 μG, E_bre=0.4 TeV, and E_cut=60 TeV. Due to the lack of radio flux measurement of the western and northern shells, we double the radio flux of the eastern shell and treat it as the flux for the whole remnant as shown by the gray dots in Figure 1.

3.7.HESSJ1534−571 (G323.7−1.0)

HESSJ1534−571 was discovered in the H.E.S.S. Galactic Plane Survey (Puehlhofer et al. 2015) and has a shell-type structure that matches the radio image of the SNR candidate G323.7−1.0, confirming its classification as an SNR (H.E.S.S. Collaboration et al. 2018d). Its radio counterpart G323.7−1.0 was first discovered in the second epoch Molonglo Galactic Plane Survey at 843 MHz, and shows a faint oval shell with a size of 51'×38' (Green et al. 2014). Based on the Σ−D relation and the CO observations, Maxted et al. (2018b) suggested that this SNR evolved in a wind-blown cavity associated with the Scutum-Crux arm gas at a distance of 3.5 kpc. At this distance, an age of 8–24 kyr and a velocity of 400–1200 km s⁻¹ are derived from the Sedov relation (Maxted et al. 2018b). In this paper, we adopt an age of 10 kyr and a velocity of 1000 km s⁻¹, which are also consistent with that estimated by scaling an RX J1713.7−3946-like SNR in the Sedov phase. Due to the lack of full measurement of its radio flux, we use 0.49 Jy to be a lower limit and 0.98 Jy to be an upper limit following Maxted et al. (2018b). At the X-ray band, based on the analysis of the Suzaku data, no significant nonthermal emission was detected toward this SNR, giving an X-ray flux upper limit of ∼2×10⁻¹¹ erg cm⁻² s⁻¹ in the 2–12 keV band (H.E.S.S. Collaboration et al. 2018d). At the GeV band, Araya (2017) analyzed 8.5 years of Fermi-LAT data and found an extended GeV source with a hard spectrum Γ≈1.34 inside the shell. For Δα=0.8, β=2.0, and infrared photons with a temperature of 20 K and an energy density of 0.8 eVcm⁻³ (Araya 2017), we obtain the following best-fit model parameters: α=1.3, W_e=6.0×10⁴⁷ erg, B_SNR=5 μG, E_bre=10 GeV, and E_cut=15 TeV. For this remnant, the electron index below the break is significantly harder than that of others, which is consistent with the fact that older SNRs trend to have harder radio spectra (Reynolds et al. 2012; Zeng et al. 2019), and may be tested with future radio observations. Given the lack of constraints from radio and X-ray observations, the SED can certainly be fitted with a simpler single power-law model. For the purpose of developing a unified model below, we adopt these broken power-law results.

In our modeling, the two parameters Δαand β are fixed artificially for all sources with a broken power-law particle distribution. Δα=0.8 can give a good fit to spectra of all sources. Of course, the best-fit value may be different for different sources. For example, Δα=0.6 for RXJ0852.0−4622 and Δα=1.0 for RXJ1713.7−3946 will give better fits. In this study, however, we try to develop a unified model. For better comparison of model parameters, it is necessary to have as few free parameters as possible. For SN1006, a single power-law model is adopted and a model with β=2 will lead to a poor fit to the X-ray spectrum. We therefore assume β=1. Moreover, as we will see below, for all sources studied here, the break energy is comparable to the age-limited maximum energy, implying that particle acceleration above the break energy already stops at present and the shape of the cutoff is dictated by the radiative energy loss process. β=2 gives a cutoff sharp enough to fit the X-ray spectra of RXJ1713.7−3946. This is not the case for SN1006, where the acceleration and energy loss timescales are comparable at the cutoff energy.

Because the radiative energy loss timescale at the break energy is much longer than the age of the corresponding SNR, according to the model proposed by Ohira et al. (2012), this break energy should be associated with the age-limited maximum energy. Then, one can obtain η_g for each source. The left panel of Figure 2 shows the dependence of the shock speed u_s, magnetic field B_SNR, and η_g on the age of the SNR. Following Ohira et al. (2012), we divide the evolution of the shock into free expansion and Sedov phases, which starts at t_sed, and assume u_s=u₀(t/t_sed)^m−1, , and η_g=η₀(t/t_sed)^λ for t>t_sed. For t_sed=200 yr, u₀=10,000 km s⁻¹, B₀=120 μG, η₀=1, m=0.4, α_B=0.9, and λ=2.0, we obtain a good fit to the evolution of u_s, B_SNR, and η_g. The latter can be obtained with Equation (2) assuming . The corresponding age-limited maximum energy, cooling-limited maximum energy, and cooling break energy are shown with the red, green, and blue solid lines in the right panel of Figure 2, respectively. The cutoff and break energies of the spectral fits are consistent with such a scenario. We therefore have a unified model for the evolution of the high-energy electron distribution accelerated by their shocks. The decrease of α and E_bre with the SNR age t is also consistent with the finding of Zeng et al. (2019). Since SNRs studied here likely have different kinds of progenitors, it is not surprising that parameters of individual sources may have significant deviations from the unified model. G150.3+4.5 and G323.7−1.0 appear to be evolving in very low density environments with very low magnetic fields, implying very inefficient particle acceleration. Moreover, escape processes may play an important role in reducing the energy content of high-energy electrons trapped in relatively order SNRs.

In the Sedov phase, α_B=0.9 implies that the magnetic field scales as which is predicted by some nonlinear diffusive shock acceleration theory (e.g., Bell 2004) and is supported by observation data including SN 1993J (Vink 2008). If the magnetic field amplification is dominated by turbulence stretching into the downstream, the magnetic field should be proportional to u_s. We also plot this case (α_B = 0.6) with dashed lines in Figure 2 with B₀=70 μG.

As the youngest SNR in our sample, SN1006 is different from others for their single power-law distribution. According to the unified model, this source is going through a particular stage of the shock evolution when the acceleration and energy loss timescale at the cutoff energy is comparable to its age, justifying the single power-law distribution. After this stage, the maximum energy to which the shock can accelerate particles is age limited, corresponding to the break energy of our spectral fit. Electrons above this break energy then only experience radiative energy loss. Therefore, the maximum energy that the shock of an SNR can accelerate electrons to is comparable to the cutoff energy of SN1006. However, the cutoff energies of other SNRs in our sample are higher than those of SN1006. This is mainly due to the fact that the β for SN1006 is different from that for others. For a given high-energy emission, a more gradual cutoff of the electron distribution leads to a lower value for the cutoff energy. We also notice that for a few sources, the energy loss timescale at the cutoff energy is shorter than the age of the corresponding SNR. These high-energy electrons likely come from upstream of the shock where the magnetic field is weak and enter the main emission region in the downstream as the remnant expands. The dotted line in the right panel of Figure 2 shows the electron energy where the synchrotron energy loss timescale in a magnetic field of 8 μG is equal to the age t.

From α_B and λ, one can obtain the time evolution of the break energy . The maximum energy that the shock can accelerate particles to therefore decreases very rapidly as the shock slows down. Such a sharp decrease is closely related to the decrease of the magnetic field and the increase of η_g. Together with the injection rate into the shock acceleration process and the escape rate from SNRs, it determines Δα. For Δα=0.8, we find that the normalization of electrons below the break energy should scale as , implying that the injection rate into the shock acceleration process should be proportional to t^1.5 in the absence of escape. Our spectral fit shows that W_e is more or less independent of t (Table 2). The particle escape process therefore must play a role in determining W_e. Since W_e is almost independent of t, a broken power-law distribution of high-energy electrons may be realized with an energy independent escape time and an injection rate inversely proportional to t. Then the ion injection needs to be different from electrons to give rise to very efficient low-energy ion acceleration in old SNRs (Zhang & Liu 2018).