Cosmic ray protons and electrons from supernova remnants

The spectrum of cosmic ray protons and electrons released by supernova remnants throughout their evolution is poorly known, because of the difficulty in accounting for particle escape and confinement in the downstream of a shock front, where both adiabatic and radiative losses are present. Here we calculate the spectrum of cosmic ray protons released during the evolution of supernovae of different types, accounting for the escape from upstream and for adiabatic losses of particles advected downstream of the shock and liberated at later times. The same calculation is carried out for electrons. The magnetic field in the post-shock region is calculated by using an analytic treatment of the magnetic field amplification due to non--resonant and resonant streaming instability and their saturation. We find that when the field is the result of the growth of the cosmic-ray--driven non--resonant instability alone, the spectrum of electrons and protons released by a supernova remnant are indeed different, but such a difference becomes appreciable only at energies $\gtrsim 100-1000$ GeV, while observations of the electron spectrum require such a difference to be present at energies as low as $\sim 10$ GeV. An effect at such low energies requires substantial magnetic field amplification in the late stages of the supernova remnant evolution (shock velocity $\ll 1000$ km/s), perhaps not due to streaming instability but hydrodynamical processes. We comment on the feasibility of such conditions and speculate on the possibility that the difference in spectral shape between electrons and protons may reflect either some unknown acceleration effect, or additional energy losses in cocoons around the sources.


Introduction
Although there is strong evidence that particle acceleration takes place in supernova remnants (SNRs), it is still debated whether these objects can be the sources of all Galactic cosmic rays (CRs).The problem has various levels: First, the processes of particle acceleration and magnetic field amplification in an individual SNR depend on the type of explosion and the type of environment where it occurs.Second, the spectrum of particles accelerated at the shock and the one that the SNR releases in the surrounding interstellar medium (ISM) are, in general, quite different.The latter is typically made of two components, the one that escapes the remnant at any given time from the upstream region and the one that is advected downstream and is eventually liberated when the SNR dissipates.The particles trapped in the downstream region are affected by losses, which in general act differently on protons and electrons.Even this simple line of thought leads to two conclusions: (1) Although the instantaneous spectrum of particles accelerated at the SNR shock is a power law in momentum, the released spectrum does not need to be; (2) the spectra of protons and electrons from an SNR can, in general, be different.It is worth keeping in mind that nonlinear effects might lead to a slight deviation from the perfect power laws predicted in linear theory (Reynolds & Ellison 1992;Malkov & Drury 2001).
As previously discussed by Cristofari et al. (2020), the spectrum of CR protons from different types of SNRs is hardly a pure power law, and it extends to a maximum momentum that depends rather critically on the type of SNR and on the environment in which the explosion takes place.The structure in the spectrum is due to the contribution of different periods during the evolution of the remnant as well as the overlap of the

