# On “box” models of shock acceleration and electron synchrotron spectra

## 1 Introduction

The EGRET detection aboard the Compton Gamma Ray Observatory of at least two supernova remnants (Esposito esposito96 (1996)) and more than fifty active galactic nuclei (Thomson et al. thompson95 (1995)) has given strong evidence of particle acceleration in these objects. This evidence is strengthened even more by the detection of SN 1006 (Tanimori tanimori98 (1998)) and the BL Lac objects Mkn 421 (Punch et al. punch92 (1992)) and Mkn 501 (Quinn et al. quinn96 (1996)) by ground based Cherenkov detectors at TeV energies.

A particularly attractive mechanism for producing the required radiating high energy particles is the diffusive shock acceleration scheme, which has already been put forward to predict TeV radiation from supernova remnants (Drury, Aharonian & Völk drury94 (1994), Mastichiadis masti96 (1996)) or explain the observed flaring behaviour in X-rays and TeV $`\gamma `$-rays from active galactic nuclei (Kirk, Rieger & Mastichiadis kirk98 (1998)). This scheme was originally proposed as the mechanism responsible for producing the nuclear cosmic ray component in shock waves associated with supernova remnants (Krymsky krymsky77 (1977), Axford axford81 (1981)). Based on this picture many authors (Bogdan & Völk bogdan83 (1983), Moraal & Axford moraal83 (1983), Lagage & Cesarsky lagage83 (1983), Schlickeiser schlickeiser84 (1984), Völk & Biermann volk88 (1988), Ball & Kirk ball92 (1992), Protheroe & Stanev protheroe98 (1998)) have used, under various guises, a simplified but physically intuitive treatment of shock acceleration, sometimes referred to as a “box” model.

In this paper we examine the underlying assumptions of the “box” model (§2) and we present an alternative more physical version of it (§3). We then include synchrotron and inverse Compton losses as a means of spectal modification and we determine the conditions under which “pile-ups” can occur in shock accelerated spectra (§4). The “box” model can also be extended to include the nonlinear effect of the particle pressure on the background flow (§5).

## 2 The “box” model of diffusive shock acceleration

The main features of the “box” model, as presented in the literature (see references above) and exemplified by Protheroe and Stanev (protheroe98 (1998))) can be summarised as follows. The particles being accelerated (and thus “inside the box”) have differential energy spectrum $`N(E)`$ and are gaining energy at rate $`r_{\mathrm{acc}}E`$ but simultaneously escape from the acceleration box at rate $`r_{\mathrm{esc}}`$. Conservation of particles then requires

$$\frac{N}{t}+\frac{}{E}\left(r_{\mathrm{acc}}EN\right)=Qr_{\mathrm{esc}}N$$
(1)   
where $`Q(E)`$ is a source term combining advection of particles into the box and direct injection inside the box.

In essence this approach tries to reduce the entire acceleration physics to a “black box” characterised simply by just two rates, $`r_{\mathrm{esc}}`$ and $`r_{\mathrm{acc}}`$. These rates have of course to be taken from more detailed theories of shock acceleration (eg Drury drury91 (1991)). A minor reformulation of the above equation into characteristic form,

$$\frac{N}{t}+r_{\mathrm{acc}}E\frac{N}{E}=QN\left(r_{\mathrm{esc}}+r_{\mathrm{acc}}+E\frac{r_{\mathrm{acc}}}{E}\right)$$
(2)   
is useful in revealing the character of the description. This is equivalent to the ordinary differential equation,

$$\frac{dN}{dt}=QN\left(r_{\mathrm{esc}}+r_{\mathrm{acc}}+E\frac{r_{\mathrm{acc}}}{E}\right)$$
(3)   
on the family of characteristic curves described by

$$\frac{dE}{dt}=r_{\mathrm{acc}}E$$
(4)   
giving the formal solution,

