{'fraction_non_alphanumeric': 0.03419426227296035, 'fraction_numerical': 0.010930289045421423, 'mean_word_length': 4.459555031678208, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 10, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Transient Emission From Dissipative Fronts in Magnetized, Relativistic Outflows. I. Gamma-Ray Flares

The transient emission produced behind internal shocks that are driven by overtaking collisions of a magnetized, relativistic outflow is considered. A self-consistent model capable of describing the structure and dynamics of the shocks and the time evolution of the pair-and gamma-ray distribution functions is developed and applied to gamma-ray flares in blazars, in the case in which gamma-ray production is dominated by inverse Compton scattering of external radiation (ERC). The dependence of the flare properties on magnetic field dissipation rate, intensity of ambient radiation, and the thickness of expelled fluid slabs is analyzed. It is shown that i) the type of gamma-ray flare produced by the model is determined by the ratio of the thickness of ejected fluid slab and the gradient length scale of ambient radiation intensity, ii) the radiative efficiency depends sensitively on the opacity contributed by the background radiation, owing to a radiative feedback, and is typically very high for parameters characteristic to the powerful blazars, and iii) the emitted flux is strongly suppressed at energies for which the pair-production optical depth is initially larger than unity; the time lag and flare duration in this energy range increase with increasing gamma-ray energy. At lower energies, flaring at different gamma-ray bands occurs roughly simultaneously, but with possibly different amplitudes. Some observational consequences are discussed.


Introduction
The discovery of strong, variable gamma-ray sources identified with blazars by CGRO has motivated many theoretical investigations concerning the emission from relativistic jets, and several models of gamma-ray blazars have recently been developed (e.g., Burns & Lovelace 1982;Dermer & Schlickeiser 1993;Bloom & Marscher 1993;Manneheim 1993;Sikora, et al. 1994;Blandford & Levinson 1995 [BL95]; Ghisellini & Madau 1996). However, most of these efforts have been devoted to examine spectral properties of blazars using steady-state models, which are most suitable for exploring quiescent states. Although some progress in identifying the emission mechanisms in individual sources (e.g., Sambruna et al. 1997;Sikora et al. 1997), and constraining the structure and dynamics of relativistic jets on small scales (Levinson 1996) has been made, it has became widely recognized that quantitative analysis of the growing body of variability data is crucial for advancing our understanding of the nature of relativistic jets further. Recent observational efforts (e.g., Reich et al. 1993;Maraschi et al. 1994;Wagner 1996;Buckley et al. 1996;Takahashi et al. 1996 ;Mattox et al. 1997;Wehrle et al. 1997; for a recent review see Ulrich 1997) to characterize the transient emission in blazars motivate detailed analysis of various variability mechanisms, and the development of time dependent models. Different processes may lead to time variability of blazar emission, including sudden changes in particle injection rate and/or magnetic field, changes in the bulk speed, and temporal variations of the intensity of background radiation in ERC models. This paper investigates the possibility that gamma-ray flares observed in blazars are produced by internal shocks propagating in a magnetized, relativistic jet. The work presented in this paper extends an earlier work by Romanova & Lovelace (1997;hereafter RL97). An outline of dissipative fronts is given in §2. In §3 we present the model and derive the basic equations. The results are described in §4. We conclude in §5.