Magnetic field amplification in SNRs
In this section we summarize our current understanding of CR-induced magnetic field amplification at SNR shocks.Magnetic field perturbations can be produced in the shock proximity due to a variety of processes, but only some of them lead to important effects in terms of particle scattering.For instance, the propagation of a shock front in a medium with density inhomogeneities can excite a Richtmeier-Meshkov instability (Giacalone & Jokipii 2007) that leads to the growth of perturbations downstream of the shock on a timescale of order /v A , where is the spatial scale of the density inhomogeneities upstream and v A the Alfvén speed.Although such a magnetic field may be important in terms of determining the morphology of synchrotron emission from a remnant and the strength of the downstream magnetic field, it does not appreciably affect the diffusion time of accelerated particles upstream, and hence it does not lead to a substantial increase in the maximum energy that can be reached through diffusive shock acceleration (DSA).A fundamental step forward in the investigation of the interaction between CRs and the surrounding medium has been made with the discovery of the nonresonant hybrid instability (Bell 2004), which is expected to be excited upstream of a shock due to the accelerated particles themselves.This is a current-driven instability, excited both by CR particles that are leaving the acceleration region and by CRs diffusively confined in the vicinity of the shock front.If the density of CRs with momentum >p at the shock is n CR (>p), the electric current that these particles carry is J CR (>p) = eD(∂n CR /∂z) shock = ev sh n CR (>p) and it extends over a precursor distance ∼D(p)/v sh .On the other hand, the density of CRs escaping toward the upstream infinity is limited to the highest energy particles and can be estimated as ∼n CR (>p)(v sh /c).Hence, the corresponding current is J CR (>p) ≈ en CR (v sh /c)c = en CR v sh , numerically equivalent to that estimated above for the same momentum despite the fact that escaping particles are assumed to be streaming away ballistically (namely, moving at roughly the speed of light, c).If the differential spectrum of accelerated particles at the shock is , the normalization can be easily found by requiring that the CR pressure be a fraction, ξ CR , of the ram pressure at the shock location, ρv 2 sh : 1 3 p min /m p c dx x 4−α /(1 + x 2 ) 1/2 .We notice that the normalization constant defined in this way depends very weakly on the minimum and maximum momenta, p min and p max , provided 4 ≤ α < 5, as expected for particle acceleration by DSA.In particular, for α = 4 one has I(α) ≈ ln(p max /m p c).The spectrum of accelerated particles may be outside this range only if nonlinear effects due to the CR pressure lead to the formation of a precursor upstream, which in turn may lead to spectra harder than p −4 (e.g., Jones et al. 2001;Malkov & Drury 2001).However, in practice, strong spectral modification should not be expected, because of numerous other effects (Berezhko & Ellison 1999;Caprioli et al. 2009a).Moreover, recent self-consistent kinetic simulations of strong shocks suggest that the formation of a shock post-cursor naturally leads to spectra with α between 4 and 5 (Caprioli et al. 2020), consistent with γ-ray observations of SNRs (Caprioli 2011).In general, the slope depends on the shock Mach number (if the shock is not strong), but even at the end of the Sedov-Taylor (ST) phase the Mach number remains much larger than unity (see Sect. 4).
The current carried by accelerated particles with momentum >p is Here we especially focus on the scenario where nonresonant hybrid instability is excited by escaping particles since this channel leads to the formation of magnetic perturbations far upstream.As discussed by Bell (2004), the fastest growing mode is associated with a wavenumber, k max , and can be written as γ max = k max v A , where v A = B 0 / 4πρ is the Alfvén speed in the unperturbed field, B 0 .The wavenumber k max is determined based on the condition: and the excitation of the instability occurs only if k max > r L (p) = pc/eB 0 , which translates to the following constraint: In other words, the instability is excited if the energy density in the form of escaping particles is larger than that of the preexisting magnetic field.As discussed in many previous works, there are different approaches to the saturation of the instability.The most intuitive one, based on comparing the plasma displacement due to the J CR × δB/c force with the Larmor radius of particles in the amplified field, δB, leads to the conclusion that the field stops growing when the energy density in the form of escaping particles equals that in the amplified field: At saturation, the spatial size of the perturbations becomes comparable with the Larmor radius in the amplified field, such that the current is now disrupted because of efficient CR scattering off the self-generated perturbations, thereby causing the drive for magnetic field amplification to stop.Some comments on Eq. ( 5) are in order: (1) From the point of view of the scattering of particles with momentum p, the diffusion coefficient is D(p) ∝ p/δB.If the spectrum of accelerated particles is ∼p −4 , then the diffusion coefficient is Bohm-like (linear in momentum).In a general case, α = 4 + , the diffusion coefficient turns out to be D(p) ∼ p 1+ /2 .Since typically ∼ 0 ÷ 0.3, the expected deviations from Bohm-like behavior are small.(2) With the exception of the case = 0, Eq. ( 5) shows a weak dependence upon the momentum p where the current is calculated.(3) In terms of the magnetic field immediately upstream of the shock, all particles should be included in estimating the magnetic field in Eq. (5).For α = 4 ( = 0), the resulting field does not depend on this choice, but for > 0 the  Vink (2012).The lines refer to the growth of nonresonant modes with an injection spectrum with slope α = 4 (solid red) and α = 4.3 (dotted red) and to the recipe adopted by Diesing & Caprioli (2019;dashed green), which is independent of the spectrum of accelerated particles.
resulting magnetic field shows a weak dependence on the minimum momentum.Assuming a minimum momentum χm p v sh , with χ > 1, one has: If the turbulent field upstream of the shock is roughly isotropic and the perpendicular components are compressed at the shock, with compression factor r, the mean value of the compressed downstream magnetic field is B 2 ≈ B 1 (1 + 2r 2 )/3.For r = 4, the compression factor is √ 11.It follows that: For the parameter values appropriate for the Tycho SNR, this expression returns a downstream magnetic field that ranges between 150 µG (α = 4) and 400 µG (α = 4.3), in good agreement with the value inferred from multiwavelength observations (e.g., Morlino & Caprioli 2012).
In Fig. 1 we compare the quantity in Eq. ( 7) with the corresponding value measured in several SNRs as reported by Vink (2012).The solid and dotted lines show the result of Eq. ( 7) for the cases α = 4.3 and α = 4, respectively, assuming χ = 3.When the amplified field upstream becomes smaller than B 0 , the magnetic field relevant for synchrotron losses becomes of order √ 11B 0 .For typical values of parameters in the ISM, this is reflected as B 2 2 /8πρ ∼ 2 × 10 12 erg g −1 in Fig. 1.This condition also identifies the end of the stage where nonresonant modes get excited (see also Amato & Blasi 2009) and typically occurs when v sh 1000 km s −1 .
Equation ( 5) clearly shows how, for this type of instability, when saturation is reached the quantity δB 2 /8πρ scales with the shock velocity as ∝v 3 sh (∝v 2.7 sh for α = 4.3) and is independent of the strength B 0 of the preexisting magnetic field.Other recipes for saturation (see for instance the one proposed by Riquelme & Spitkovsky 2009) suggest some weak dependence of the magnetic field at saturation on B 0 .These recipes lead to a somewhat less amplified magnetic field, a point to keep in mind in the discussion below.
A62, page 3 of 11 The dashed green curve in Fig. 1 represents the same quantity calculated following the recipe originally put forward by Caprioli (2011) and used in Morlino & Caprioli (2012).A few comments are in order concerning the genesis of this approach.Resonant streaming instability was originally included in nonlinear theories of DSA by Amato & Blasi (2006), where the equation for the growth of these modes was solved analytically.The saturation level derived from this simple approach can be written as: where M A = v sh /v A,0 is the Alfvén Mach number calculated with respect to the Alfvén speed in the field B 0 , v A,0 = B 0 / 4πρ.It follows that the compressed field downstream reads: It is worth noting that, since M A = v sh /v A,0 , effectively B 2 2 ∝ v sh .In Fig. 1, this is shown as a dashed purple line.We should stress that this amplified field should be used only when the Bell modes are not allowed to grow (low shock velocity).In the opposite condition, the growth rate of the resonant instability is proportional to n 1/2 CR (>p) and leads to a lower magnetic field than the nonresonant instability, and usually δB/B 0 1 (Zweibel 1979;Amato & Blasi 2009).Morlino & Caprioli (2012) proposed a formal modification of Eq. ( 9) so as to mimic the onset of nonlinear effects: The proposed modification consists in attempting to interpret the Mach number in Eq. ( 8) as the one calculated with respect to the Alfvén speed in the amplified field, δB.This assumption results in a modified expression for the amplified field that, after compression at the shock, reads: This quantity is plotted in Fig. 1 for r = 4 (dashed green line).
The most distinctive characteristic of this trend is the scaling with v 2 sh , which is quite different from the ∝v 3 sh typical of nonresonant hybrid modes (see Eq. ( 7) for α = 4).Equation ( 7) results in a much larger magnetic energy density downstream at late times, when the shock velocity drops.
Though being a viable phenomenological approach, this is based on a few assumptions that are worth keeping in mind: (1) Although based on a perturbative approach, Eq. ( 8) leads to B 1 /B 0 1 for high shock speeds, which is unphysical if modes are strictly resonant.(2) As pointed out above, the growth rate used to describe resonant modes, and which led to Eq. ( 8), is only valid in the case in which the energy density of particles is much smaller than B 2 0 /4π, which is certainly not the case for fast shocks (in fact, this is the very regime where nonresonant hybrid modes grow).This was discussed at length by Amato & Blasi (2009, see also Blasi 2013).
To add to this rather confused situation, hybrid kinetic simulations of particle acceleration at shocks lead to a prescription for the amplified field that reads (Caprioli & Spitkovsky 2014): which is formally similar to Eq. ( 8) but with a different numerical factor.It is not easy to encapsulate these simulation results in the above theory: The simulations were run for a range of parameters in which Bell modes would grow and are seen to be growing.On the other hand, these simulations are nonrelativistic, which implies that the anisotropy of the accelerated particles (which for an astrophysical shock is ∼v sh /c ∼ 10 −2 ÷ 10 −3 ) is ∼v sh /v ∼ 0.1.This might lead to a larger fraction of the CR energy being channeled into CRs compared with an SNR shock.Moreover, if accelerated particles are forced to remain nonrelativistic, the growth of the Bell instability can be shown to saturate to a value that is larger than the one in Eq. ( 6).In any case, no dependence upon the background field B 0 is expected, based on semi-analytical arguments, unless more complex phenomena come into play.
Finally, the magnetic field energy density in Eq. ( 11) shows a scaling B 2 2 /ρ ∝ v sh , but the nonrelativistic simulations discussed by Caprioli & Spitkovsky (2014) were run for a given shock velocity and changing B 0 so as to achieve different Alfvén Mach numbers.This procedure is optimal for unraveling the scaling of B 2 with v sh /v A but would not reveal the additional scaling with the actual CR speed, c/v A .These aspects definitely deserve further investigation.
In the following we focus on the investigation of the two prescriptions for magnetic field amplification described in Eqs. ( 7) and ( 10).We also comment on the case in which the field is amplified through nonresonant streaming instability when the appropriate condition is satisfied and by resonant streaming instability (Eq.( 8)) at later times.
A quick inspection of Fig. 1 shows that data points are too sparse to allow us to infer a clear dependence of the amplified field on shock speed, v sh , although there is a mild preference for a ∝v 3 sh trend.As one can see from Eq. ( 7), the nonresonant hybrid modes are expected to lead to B 2 2 /8πρ ∼ v 7−α sh , which compares well with the observed trend for α = 4÷4.3.As discussed below, most of the modification of the electron spectrum due to energy losses occurs at late times (low shock speeds).Hence, for the purpose of calculating the difference in spectrum between electrons and protons, it is of critical importance to understand what the downstream magnetic field for older SNRs is.
The importance of the nonresonant hybrid instability for particle scattering is well known and will be briefly summarized here only for the sake of completeness.Following Bell et al. (2013) and Cristofari et al. (2020), the maximum energy of protons can be estimated by requiring the condition t 0 dt γ max (t ) 5, which leads to: Introducing the expression for the total magnetic field upstream, B 1 , one can rewrite this expression as where Ωc = eB 1 /m p c is the cyclotron frequency and ṽA = B 1 / 4πρ is the Alfvén speed, both calculated in the amplified magnetic field, B 1 .It is worth pointing out that for α = 4.3 the total magnetic field at the shock is slightly larger than for α = 4, but the corresponding p max is lower, because there is less power available at the scales resonant with particles with momentum p max .The maximum momentum of electrons is computed by equating the acceleration time to the minimum between the synchrotron time and the age of the SNR (see, e.g., Blasi 2010).