$`N(E,t)={\displaystyle _0^t}`$ $`Q`$ $`(t^{},E^{})`$ (5)   
$`exp\left[{\displaystyle _t^{}^t}\left(r_{\mathrm{esc}}+r_{\mathrm{acc}}+E{\displaystyle \frac{r_{\mathrm{acc}}}{E}}\right)𝑑t^{\prime \prime }\right]dt^{}.`$   
The number of particles at energy $`E`$ and time $`t`$ in the “box” is given simply by an exponentially weighted integral over the injection rate at earlier times and lower energies. Of particular interest is the steady solution at energies above those where injection is occuring which is easily seen to be a power-law with exponent

$$\frac{\mathrm{ln}N}{\mathrm{ln}E}=\left(1+\frac{r_{\mathrm{esc}}}{r_{\mathrm{acc}}}+\frac{\mathrm{ln}r_{\mathrm{acc}}}{\mathrm{ln}E}\right).$$
(6)   
At first sight (to one familiar with shock acceleration theory) it appears odd that the exponent depends not just on the ratio of $`r_{\mathrm{esc}}`$ to $`r_{\mathrm{acc}}`$ but also on the energy dependence of $`r_{\mathrm{acc}}`$. However, as remarked by Protheroe and Stanev, the physically important quantity is not the spectrum of particles inside the fictitious acceleration “box” but the escaping flux of accelerated particles $`r_{\mathrm{esc}}N`$ and this is a power-law of exponent

$$\frac{\mathrm{ln}(r_{\mathrm{esc}}N)}{\mathrm{ln}E}=\left(1+\frac{r_{\mathrm{esc}}}{r_{\mathrm{acc}}}+\frac{\mathrm{ln}r_{\mathrm{acc}}}{\mathrm{ln}E}\frac{\mathrm{ln}r_{\mathrm{esc}}}{\mathrm{ln}E}\right).$$
(7)   
Thus provided the ratio of $`r_{\mathrm{acc}}`$ to $`r_{\mathrm{esc}}`$ is fixed, the power-law exponent of the spectrum of accelerated particles escaping from the accelerator is determined only by this ratio whatever the energy dependence of the two rates.

## 3 Physical interpretation of the box model

We prefer a very similar, but more physical, picture of shock acceleration which has the advantage of being more closely linked to the conventional theory. For this reason we also choose to work in terms of particle momentum $`p`$ and the distribution function $`f(p)`$ rather than $`E`$ and $`N(E)`$.

The fundamental assumption of diffusive shock acceleration theory is that the charged particles being accelerated are scattered by magnetic structures advected by the bulk plasma flow and that, at least to a first approximation, in a frame moving with these structures the scattering changes the direction of a particle’s motion, but not the magnitude of its velocity, energy or momentum. If we measure $`p`$, the magnitude of the particle’s momentum, in this frame, it is not changed by the scattering and the angular distribution is driven to being very close to isotropic. However if a particle crosses a shock front, where the bulk plasma velocity changes abruptly, then the reference frame used to measure $`p`$ changes and thus $`p`$ itself changes slightly. If we have an almost isotropic distribution $`f(p)`$ at the shock front where the frame velocity changes from $`𝐔_\mathrm{𝟏}`$ to $`𝐔_\mathrm{𝟐}`$, then it is easy to calculate that there is a flux of particles upwards in momentum associated with the shock crossings of

$`\mathrm{\Phi }(p,t)`$ $`=`$ $`{\displaystyle p\frac{𝐯(𝐔_\mathrm{𝟏}𝐔_\mathrm{𝟐})}{𝐯^\mathrm{𝟐}}p^2f(p,t)𝐯𝐧𝑑\mathrm{\Omega }}`$ (8)   
$`=`$ $`{\displaystyle \frac{4\pi p^3}{3}}f(p,t)𝐧(𝐔_\mathrm{𝟏}𝐔_\mathrm{𝟐})`$   
where $`𝐧`$ is the unit shock normal and the integration is over all directions of the velocity vector $`𝐯`$. Notice that this flux is localised in space at the shock front and is strictly positive for a compressive shock structure.