Dissipative fronts in magnetized, relativistic outflows
Temporal fluctuations in the parameters of a MHD outflow lead to excitation of waves that steepen into shocks at some distance from the fluid ejection point, and the ultimate formation of a front consisting of a pair of shocks and a contact discontinuity across which the total pressure (kinetic plus magnetic) is conserved. The front propagates at a speed intermediate between that of the colliding streams and expands at a rate proportional to the relative velocity between the two shocks. Consequently, there is a net energy flow into the front which is balanced by either adiabatic cooling, owing to the front expansion, or radiative cooling, depending upon the conditions in the front. In the case of Poynting flux dominated outflows, rapid magnetic filed dissipation is required in order for a significant fraction of the outflow energy to be converted into kinetic energy behind the shocks. Efficient magnetic field dissipation may be achieved in strongly turbulent MHD flows (Thompson 1997) (an ordered magnetic field component and an associated Poynting flux can still be defined), or in the case of a sudden reversal of magnetic field lines (Romanova & Lovelace 1992).
We envision that the Lorentz factor of a fluid expelled from the central engine increases suddenly from Γ + to Γ − > Γ + , and denote by τ acc the characteristic time change of the outflow parameters (of order the dynamical time in the injection region). The ejection of the fast fluid is assumed to persist for time τ inj , after which it abruptly terminates, and the outflow Lorentz factor drops to lower values. As described above, this would lead to the creation of forward and reverse shocks that will propagate across the slow and fast fluid slabs, respectively. (The drop in Lorentz factor after time τ inj will lead also to the formation of a rare-faction wave behind the rapid slab.) The position r o at which the front is created has been calculated by Levinson & van Putten (1997) using simple-wave analysis. As they have shown, in the limit in which the disturbance speed with respect to the rest frame of the boundary at which the fluid is injected (henceforth referred to as injection frame; this might be e.g., the rest frame of the central engine ) is highly relativistic, r o = κcτ acc Γ 2 + Γ 2 A , where Γ A is the Lorentz factor associated with the Alfven speed with respect to the fluid rest frame, and κ is a numerical factor that depends on the details of fluid injection (κ = 1/3 when the 4-velocity of the expelled fluid increases according to u = sinh[λ + t/τ ]; t ≥ 0). In the case of extragalactic jets, Γ + is typically of order a few. Γ A may range from ∼1, in the case of a weakly magnetized flow, to ≥ 10 beneath the annihilation radius in Poynting flux jets (BL95; Levinson 1996). Thus, if we associate cτ acc with the gravitational radius of the putative black hole, we anticipate r o to lie in the range between 10 15 and 10 18 cm in blazars.
The creation of fronts in non-relativistic, hydrodynamic flows has been considered by Raga et al. (1990) in the context of HH objects. The basic idea has later been generalized and applied to Poynting flux jets in blazars by RL97. These authors proposed a scenario in which the gamma-ray and synchrotron outbursts often seen in blazars are ascribed to cooling of relativistic electrons accelerated inside highly dissipative fronts propagating in a magnetically dominated jet. The physical processes included in their model are inverse Compton scattering of external radiation, synchrotron and SSC emission. In their treatment, however, the spectrum of emitting particles is assumed a priori, and is characterized by two break energies which are allowed to evolve with time. The structure of the front is not calculated self-consistently, but rather the front is assumed to expand freely with the speed of sound. As we show, the front structure can be significantly altered by radiative losses. Moreover, their model does not account for absorption of gamma-rays escaping the front by pair production on ambient photons ahead of the front. As shown below, this process may considerably affect the emitted spectrum, particularly in the powerful blazars.
In this paper we generalize RL97 treatment in several important ways. Firstly, the time evolution of the energy distribution of pairs and gamma-rays in the front is calculated self-consistently, by numerically integrating the coupled kinetic and MHD equations. Secondly, the effect of radiative losses on the evolution of the front structure is taken into account. Thirdly, pair cascades inside and ahead of the front are accounted for properly. Fourthly, finite length outbursts are considered. And finally, we use different initial conditions. Our model involves the following key parameters: i) the dissipation rate of magnetic energy, ii) the maximum injection energy of electrons, denoted by E emax , iii) the fraction η of dissipation energy that is tapped to the injection of electrons to E emax , and iv) the ratio of the thickness of ejected fluid slab and the gradient length scale of the intensity of background radiation. We suppose that on scales of interest the outflow material is dominated by e ± plasma, and ignore possible effects associated with entrainment of ions (which are likely to be present in the MHD wind confining the jet) on the front dynamics. In the present work, we consider only gamma-ray flares produced by the ERC mechanism. The extension of this work to encompass synchrotron and SSC flares will be presented in a follow-up paper.