Cumulative spectra of CRs at SNRs
The contribution of CRs from an individual SNR can be written as the sum of two components: (1) particles accelerated at the SNR shock and trapped downstream of the expanding shell until the time when the shock dissipates away and the particles are released into the ISM after having suffered the effect of losses, N loss , and (2) particles accelerated at the shock up to the highest energy that can be achieved at that time and leaving the acceleration region moving upstream, N esc .Although the latter contribution is strongly peaked around the maximum energy reached at that time, the integration over the whole expansion history of the remnant leads to a continuous spectrum of CRs released into the ISM (Caprioli et al. 2009b).
The spectrum of particles accelerated at a strong shock in the test-particle limit (Caprioli et al. 2009b), assuming a free escape boundary condition at some location upstream, reads: For p p max (t), f p (p) ∝ (p/p inj ) −3r/(r−1) , and for p p max (t), f p (p) ∝ exp − p p max (t) .Working under the usual assumption that a fraction of the ram pressure of the SNR shock is converted into CRs, it is easy to write the spectrum of CRs accelerated at the SNR shock at any given time.In a given time interval, dt, the number of particles of momentum p dn acc (p, t) accelerated at the shock reads: where v sh /r is the velocity downstream of the SNR shock and r is the compression factor at the shock.
If accelerated particles could escape the acceleration region immediately after penetrating downstream, the total number of particles integrated over the entire active lifetime of an SNR, from t 0 , typically the beginning of the free expansion phase, to T SN , the end of the ST phase, would read: This quantity could be interpreted as the flux of CR particles contributed by an SNR in the absence of adiabatic energy losses.
The departures from such a spectrum will be used later as an index of the effect of losses on the spectrum of CRs from each SNR.In principle, if the SNR shell is broken in some locations where escape is allowed, the time-integrated CR flux would be somewhat similar to that estimated using Eq. ( 16).
The spectrum of electrons at the SNR shock is calculated as in Morlino et al. (2009), using the approximated expression proposed by Zirakashvili & Aharonian (2007): where K ep is the electron-to-proton ratio, typically in the range 10 −4 −10 −2 .The cumulative spectrum of electrons, N e acc , can then be calculated as in Eq. ( 16), using f e instead of f p .
Below, we calculate the spectrum of protons and electrons that leaves an SNR after accounting for the effect of adiabatic and radiative energy losses.