This spatially localised flux in momentum space is the essential mechanism of shock acceleration and in our description replaces the acceleration rate $`r_{\mathrm{acc}}`$. The other key element of course is the loss of particles from the shock by advection downstream. We note that the particles interacting with the shock are those located within about one diffusion length of the shock. Particles penetrate upstream a distance of order $`L_1=𝐧𝐊_\mathrm{𝟏}𝐧/𝐧𝐔_\mathrm{𝟏}`$ where $`𝐊`$ is the diffusion tensor and the probability of a downstream particle returning to the shock decreases exponentially with a scale length of $`L_2=𝐧𝐊_\mathrm{𝟐}𝐧/𝐧𝐔_\mathrm{𝟐}`$. Thus in our picture we have an energy dependent acceleration region extending a distance $`L_1`$ upstream and $`L_2`$ downstream. The total size of the box is then $`L(p)L_1(p)+L_2(p)`$. Particles are swept out of this region by the downstream flow at a bulk velocity $`𝐧𝐔_\mathrm{𝟐}`$.

Conservation of particles then leads to the following approximate description of the acceleration,

$$\frac{}{t}\left[4\pi p^2fL\right]+\frac{\mathrm{\Phi }}{p}=Q𝐧𝐔_\mathrm{𝟐}4\pi p^2f,$$
(9)   
that is the time rate of change of the number of particles involved in the acceleration at momentum $`p`$ plus the divergence in the accelerated momentum flux equals the source minus the flux carried out of the back of the region by the downstream flow. The main approximation here is the assumption that the same $`f(p,t)`$ can be used in all three terms where it occurs. In fact in the acceleration flux it is the local distribution at the shock front, in the total number it is a volume averaged value, and in the loss term it is the downstream distribution which matters. Diffusion theory shows that in the steady state all three are equal, but this need not be the case in more elaborate transport models (Kirk, Duffy & Gallant kirk96 (1996)).

Substituting for $`\mathrm{\Phi }`$ and simplifying we get the equation

$$L\frac{f}{t}+𝐧𝐔_\mathrm{𝟏}f(p)+\frac{1}{3}𝐧\left(𝐔_\mathrm{𝟏}𝐔_\mathrm{𝟐}\right)p\frac{f}{p}=\frac{Q}{4\pi p^2}$$
(10)   
which is our version of the “Box” equation. Note that this, as is readily seen, gives the well known standard results for the steady-state spectrum and the acceleration time-scale. In fact our description is mathematically equivalent to that of Protheroe and Stanev as is easily seen by noting that

$$r_{\mathrm{acc}}=\frac{𝐧(𝐔_\mathrm{𝟏}𝐔_\mathrm{𝟐})}{3L},r_{\mathrm{esc}}=\frac{𝐧𝐔_\mathrm{𝟐}}{L},N=4\pi p^2fL.$$
(11)   
However our version has more physical content, in particular the two rates are derived and not inserted by hand. It is also important to note that in our picture the size of the “box” depends on the particle energy.

## 4 Inclusion of additional loss processes

In itself the “box” model would be of little interest beyond providing a simple “derivation” of the acceleration time scale. Its main interest is as a potential tool for investigating the effect of additional loss processes on shock acceleration spectra. One of the first such studies was that of Webb, Drury and Biermann (webb84 (1984)) where the important question of the effect of synchrotron losses was investigated (see also Bregman et al. bregman81 (1981)). An interesting question is whether or not a “pile-up” occurs in the accelerated particle spectrum at the energy where the synchrotron losses balance the acceleration. Webb, Drury and Biermann (webb84 (1984)) found that pile-ups only occured if the spectrum in the absence of synchrotron losses (or equivalently at low energies where the synchrotron losses are insignificant) was harder than $`fp^4`$. However Protheroe and Stanev obtain pile-ups for spectra as soft as $`fp^{4.2}`$.