Front dynamics
The stress-energy tensor of a magnetized flow can be written as a sum of two contributions: the energy-momentum tensor associated with the electromagnetic field, T µν em , and the energy-momentum tensor associated with the system of particles, T µν s . T µν s and the particle flux, N µ , are related to the electron distribution function, defined in §4.2 below through, and where n is the proper density and u α the fluid 4-velocity. In terms of the proper gas pressure, p, specific enthalpy, h, magnetic induction, b α =F αβ u β , withF αβ being the dual electromagnetic tensor, and rest mass density ρ = m e n, the stress-energy tensor takes the form Here g αβ denotes the metric tensor, h * = h + b 2 /(ρc 2 ), p * = p + b 2 /2, and b 2 = b α b α .
Let Q and S ν denote the rate of change of the electron density due to e ± pair creation and annihilation, and the rate of change of the ν component of the stress-energy tensor resulting from radiative losses, respectively, and Q β b denotes the source term associated with magnetic field dissipation. The fluid equations can then be expressed as, The terms Q and S µ are expressed in §4.3 below in terms of the Boltzmann collision operators. The evaluation of these source terms at any given time involves the integration of the kinetic equations (eqs. [14] and [15]) for the pairs and photons.
We consider the evolution of a front produced by a collision of fast and slow fluids having velocities β − and β + , and corresponding Lorentz factors Γ − and Γ + , respectively. We suppose that each shock decays after it crosses the corresponding fluid slab, at which point energy deposition into that slab ceases. At times shorter than the shock crossing times, the set of equations (4) can be integrated over a portion of space-time containing the front, using the generalized Gauss's theorem, to obtain a set of ordinary differential equations governing the time evolution of the front quantities. The derivation of the front equations is given in the appendix. In what follows, subscript minus (plus) refers to quantities leftward (rightward) of the contact discontinuity, and subscript f denotes quantities inside the front. For simplicity, we shall restrict our analysis to the case in which the colliding fluids have the same proper density, pressure, and magnetic field, viz., n + = n − , p + = p − , b+ = b−. In this case a symmetric front will be created, i.e., The problem is then characterized by seven independent variables: where β s± are the corresponding shock velocities, and β c is the front velocity (i.e., the velocity of the contact discontinuity). Since the proper energy density of the fast and slow fluids is the same, the energy flux of the fast outflow, as measured in the injection frame, is larger than that of the slow outflow by roughly a factor of (Γ − /Γ + ) 2 . Thus, this case corresponds to a situation wherein the central engine undergoes an outburst during which energy is being injected into the system. An alternative possibility, which will not be considered here, is that the fluid parameters change in such a way as to sustain the associated energy flux unchanged. This would lead to the creation of an asymmetric front, which in some circumstances can give rise to a somewhat higher radiative efficiency. To simplify the analysis further, we suppose that the proper density, pressure and magnetic field in the front are homogeneous, that d + = d − = d, and that the magnetic field dissipation rate is proportional to the injection frame magnetic field, B f , and inversely proportional to the front length ∆X, viz., The latter assumption is not crucial, but is convenient since it admits steady-state solutions in the adiabatic case (i.e., in the absence of radiative losses) which simplifies the choice of initial conditions (cf §3.5). As we have verified, the essential features of the solution are insensitive to the choice of the form of the term associated with magnetic field dissipation. With the above simplifications eqs. (A7) reduce to where x o = ct, t being the time measured in the injection frame, and where the shock and front velocities obey the algebraic equations, The last two terms on the RHS of equations (5)-(8) account for the flux of the corresponding quantity incident into the front through the shock surfaces, whereas the second term on the RHS describes the rate of change of the corresponding front quantity due to the front expansion. Finally, the expansion of the front is governed by the equation, and the crossing times of the forward and reverse shocks, denoted by t R± , are given implicitly by As discussed earlier, the above equations are valid only for times shorter than min {t R+ ,t R− }. The equations describing the evolution of the system at later times can be derived in a similar manner. Formally, when t > t R+ we replace β s+ by β c in eqs. (5)-(11) , and set N + , B + , and T µν In doing so we ignore the rare-faction waves that will be produced at the edges of the slabs due to the pressure gradient there.