Adiabatic and radiative losses of CRs downstream of the shock
As the SNR shock expands, the particles produced at the shock, and trapped inside the SNR downstream of the shock, suffer adiabatic and radiative losses.The latter are dominated by the emission of synchrotron photons, while inverse-Compton scattering, although included in our calculation, is typically negligible.In fact, the rate of inverse-Compton losses on cosmic microwave background photons is the same as that of synchrotron losses in a ∼3 µG magnetic field, and post-shock fields are much larger due to amplification and compression.On the other hand, inverse-Compton losses might become important for SNRs located in star-forming regions if the energy density of optical-to-infrared photons exceeds ∼100 eV cm −3 .The number of particles (protons or electrons) that are liberated by an individual SNR at the end of the evolution can be easily written in terms of the conservation of the total number of particles.In fact, the number of particles with momentum p at the end of the SN evolution (t = T SN ) is the result of all the particles produced at earlier times (t < T SN ) with momentum p > p, such that in a time T SN the momentum has degraded down to p.We can then write: The change in momentum of a particle injected at time t with momentum p due to losses is: where σ T is the Thompson cross section and L accounts for adiabatic energy losses in terms of the change in volume between the two times, t and t: where ρ down is the density downstream of the shock.If the expansion is adiabatic, ρ down ∝ P 1/γ ∝ (ρv 2 sh (t)) 1/γ , where ρ is the gas density upstream of the shock.For protons, synchrotron losses are negligible, while for electrons both adiabatic and radiative losses are important.