It is relatively straightforward to include losses of the synchrotron or inverse Compton type (Thomson regime) in the model. These generate a downward flux in momentum space, but one which is distributed throughout the acceleration region. Combined with the fact that the size of the “box” or region normally increases with energy this also gives an additional loss process because particles can now fall through the back of the “box” as well as being advected out of it (see Fig. 1). Note that particles which fall through the front of the box are advected back into the acceleration region and thus this process does not work upstream.

If the loss rate is $`\dot{p}=\alpha p^2`$ the basic equation becomes

$`{\displaystyle \frac{}{t}}\left[4\pi p^2fL\right]+{\displaystyle \frac{}{p}}\left[\mathrm{\Phi }4\pi p^4f(p)\alpha L\right]=`$   
$`QU_24\pi p^2f(p)4\pi \alpha p^4f(p){\displaystyle \frac{dL_2}{dp}}`$ (12)   
This equation is easily generalised to the case of different loss rates upstream and downstream. Simplifying equation (12) gives

$`L{\displaystyle \frac{f}{t}}+p{\displaystyle \frac{f}{p}}\left[{\displaystyle \frac{U_1U_2}{3}}\alpha pL\right]+`$   
$`f\left[U_14\alpha pL\alpha p^2{\displaystyle \frac{dL_1}{dp}}\right]={\displaystyle \frac{Q}{4\pi p^2}}.`$ (13)   
Note that for convenience we have dropped the explicit vector (and tensor) notation; all non-scalar quantities are to be interpreted as normal components, that is $`U_2`$ means $`𝐧𝐔_\mathrm{𝟐}`$ etc. Note also that our model differs from that of Protheroe and Stanev in that they do not allow for the extra loss process resulting from the energy dependence of the “box” size.

In the steady state and away from the source region this gives immediately the remarkably simple result for the logarithmic slope of the spectrum,

$$\frac{\mathrm{ln}f}{\mathrm{ln}p}=3\frac{U_14\alpha pL\alpha p^2{\displaystyle \frac{dL_1}{dp}}}{U_1U_23\alpha pL}.$$
(14)   
Note that at small values of $`p`$ we recover the standard result, that the power-law exponent is $`3U_1/(U_1U_2)`$.

Under normal circumstances both $`L_1`$ and $`L_2`$ are monotonically increasing functions of $`p`$. Thus both the numerator and denominator of the above expression, regarded as functions of $`p`$, have single zeroes at which they change sign. The denominator goes to zero at the critical momentum

$$p^{}=\frac{U_1U_2}{3\alpha L}$$
(15)   
where the losses exactly balance the acceleration. If the numerator at this point is negative, the slope goes to $`\mathrm{}`$ and there is no pile-up. However the slope goes to $`+\mathrm{}`$ and a pile-up occurs if

$$U_14U_2+3\alpha p^2\frac{dL_1}{dp}>0\text{at}p=p^{}.$$
(16)   
In the early analytic work of Webb et al the diffusion coefficient was taken to be constant, so that $`dL_1/dp=0`$ and this condition reduces to $`U_1>4U_2`$ in agreement with their results. However if, as in the work of Protheroe and Stanev, the diffusion coefficient is an increasing function of energy or momentum, the condition becomes less restrictive. For a power-law dependence of the form $`Kp^\delta `$ the condition for a pile-up to occur reduces to

$$U_14U_2+\delta \left(U_1U_2\right)\frac{L_1}{L_1+L_2}>0$$
(17)   
(The equivalent criterion for the model used by Protheroe and Stanev is slightly different, namely

$$U_14U_2+\delta \left(U_1U_2\right)>0$$
(18)   
because of their neglect of the additional loss process.)