Kinetic equations
Let f e (p e , x) and f γ (p γ , x) denote the distribution functions of electrons (we do not distinguish electrons from positrons, designating both by subscript e) and gamma-rays as measured in the injection frame, respectively. The corresponding Boltzmann equations can be written as and The operators C e and C γ on the RHS of the above equations represent the interaction between pairs and gamma-rays, and are given explicitly in BL95. They are restricted to the condition by virtue of energy conservation. The injection operator, C inj , represents the change in energy and momentum of electrons (positrons) due to their interaction with the large scale electromagnetic field; that is, as a result of shock acceleration, magnetic reconnection, or stochastic acceleration owing to absorption of nonlinear plasma waves. Since the injection operator preserves the total electron number density, it must satisfy, The angular distribution of pairs and gamma-rays is expected to be strongly beamed along the direction of propagation of the front in the injection frame. Thus, we can approximate the distribution functions , with µ α being the cosine of the angle between the momentum of a species α and the flow velocity, and n α (E α ) the corresponding number density of a species α per unit energy. Integrating eqs. (14) and (15) over a region encompassing the front in the manner described in the appendix, assuming again that the distribution functions are homogeneous inside the front, and averaging over µ α yields to a good approximation, Now, the electron density in the fluids exterior to the front, viz., n e+ and ne−, are given as input, and can be absorbed into the definition of the injection operator. For convenience we redefine the injection operator as follows: This condition simply reflects the fact that particle acceleration preserves the total number density of particles incident into the front. A second constraint onC inj is obtained by multiplying eq. (18) by E e , and then integrating over dE e : where T oo sf is given by eq. (2). Note that the RHS of the last equation equals the fraction of total energy flux dissipated inside the front.
We must also determine n γ+ and n γ− . Under the approximation of perfect beaming invoked above n γ− = 0 and n γ+ = n γf , since the only source of beamed radiation is the front itself. With the above results eq. (19) simplifies to, Again, we assume that for t > t R± the time change of the distribution functions is dictated by eqs. (18) and (22) with β s− and β s+ replaced by β c and n e± set to zero.
The emergent gamma-ray spectrum depends on the pair production opacity contributed by the ambient radiation field ahead of the front. Gamma-rays for which the pair production optical depth to infinity largely exceeds unity, will be converted into lower energy pairs and gamma-ray via pair cascades upstream the forward shock. (We note that the escape probability of gamma-rays at a given energy from the front may be considerably different than the probability that they will escape the system to infinity. In fact, the pair production optical depth of the front itself will be smaller than unity as long as its axial length is smaller than the corresponding gamma-spheric radius.) This leads, as shown below, to a strong suppression of the emitted flux at the corresponding energies. The calculation of the emitted spectrum proceeds as follows: At every time step we integrate the equations describing the spatial evolution of the cascade, starting at the current position of the forward shock, r(t) = r o + β s+ dt, and subject to the boundary conditions n e [r(t), E e ] = n e+ , n γ [r(t), E γ ] = n γf (t, E γ ). Here n γ (r, E γ ) is the number density per unit energy of escaping gamma-rays ahead of the forward shock. The assumption underline these calculations is that the cascade develops on a sufficiently short time scale to render any retardation effects negligible. We also ignore possible alterations of the front structure due to momentum exchange between emitted gamma-rays and the slow outflow (i.e., upstream of the forward shock). The inclusion of this effect complicates the numerics substantially and is beyond the scope of our analysis. We anticipate, though, that such alterations of the front structure will not affect the emission characteristics significantly.