Escaping particles
The flux of particles escaping the accelerator from the upstream region can be written following Caprioli et al. (2009b, see also Cristofari et al. 2020): where we introduced the function which describes the spectral shape of the particles that can escape upstream at a given time.For protons, the function G(p, t) is strongly peaked around the maximum momentum of protons, p max (t), at the given time, t.For electrons, when the maximum A62, page 5 of 11 energy is determined by diffusion, the meaning of this function is the same as for protons.When the magnetic field is large enough that synchrotron losses dominate the maximum energy of accelerated electrons (namely for young SNRs), the function G is vanishingly small and no escape upstream is possible.

Evolution of the shock in the circumstellar environment
The evolution in time of SNR shocks directly impacts the spectra of particles injected by these sources into the ISM, when integrated on the SNR life span.In this section we briefly summarize the calculations adopted here for the description of such evolution, which is discussed in more detail in Cristofari et al. (2013).
For the sake of calculating the contribution to the CR flux, SNe are broadly classified into two groups, depending on the mechanism triggering the explosion: thermonuclear SNe (type Ia) and core-collapse SNe (type II).In addition, we consider a peculiar type of very energetic core-collapse SNe, which we call type II*.This type is introduced here just to illustrate the wide range of physical parameters that apply to type II SNe and to demonstrate that these energetic events can accelerate particles up to the petaelectronvolt range (when E SN and Ṁ are sufficiently high and M ej is sufficiently low; see the discussion on the parameter space presented by Cristofari et al. 2020).Schematically, SNR shocks from type Ia SNe expand in a uniform ISM.The evolution of the shock radius, R sh , and shock velocity, v sh , in the ejecta-dominated phase and the ST phase are well described by self-similar solutions (see, e.g., Chevalier 1982;Truelove & McKee 1999;Ptuskin & Zirakashvili 2005).
Here we rely on the approach presented in Truelove & McKee (1999, 2000), adopting the same formalism as in Cardillo et al. (2015): where , λ FE = 4/7, λ ST = 2/5, and . This description holds until the end of the ST phase (i.e., the beginning of the radiative phase) of SNR evolution.This transition typically occurs when the age of the SNR becomes on the order of the cooling time: t cool ≈ 10 3 n 0 1 cm −3 −1 v sh (t) 10 8 cm s −1 3 kyr (Blondin et al. 1998).In this work, we assume E SN = 10 51 erg, M ej = 1 M , and n 0 = 1 cm −3 as reference values for type Ia SNe.It is worth recalling that in the presence of efficient CR acceleration at the SNR shock, as discussed by Diesing & Caprioli (2018), the evolution in time of the shock in the final stages of the evolution may be affected rather remarkably by the CR pressure, and the beginning of the radiative phase may be delayed.We do not include these effects here.
In general, type II SNR shocks expand in a complex medium structured by the evolution of the massive progenitor star before the explosion of the SN.Throughout its main sequence, the stellar wind produces a low-density hot-temperature bubble that is in pressure equilibrium with the ISM.Later, when the massive star reaches the final stages of stellar evolution, typically entering the red super giant (RSG) stage, a slow dense wind forms.
Therefore, after the SN explosion, the shock expands through a dense wind and then through a low density bubble before finally reaching the ISM.The density of the dense wind created by the RSG is typically n w = Ṁ/(4πmu w r 2 ), were Ṁ ∼ 10 −5 M yr −1 is the mass-loss rate, u w ∼ 10 6 cm s −1 the velocity of the wind, and m = m p (1 + 4 f He )/(1 + f He ) ∼ 1.27m p the mean mass of the ISM per hydrogen nucleus (Weaver et al. 1977).The density of the low density bubble is n b = 0.01 (L 6  36 n 19 0 t −22 Myr ) 1/35 cm −3 , with t Myr the duration of the main sequence, which is on the order of a megayear (Longair 1994).The transition between the dense RSG wind and the low density cavity, r 1 , is set by equating the RSG wind pressure to the thermal pressure of the hot cavity: where k is the Boltzmann constant.The radius of the hot bubble is r b = 27(L 36 /1.27n 0 ) 1/5 t 3/5 Mpc pc, where L 36 is the main sequence star wind power in units of 10 36 erg s −1 .In the case of type II SNe, the evolution in the structured medium described above does not allow for direct self-similar solutions.It is, however, possible to work under the "thin-shell" approximation, considering that the swept-up gas is located in a think layer behind the shock wave (see, e.g., Ostriker & McKee 1988;Bisnovatyi-Kogan & Silich 1995).In the case of spherically symmetric distribution, Ptuskin & Zirakashvili (2005) obtained: and For a typical type II SN, we assume E SN = 10 51 erg, M ej = 5 M , and Ṁ = 10 −5 M yr −1 .For a type II*, we assume E SN = 5 × 10 51 erg, M ej = 1 M , and Ṁ = 10 −4 M yr −1 .