For the case where $`L_1/L_2=U_2/U_1`$ and with $`\delta =1`$ this condition predicts that shocks with compression ratios greater than about $`r=3.45`$ will produce pile-ups while weaker shocks will not. In Figures 1 and 2 we plot the particle spectra up to $`p^{}`$ for a range of values of $`\delta `$ and with $`r=4`$ and $`r=3`$ respectively.

Thus there is no contradiction between the (exact) results of Webb et al and those of Protheroe and Stanev; the apparent differences can be attributed to the energy dependence of the diffusion coefficient. Indeed, looking at the results presented by Protheroe and Stanev, it is clear that the pile-ups they obtain are less pronounced for those cases with a weaker energy dependence.

## 5 Nonlinear effects

At the phenomenological and simplified level of the “box” models it is possible to allow for nonlinear effects by replacing the upstream velocity with an effective momentum-dependent velocity $`U_1(p)`$, reflecting the existence of an extended upstream shock precursor region sampled on different length scales by particles of different energies. Higher energy particles, with larger diffusion length scales, sample more of the shock transition and have larger effective values of $`U_1(p)`$; thus $`U_1(p)`$ must be a monotonically increasing function of $`p`$. Repeating the above analysis with a momentum-dependent $`U_1`$ the logarithmic slope of the spectrum is in this case

$$\frac{\mathrm{ln}f}{\mathrm{ln}p}=3\frac{U_14\alpha pL+{\displaystyle \frac{p}{3}}{\displaystyle \frac{dU_1}{dp}}\alpha p^2{\displaystyle \frac{dL_1}{dp}}}{U_1U_23\alpha pL}$$
(19)   
with a pile-up criterion of,

$$U_1(p)4U_2p\frac{dU_1}{dp}+3\alpha _1p^2\frac{dL_1}{dp}>0\text{at}p=p^{}$$
(20)   
We see that whether or not the nonlinear effects assist the formation of pile-ups depends critically on how fast they make the effective upstream velocity vary as a function of $`p`$. By making $`U_1(p^{})`$ larger they make it easier for pile-ups to occur. On the other hand, if the variation is more rapid than $`U_1p`$, the derivative term dominates and inhibits the formation of pile-ups.

In most cases the shock modification will be produced by the reaction of accelerated ions, and the electrons can be treated as test-particles with a prescribed $`U_1(p)`$. However in a pair plasma, or if one applies the “box” model to the ions themselves, the effective upstream velocity has to be related to the pressure of the accelerated particles in a self-consistent way. We require in the ”box” model a condition which describes the reaction of the accelerated particles on the flow. Throughout the upstream precursor and in the steady case both the mass flux, $`A\rho U`$, and the momentum flux, $`AU+P_C`$ are conserved. Here $`P_C`$ is the pressure contained in energetic particles and the gas pressure is assumed to be negligible upstream. At a distance $`L_1(p)`$ upstream only particles with momenta greater than $`p`$ remain in the acceleration region. This suggests that in the ”box” model the reaction of the particles on the flow is described by the momentum flux conservation law

$$AU_1(p)+_p^{p_{\mathrm{max}}}4\pi p^2f\frac{pv}{3}𝑑p=\mathrm{constant}$$
(21)   
where $`p_{\mathrm{max}}`$ is the highest momentum particle in the system and $`v`$ is the particle velocity corresponding to momentum $`p`$. Differentiating with respect to $`p`$ gives

$$A\frac{dU_1(p)}{dp}=4\pi p^2f(p)\frac{pv}{3}.$$
(22)   
With no losses and for $`U_1(p)U_2`$ we can now recover Malkov’s spectral universality result for strong modified shocks (Malkov, malkov98 (1998)). In the limit of $`U_2=0`$ and $`\alpha =0`$ the conservation equation reduces to the requirement than the upward flux in momentum space be constant (equation (9)),

$$\mathrm{\Phi }=\frac{4\pi p^3}{3}f(p)U_1(p)=\mathrm{\Phi }_0.$$
(23)   
When combined with equation (22) this gives