Determination of the source functions
The source term describing the rate of change of the electron density in the front due to formation of pair cascades, is determined by integrating eq. (18) over energy, and using eqs. (20) and (1). One then obtains, where κ pp (x, E γ ) is the pair production opacity, and the second equality has been obtained using eq. (4.12) of BL95. Likewise, the radiative energy loss rate is obtained by taking the first moment of eq. (18), i.e., multiplying the equation by E e and then integrating over dE e , and by using eqs. (21) and (2): Under the assumption of perfect beaming the momentum loss rate equals the energy loss rate, viz., S x = cS o . By employing eqs. (16) and (22) we can rewrite the rates associated with radiative losses in the form where T oo γ = E γ n γ d ln E γ is the gamma-ray energy density. Note that the latter equation is essentially the energy equation for the system of gamma-rays.

The injection operator
As stated above, the injection operator represents the energy redistribution of electrons (positrons) resulting from their interaction with the electromagnetic field, e.g., due to Fermi acceleration, magnetic reconnection, or absorption of nonlinear plasma waves excited by the shock. The physics of these processes is yet to be understood better before C inj can be derived from first principles. Here we settle for a simple prescription in which a small fraction of the pairs incident into the front are injected to some maximum energy, assumed to be fixed in the observer frame, and the rest are redistributed at much lower energies. The maximum injection energy, E emax , and the fraction η of the dissipation energy that is being injected to E emax are treated as free parameters. To be concrete, we choosē with n 1 , n 2 and E o determined from eqs. (20) and (21), and the requirement that n 2 E emax comprises a fraction η of the total dissipation energy. Other forms ofC inj can be readily derived, for example, a thermal distribution with a power law tail. We find though that the results are insensitive to the exact form of the injection operator provided that particle acceleration is sufficiently efficient, in the sense that a significant fraction of the injection energy is carried by electrons having energies near E emax . It should be noted that the assumption that E emax is fixed in the observer frame is unrealistic. A more realistic prescription would be to take E emax to be fixed in the rest frame of the front. Unfortunately, this complicates the numerics considerably. However, variations of E emax due to the changing front velocity should not affect significantly the evolution of the gamma-ray distribution function at sufficiently low energies.

Initial condition
In order to integrate the above equations, one must specify the initial position and structure of the front, and the initial energy distribution of electrons. A simple and convenient choice would be to use the structure of an adiabatic front. Below, we employ the front structure calculated by Levinson & Van Putten (1997). Specifically, we first solve algebraically eqs. (5)-(11) in the absence of radiative losses (i.e., with Q = S o = S x = α b = 0), and assuming steady state (∂/∂x o = 0), for the chosen input parameters. The front quantities thereby obtained then serve as initial values for the integration of the radiative front equations. This choice is viable if the creation of the front occurs on a timescale much shorter than the radiative cooling time. An initial electron distribution subject to the restrictions that the number density and total energy (temperature) equal those computed for the adiabatic front is also specified. The initial thickness of the front, ∆X(t = 0), is taken to be small enough, so that the energy loss rate due to cooling of the initial population of electrons is well below the rate of energy deposition in the front. This renders the results highly insensitive to the choice of initial electron distribution.