Results
In this section we discuss our results in terms of the spectrum of protons and electrons injected into the ISM by different types of SN explosions.As discussed above, the CR spectrum released by an individual source is generally the sum of (1) the contribution due to particles that escaped upstream at any given time and (2) the contribution of CRs trapped downstream, where energy losses are in action, and liberated in the final stages of the SN evolution, which we assume coincide with the beginning of the radiative phase.While the proton escape flux is always present, the electron escape flux vanishes when radiative losses limit the maximum energy.It is this phase of radiative and adiabatic losses that can potentially make the overall spectrum of electrons different from that of protons if radiative losses are sufficiently severe.
The spectrum of protons liberated by type Ia, II, and II* SNRs is shown in Fig. 2, with the corresponding slopes shown in the bottom panels.As discussed by Cristofari et al. (2020), the highest energies are reached at very early times, but they do not reflect the equally large energies in the final spectrum of protons because the amount of mass that the SN shock processes at such early times is very small.In fact, the effective maximum energy is the one reached at the beginning of the ST phase.The easiest case is that of type Ia SNRs, where the explosion occurs in  (thin lines) if they were instantaneously liberated into the ISM (broken shell assumption).The dashed curves illustrate the effect of adiabatic losses in the downstream region, while the dotted lines refer to the escape flux from the upstream region.In the bottom part of each panel we also show the local slope of the spectrum q(p) at a given momentum.
the normal ISM, with a spatially constant gas density and background magnetic field.For type Ia SNRs the effective maximum energy is a few tens of TeV (left panel of Fig. 2).There is an additional spectral steepening at somewhat lower energies due to the temporal evolution of the maximum energy.More specifically, the steepening occurs at the maximum energy reached at the end of the ST phase, typically a few TeV.The flux of escaping CR protons starts at about the same energy, as is clearly visible in Fig. 2. For a strong shock, such as the one expected for a young SNR expanding in the normal ISM, the spectrum of accelerated particles at the shock location has a slope very close to 4 (thick lines in Fig. 2).Nevertheless, as recently discussed by Caprioli et al. (2020), the spectrum can be steeper if the finite velocity of scattering centers in the downstream plasma is taken into account.For this reason, in Fig. 2 we also show the case α = 4.3 (thin lines).In all cases of interest, the spectra of CR protons that are injected into the ISM (as the sum of the two contributions) are quite close to the spectrum at the shock in terms of slope, with the exception of the highest energies, as discussed above.
For type II SNRs, the spectrum of CR protons is shown in the middle panel of Fig. 2. For the sake of making a fair comparison between the three types of SN explosions, here we used an acceleration efficiency of ξ CR = 0.1 for all of them.As discussed by Cristofari et al. (2020), because of the different rates of occurrence of these events in the Galaxy, for type II SNRs the efficiency is required to be somewhat lower than for type Ia, which is also reflected in a lower value of the maximum energy of particles accelerated at the shock (see Eq. ( 12)).Despite this bias, the maximum achievable energy for type II SNRs remains on the order of ∼10 5 GeV and falls short of the knee by a large amount, as already pointed out by Cristofari et al. (2020).
Only when parameters are pushed to the extreme (what we have called here type II* SNRs) can the maximum energy reach the knee, as shown in the right plot of Fig. 2. As already pointed out by Caprioli et al. (2009b), the superposition of the escape flux from the different stages of shock evolution in the complex environment around these SNRs may lead to the appearance of bumps in the overall CR spectrum that might be related to the feature recently measured by DAMPE in the 10−100 TeV region of the proton spectrum (An et al. 2019).
The corresponding spectra of electrons injected by SNRs of different types into the ISM are shown in Fig. 3.The thick and thin curves refer to α = 4 and α = 4.3, respectively.The dashdotted line identifies the spectrum of particles accelerated at the shock, as if they were immediately liberated into the ISM, without energy losses.The solid lines are the spectra of electrons liberated into the ISM after adiabatic and synchrotron losses downstream of the shock, while the upstream escape flux, limited to the times when the maximum energy of electrons is not determined by energy losses, is shown in the form of dotted lines.If the SNR shell were broken or if confinement in the downstream region were energy-dependent (e.g., due to turbulence damping), the actual contribution would lie between the dashdotted and solid lines.
The rate of synchrotron losses is larger when the condition for the growth of the magnetic field through the excitation of the nonresonant hybrid instability is fulfilled.As discussed in Sect.2, B 2 2 /ρ ∝ v 7−α sh for this instability, and hence the mechanism becomes less effective or even ineffective in the late stages of SNR evolution; these stages are, however, crucial for the production of low energy electrons.As a consequence, the effect of radiative energy losses is only important at energies at or above teraelectronvolt levels, while it is minor at lower energies, as clearly shown in Fig. 3, independent of the type of SNR considered.It follows that if the magnetic field amplification is mainly due to the excitation of Bell modes, then losses cannot be the main reason for a difference in the spectra of protons and electrons in SNRs.This result is illustrated even more clearly in Fig. 4, where we show the difference in slope between protons and electrons as a function of energy for type Ia (left), type II (middle), and type II* (right) SNRs.The thick and thin curves refer to the cases α = 4 and α = 4.3, respectively.The dashed lines were obtained using the growth of nonresonant hybrid modes to calculate the magnetic field, while the dotted lines are based on Eq. ( 10), the same as in Diesing & Caprioli (2019).The reason for the very different results is that the two recipes lead to quite different predictions for the downstream magnetic field at late times, which are the most important times for determining the electrons' spectrum: The nonresonant instability leads to an intense field in the early phases when the shock is faster, but the instability virtually shuts off at later times when the shock velocity drops below ∼1000 km s −1 .On the other hand, the prescription used by Diesing & Caprioli (2019) leads to larger magnetic fields being retained even in late stages, such that the acceleration of electrons remains loss-dominated at such late times.This is reflected in steeper electron spectra.The main lesson to be learned from this exercise is that the spectrum of electrons liberated by SNRs into the ISM is strongly dependent upon knowledge of the magnetic field in the phases of SNR evolution about which we know the least.
The issue boils down to whether we have observational bounds or theoretical arguments that can be used to assess the credibility of different scenarios.From the observational point of view, to our knowledge, there is no evidence in favor of or against having relatively large magnetic field amplification at late times.From the theoretical point of view, the only prejudice that we can mention is that the recipe based on Eq. ( 10) is somewhat weak, for the reasons explained above.To be somewhat more conservative and limit our attention to scenarios that are based on some sort of theoretical ground, one could say that the magnetic field downstream is provided by Eq. ( 7) when the energy density of accelerated particles satisfies the condition for the growth of the nonresonant instability, and by Eq. ( 9) when the nonresonant instability stops growing and only resonant Alfvén modes get excited.In both cases, when the amplified magnetic field becomes much smaller than B 0 , the field to be used for the purpose of calculating the electrons' energy losses reduces to B 0 .
To further address the importance of the late times in SNR evolution for shaping the electron spectrum, we considered a rather extreme situation in which particle acceleration is efficient until 80 kyr.Indeed, it has been shown that the duration of the ST phase can be substantially extended, for instance because of the CR pressure on the SNR shock (Diesing & Caprioli 2018).In Fig. 5, we illustrate the results for a type Ia SN (for a slope 4.3 of the spectrum of accelerated particles): Although the effect of losses is more pronounced (in the Bell and Bellplus-nonresonant cases), it still remains limited to high energies, while no systematic difference is visible in the range 10−1000 GeV.