$$U_1\frac{dU_1}{dp}=\frac{\mathrm{\Phi }_0}{A}v=\frac{\mathrm{\Phi }_0}{A}\frac{dT}{dp}$$
(24)   
where we have used the elementary result from relativistic kinematics that the particle velocity $`v`$ is the derivative of the kinetic energy $`T`$ with respect to momentum. Integrating for relativistic particles, $`T=pc`$, we get the fundamental self-similar asymptotic solution found by Malkov,

$$U_1=\sqrt{\frac{2c\mathrm{\Phi }_0}{A}}p^{1/2},f=\frac{3}{4\pi }\sqrt{\frac{\mathrm{\Phi }_0A}{2c}}p^{3.5}.$$
(25)   
If the electrons are test-particles in a shock strongly modified by proton acceleration, and if the Malkov scaling $`U_1p^{1/2}`$ holds even approximately, then equation (20) predicts that a strong synchrotron pile-up appears inevitable.

It is perhaps worth remarking on some peculiarities of Malkov’s solution. Formally it has $`U_2=0`$, all the kinetic energy dissipated in the “shock” is used in generating the upwards flux in momentum space $`\mathrm{\Phi }`$ and there is no downstream advection. It is not clear that a stationary solution exists in this case. The problem is that as $`U_20`$ so $`L_2\mathrm{}`$ if a diffusion model is used for the downstream propagation. The solution appears to require some form of impenetrable reflecting barrier a finite distance downstream if it is to be realised in finite time. Also, although the accelerated particle spectrum at the shock is a universal power law, none of these particles escape from the shock region. From a distance the shock appears as an almost monoenergetic source at whatever maximum energy the particles reach before escaping from the system.

The case of a synchrotron limited shock in a pure pair plasma is also interesting. Here the upper cut-off is determined not by a free escape boundary condition but by the synchrotron losses. If most of the energy dissipated in the shock is radiated this way, the shock will be very compressive and the downstream velocity $`U_2`$ negligible compared to $`U_1`$. The same caveats about time scales apply as to Malkov’s solution, but again we can, at least as a gedanken experiment, consider a cold pair plasma hitting an impenetrable and immovable boundary. In this case, if there is a steady solution, the upward flux due to the acceleration must exactly balance the synchrotron losses at all energies. In general it appears impossible to satisfy both this condition and the momentum balance condition for $`p<p^{}`$ unless the diffusion coefficient has an artificially strong momentum dependence. However a solution exists corresponding, in the box model, to a Dirac distribution at the critical momentum $`p^{}`$. This steady population of high energy electrons has enough pressure to decelerate the incoming plasma to zero velocity and radiates away all the absorbed energy as synchrotron radiation. This extreme form of pile-up may be of interest as a means of very efficiently converting the bulk kinetic energy of a cold pair plasma into soft gamma-rays.

## 6 Conclusion

A major defect of all “box” models is the basic assumption that all particles gain and loose energy at exactly the same rate. It is clear physically that there are very large fluctuations in the amount of time particles spend in the upstream and downstream regions between shock crossings, and thus correspondingly large fluctuations in the amount of energy lost. The effect of these variations will be to smear out the artificially sharp pile-ups predicted by the simple “box” models. However our results are based simply on the scaling with energy of the various gain and loss processes together with the size of the acceleration region. Thus they should be relatively robust and we expect that even if there is no sharp spike, the spectrum will show local enhancements over what it would have been in the absence of the synchrotron or IC losses in those cases where our criterion is satisfied.

## 7 Acknowledgments

This work was supported by the TMR programme of the EU under contract FMRX-CT98-0168. Part of the work was carried out while LD was Dozor Visiting Fellow at the Ben-Gurion University of the Negev; the warm hospitality of Prof M Mond and the stimulating atmosphere of the BGU is gratefully acknowledged.