Results
Equations (5)-(13), (18), (22)-(24) and (28) have been integrated numerically using the initial conditions discussed in §3.5. The source functions Q and S o have been computed using eqs. (25) and (26). We have verified that the solution is indeed highly insensitive to the choice initial electron distribution when the initial length of the front is sufficiently small, as stated in §3.5. The equations have been modified in a manner described in §3.1 after each shock crossing (first at t = min {t R+ , t R− } and again at t = max {t R+ , t R− }), and the integration continued using the modified equations. In the following examples, the Lorentz factors of the slow and fast fluids and the rest frame Alfven 4-velocity have been taken to be 5, 20, and 10, respectively, and the initial electron distribution has been taken to be Maxwellian with appropriate temperature and density. A rapid magnetic field dissipation has been invoked with α b = 0.5. The standard soft photon intensity (consisting of a broken power law with a steeper slope below about 0.5 KeV) defined in  has been adopted for the calculations. To be concrete, we used the following form for the intensity of external radiation: Here E o is the break energy and is taken to be 0.5 keV, and ǫ << 1 is the fraction of nuclear luminosity that is reprocessed or scattered by surrounding gas. As a check, we integrated eq. (27) to obtain T oo γ and compared the result with total flux emitted (i.e., E γ n γ d ln E γ ) at every time step. The agreement was typically better than 2%.
To illustrate the dynamics of the system we first consider the limit of continuous ejection; that is, d/r o → ∞. The time evolution of the front quantities is depicted in fig. 1, where the 4-velocity of the contact discontinuity and the two shocks (upper left panel), the proper density (upper right panel), the total pressure (bottom left panel), and the proper magnetic pressure (bottom right panel) are plotted against log(ct/r o ), t being the injection frame time, for d/r o = 100, (ǫL s ) 45 /r 16 = 1 (solid lines) and (ǫL s ) 45 /r 16 = 10 −2 (dashed lines). Here (ǫL s ) 45 is the scattered luminosity in units of 10 45 erg s −1 , and r 16 = r o /(10 16 cm). The corresponding radiative efficiency, defined as the fraction of dissipation power radiated by the front, and given explicitly by, is exhibited in fig. 2. As seen from figs. 1 and 2, the front velocity and expansion rate, initially equal those of an adiabatic front, decrease as the emitted energy flux rises. The rest mass density and total pressure increase correspondingly. After the peak emission is reached, the front starts accelerating and adiabatic cooling becomes gradually more important until, ultimately, the initial structure and velocity of the front are restored. The increase in number density is partly due to pair production and partly due to enhanced compression resulting from the drop in expansion rate. The change in total pressure is primarily due to magnetic field compression. Quite generally, we find that when t R± > l/c, the maximum of the flux is reached at a distance l from the creation radius r o , where l is the gradient length scale of the intensity of ambient radiation at r o (l ∼ r o in this example), and the duration of the flare corresponds to several times l (see below). When t R± < l/c the timescale of the flare is determined by the shock crossing times. The evolution of the system in this case is essentially the same as described above for times shorter than t R± . At later times, however, energy supply to the hot outflow terminates, and the outflow begins to cool radiatively and decelerate until radiative losses become small. The radiative efficiency, namely the ratio of radiated and dissipated energy fluxes increases with increasing values of (ǫL s ) 45 /r 16 , as seen from fig. 2. It approaches 80 percent near the peak in this example for (ǫL s ) 45 /r 16 = 1, and is generally around this value for typical parameters of gamma-ray blazars.
The rate of energy dissipation in the front is presented in fig. 3, and it is seen that it increases as the radiated flux increases. The reason is that radiative losses lead to deceleration of the front and shocks (see fig. 1) and, consequently, enhancement in the rate of energy deposition behind the shocks, particularly the reverse shock, as can be inferred from eq. (7). This positive feedback has important implications for the emitted spectrum, particularly in the presence of large pair production opacity, that are discussed below. Fig. 4 depicts typical light curves computed for d/r o > 1 and different values of (ǫL s ) 45 /r 16 . In this figure the total apparent gamma-ray luminosity is plotted against the time measured by a distant observer, where β s+ is the velocity of the forward shock, and θ is the angle to the line of sight (θ = 0 in this example). The strong dependence of the flare intensity on the parameter (ǫL s ) 45 /r 16 is evident. This is one consequence of the positive feedback mentioned above. The decay of the flare in this regime depends on the radial variation of the intensity of background radiation and the expansion rate of the front, and is typically longer and more gradual than the rise. Further, the rise time increases slightly with increasing L 45 /r 16 . The reason is that larger opacity results in smaller shock and front velocities during peak emission (see fig. 1) and, therefore, smaller beaming factor. From fig. 4. it is seen that the rise time is of order (l/2c)Γ −2 s+ (t peak ) ≃ 10 −2 r o /c, where Γ s+ (t peak ) is the Lorentz factor of the forward shock near maximum flux, and that the flare duration is several times longer. For r o between 10 15 and 10 18 cm (cf. §2) this corresponds to a rise time in the range between several minutes and several weeks. However, gamma-rays of a given energy cannot escape from a radius smaller than their gamma-spheric radius, so that the duration of a flare in a given band is also limited by the pair production opacity. In general, we find that for parameters typical to the powerful blazars the flare duration in the EGRET band can range from several hours to several weeks for small viewing angles (i.e., θ < Γ −1 ), when d/r o > 1. The shape of the flare is substantially altered when d/r o becomes sufficiently small. This is demonstrated in Fig. 5 where light curves computed for L 45 /r 16 = 1 and different values of d/r o are displayed. As seen, for values of d/r o smaller than unity the time scale of the flare, as measured by a distant observer, is of order d/c; the light curves exhibit a steep decline at times t obs > d/c. For sufficiently thin slabs (d/r o < 2 × 10 −2 in this example) the decay time of the flare is comparable to the rise time or even slightly shorter.
In order to examine the relationship between the emergent fluxes in different high-energy bands predicted by our model, we divided the energy interval into several subintervals and followed the time evolution of the flux in each band. Quite generally we find that i) only a relatively small fraction of the radiated energy is emitted above the gamma-spheric energy at which the pair-production optical depth at r ∼ r o + l (the front position at peak emission) equals unity, ii) the time of peak emission and flare duration increase with increasing gamma-ray energy in this energy range, and iii) below that energy variations in the fluxes at different X-and gamma-ray bands occurs roughly simultaneously, but not necessarily with the same amplitude. This behavior is also a consequence of the positive feedback discussed above. An example is given in fig. 6, where the time evolution of the energy fluxes emitted in four equally spaced subintervals in the energy interval 5 MeV to 50 GeV, for d/l > 1 is presented. Similar behavior has been found in the case d/l << 1. This result is in contrast with the simultaneous flaring predicted by the model of RL97, which does not take into account gamma-ray attenuation by pair production on external photons.
We also checked the dependence of the front dynamics and emission on the rate of magnetic field dissipation. We find that as U A increases above unity, the front expands more rapidly, adiabatic cooling becomes more important and, hence, the radiative efficiency drops. The rapid expansion results in a slower decay of the flare, but otherwise the shape of the light curve is not altered significantly.