Discussion
Supernova remnants contribute a spectrum of CRs comprising two terms: one is the flux of particles that escape from upstream (2) The maximum energies of protons as discussed above are mainly the result of the self-confinement of CRs in the shock proximity due to the excitation of the current-driven nonresonant streaming instability (see also Cristofari et al. 2020).We calculated the strength of the magnetic field downstream of the shock that would be expected based on this process as a function of the shock velocity, and we compared the results with the observed trend, as reported by Vink (2012).The scaling B 2 2 /ρ ∝ v 3 s , typical of this instability, is in good agreement with the sparse data available and provides the correct normalization.For old SNRs, with typical velocities of v s 1000 km s −1 , the nonresonant instability is quenched and the strength of the magnetic field is comparable with that in the surrounding medium, δB B 0 .A recipe for magnetic field evolution often used in the literature was put forward by Caprioli (2011) and Morlino & Caprioli (2012).For typical speeds of SNR shocks, such a recipe leads to larger values of the amplified field, especially for older SNRs.As discussed in Sect.2, this recipe has numerous caveats that are important to keep in mind.
(3) The amplified magnetic field described above plays a crucial role in shaping the spectrum of the accelerated electrons that are advected downstream because of radiative losses.In the cases we investigated, a substantial steepening of the spectrum was only obtained when the phenomenological prescription of Caprioli (2011) was adopted, which leads to results that are in qualitative agreement with those of Diesing & Caprioli (2019).When the magnetic field is assumed to be amplified according to the growth of the nonresonant (Bell 2004) and resonant streaming instability (Amato & Blasi 2006), the effect of losses on the electron spectrum typically reduces to a cutoff in the overall spectrum rather than a broad steepening.
The appearance of a steepening in the electron spectrum (with respect to the proton spectrum) depends on the time dependence of the maximum electron energy, p e max (t).In the initial stages of SNR evolution (ejecta-dominated to the early ST phase), the maximum energy is determined by energy losses: In the few hundred µG typical of an SNR at the beginning of the ST phase, synchrotron losses degrade the particle energy to ∼10 GeV by the end of the remnant life.Later times (in the ST phase) may contribute electrons with increasingly larger energies if the amplified magnetic field decreases sufficiently fast with the shock velocity, as is typically the case for all the above prescriptions.When p e max (t) becomes age-limited, as for protons, contributions from different times pile up, eventually leading to an overall cutoff in the teraelectronvolt range.In summary, sustained magnetic field amplification may provide a steepening in the electron spectrum, but specific recipes for the amplification determine whether the steepening extends over several orders of magnitude in energy or is instead limited to a narrow energy range, thereby appearing as a cutoff in the electron spectrum.Current data suggest that a global steepening must be present between ∼10 GeV and a few TeV (Evoli et al. 2020), although it is not clear that synchrotron losses are the cause of this phenomenon (Ohira et al. 2012).
This result shows that if SNRs are responsible for the observed spectrum of electrons, the late phases of their evolution, usually considered to be of little interest for CR acceleration, are in fact crucial in the sense that synchrotron energy losses may be at work for long periods.On the other hand, the streaming instabilities that are thought to be most effective in amplifying the magnetic field may not be able to provide the magnetic fields that, through synchrotron losses, would cause the required steepening.Since such stages are crucial to shaping the spectrum of electrons, it is worth asking whether there are physical phenomena that we did not take into account that might lead to larger (or smaller) fields at late times: One point to keep in mind is that the total compression factor at CR-modified shocks might be somewhat larger than 4 if CR acceleration remains efficient; for instance, in the simulations of Haggerty & Caprioli (2020), one has r ∼ 6−7.This can indeed lead to somewhat larger magnetic fields, as one can infer from Eq. ( 7).However, in order to have a sizeable effect on the spectrum of electrons, this should happen for low shock velocities, while it is typically expected to be more of a concern for fast shocks.One could also speculate that, in addition to CR-induced magnetic field amplification (which takes place upstream), the field could be further amplified downstream of the shock, perhaps through Richtmeier-Meshkov instability.As discussed above, the duration of the ST phase, which may also be controlled by the nonthermal CR pressure (Diesing & Caprioli 2018), may affect the spectrum of electrons as well.Finally, we note that if shock acceleration were efficient even when the shock reaches a Mach number of a few, the spectrum of accelerated particles, which is affected by the compression factor, would become steeper.For typical values of the parameters, when the shock reaches the radiative phase, its Mach number is still ∼40, so it is possible that the CR spectrum may still evolve during the early radiative stage.
Clearly, there are also physical processes that lead to smaller magnetic fields other than those predicted above, thereby making the spectral modification for electrons even less effective.For instance, when the shock moves in the ordinary ISM, as is typically the case in the late stages of the evolution of an SNR, neutral hydrogen induces a strong level of ion-neutral damping, which substantially limits the growth of perturbations.These effects are not included in our calculations, nor in numerical simulations of DSA.
From the discussion above, it is clear that if the difference in the spectrum of electrons and protons is to be attributed to radiative losses in the downstream region of an SNR shock, the required conditions appear to be rather at odds with those that we would typically assume to exist.In particular, some rather efficient mechanism for magnetic field amplification should come into effect for late, slow moving shocks.A dedicated effort to investigate these stages is needed.
On the other hand, if the results of our investigation are taken at face value, then the spectral shape of electrons and protons liberated into the ISM by an individual SNR should be very similar.Hence, it would follow that the observed difference should be attributed to phenomena occurring after the particles have been released into the ISM.If the diffusion coefficient describing transport in the Galaxy is the same for the two species, as one would expect, the only possibility left open is that electrons and protons may develop different spectral shapes while propagating in the neighborhood of the source due to the large perturbations induced by the escaping particles (Schroer et al. 2020).This possibility is currently being investigated.

Fig. 2 .
Fig.2.Spectra of protons produced at SNRs from type Ia (top), type II (center), and type II* (bottom) SNRs for α = 4 (thick lines) and α = 4.3 (thin lines) if they were instantaneously liberated into the ISM (broken shell assumption).The dashed curves illustrate the effect of adiabatic losses in the downstream region, while the dotted lines refer to the escape flux from the upstream region.In the bottom part of each panel we also show the local slope of the spectrum q(p) at a given momentum.

Fig. 3 .
Fig.3.Spectra electrons produced at SNRs from type Ia (top), type II (center), and type II* (bottom) progenitors for α = 4 (thick lines) and α = 4.3 (thin lines).The dash-dotted line indicates the spectrum of particles accelerated at the shock, the solid lines the spectra of electrons liberated into the ISM after losses downstream of the shock, and the dotted lines the upstream escape flux.In the bottom part of each panel we also show the slope of the spectrum q(p) at a given momentum.