Conclusions
This paper considers the production of gamma-ray flares through inverse Compton scattering of external radiation by pairs accelerated behind internal shocks, that are driven by temporal fluctuations in the parameters of a magnetized, relativistic outflow. A self-consistent model capable of describing the dynamics of the front and the time evolution of the angle averaged pair and gamma-ray distribution functions, has been developed and employed to calculate gamma-ray light curves in blazars.
The main results and conclusions of this study are: 1) The shape and timescale of the flare depend on the ratio of the thickness of ejected fluid slab, d, and the gradient length scale of the background radiation intensity at the radius of shock formation, l. When d/l is sufficiently small, such that the propagation time of the shocks across the fluid slab is shorter than l/c, the flare duration, as measured by an observer at small viewing angle, is ∼ d/c, and the shape of the light curve is roughly symmetric. For larger values of d/l the flare reaches its maximum at a distance l from the creation radius. It then decays gradually until time d/c, after which the flux declines steeply. For d/l > 1 the shape of the flare is determined essentially by the radial variations of the intensity of ambient radiation and the expansion rate of the front; typically, the rise is fast compared with the decay. The above results suggest that different types of flares will be produced in different sources, or even in the same source under different conditions. This expectation appears to be consistent with the variety of types of EGRET flares (e.g., Hartman et al. 1993;Kniffen et al. 1993;Mattox et al. 1997) observed in gamma-ray blazars. In the event of impulsive outflow ejection (i.e., ejection time of order the dynamical time, τ inj ≃ τ acc ; cf. §2) the axial length of the fluid slab can be as small as the size of the central engine. If the latter is associated with the gravitational radius of the putative black hole, r g , then flare durations as short as r g /c ∼ 50[M BH /(10 7 M ⊙ )] sec are plausible. Successive ejection of jets with short duty cycle may also lead to a substructure in the light curve of a longer duration outburst (as often seen in GRBs). In general we find that for parameters typical to extragalactic jets, namely Γ of order 10 and comoving Alfven 4-velocity between 0 and 10, the flare timescale can lie in the range between several minutes and several weeks (for small viewing angles). We note that the shape of the light curve may be significantly altered when the SSC process becomes important (cf. RL97). The study of SSC flares is left for future investigation.
2) The peak flux and flare intensity depend sensitively on the Thomson opacity near r o , specifically, on the value of L s /r o , where L s is the fraction of nuclear luminosity that is reprocessed or scattered across the jet. This behavior is a consequence of a positive feedback that gives rise to a strong enhancement of the energy deposition rate during peak emission. For parameters typical to the powerful gamma-ray blazars, the radiative efficiency, namely the fraction of dissipation energy that is radiated by the front is typically high (> 60%).
3) The amplitude of variations of the flux emitted at energies for which the pair-production optical depth is initially larger than unity is much smaller than that at lower energies. The flare, in this energy interval, propagates from low to high energies and its duration increases with energy. At lower energies, flaring in different X-and gamma-ray bands occurs roughly simultaneously, but with possibly different amplitudes, depending on the spectra of ambient radiation and injected electrons. The claimed invariance of the gamma-ray spectrum during the strong EGRET flare observed in PKS1622-297 (Mattox et al. 1997) is not in conflict with that model prediction, since the sensitivity was insufficient to resolve any time lags, if present, between the fluxes in the two energy subintervals (below and above 300 MeV) analyzed by those authors, particularly in view of the small amplitudes anticipated at high EGRET energies. Despite the small amplitudes, it may be possible to observe time lags between the emission at hard (but still well below TeV where absorption by the IR background is strong) and soft gamma-ray energies with a next generation gamma-ray telescope. The detection of such lags is an important test for this model. One implication of the above results is that in some cases strong flares may be observed in some energy bands while the simultaneous variations in other bands may appear modest or even small. The simultaneous X-ray/TeV flare reported recently for Mrk 421 and the lack of significant variations of the EGRET flux (Macomb, et al. 1995;Takahashi et al. 1996) is an example. The hard X-ray flux may be due to either inverse Compton emission of the thermal electrons or synchrotron emission of the highest energy electrons. Alternatively, this event may be explained as an SSC flare resulting from changes in the tail of the distribution of non-thermal electrons (Takahashi et al. 1996).
I thank the referee, R.V.E. Lovelace, for useful comments. Support by Alon fellowship and a TAU Research Authority grant is acknowledged.

A. Appendix: Derivation of the front equations
Equations (4) are of the form, Integrating the above equation over a region in space-time enclosed by the world lines associated with the left shock, the contact discontinuity and the lines defined by the equations t = t o ; t = t o + ∆t, where t o is some fiducial time [see fig. (7)] , and using the generalized Gauss's theorem, where d 4 Ω is a volume element and d 3 Σ µ is the surface area of the 3d surface enclosing the region of integration yields,