Issue 
A&A
Volume 609, January 2018



Article Number  A122  
Number of page(s)  8  
Section  Extragalactic astronomy  
DOI  https://doi.org/10.1051/00046361/201732000  
Published online  02 February 2018 
Recollimation shocks and radiative losses in extragalactic relativistic jets
^{1} INAF–Osservatorio Astrofisico di Torino, Strada Osservatorio 20, 10025 Pino Torinese, Italy
email: bodo@to.astro.it
^{2} INAF–Osservatorio Astronomico di Brera, via E. Bianchi 46, 23807 Merate, Italy
Received: 26 September 2017
Accepted: 17 October 2017
We present the results of stateoftheart simulations of recollimation shocks induced by the interaction of a relativistic jet with an external medium, including the effect of radiative losses of the shocked gas. Our simulations confirm that – as suggested by earlier semianalytical models – the postshock pressure loss induced by radiative losses may lead to a stationary equilibrium state characterized by a very strong focusing of the flow, with the formation of quite narrow nozzles, with crosssectional radii as small as 10^{3} times the length scale of the jet. We also study the timedependent evolution of the jet structure induced by a density perturbation injected at the flow base. The setup and the results of the simulations are particularly relevant for the interpretation of the observed rapid variability of the γray emission associated to flat spectrum radio quasars. In particular, the combined effects of jet focusing and Doppler beaming on the observed radiation make it possible to explain the subhour flaring events such as that observed in the flat specrum radio quasar PKS 1222+216 by MAGIC.
Key words: galaxies: jets / galaxies: active / radiation mechanisms: nonthermal / shock waves
© ESO, 2018
1. Introduction
Extragalactic relativistic jets – collimated outflows of plasma expelled by supermassive black holes residing in central regions of active galactic nuclei and traveling for hundreds of kpc – are among the most fascinating astrophysical structures. Their phenomenology is regulated by complex timedependent physical processes involving magnetic fields and plasma. Shocks and/or magnetic reconnection sites are able to accelerate charged particles to ultrarelativistic energies, whose presence is flagged by intense nonthermal emission over the entire electromagnetic spectrum, up to the TeV band (e.g., Romero et al. 2017, for a review). Relativistic jets are best studied in blazars (see e.g., Urry & Padovani 1995). The jet of these sources is closely aligned toward the line of sight and, thanks to this favorable geometry, their nonthermal emission is strongly amplified by the relativistic beaming effects.
The observation of the emission from blazar jets at the highest energies, accessible thanks to space (Fermi) and groundbased (Cherenkov arrays) instruments, is revealing numerous unexpected features. One of the most intriguing and challenging aspects is the ultrafast variability (doubling time ≲ 1 h, down to a few minutes) detected in several blazars (Aharonian et al. 2007; Albert et al. 2007; Arlen et al. 2013; Aleksić et al. 2011, 2014). These very small timescales are often smaller than the lightcrossing time of the horizon of the supermassive black hole (e.g., Vovk & Babić 2015) and their interpretation requires extreme physical conditions (e.g., Begelman et al. 2008). In particular, very rapid variability events detected in the most powerful blazars, the flat spectrum radio quasars (FSRQ), are the most challenging to interpret. FSRQ display the classical broad emission lines of quasars, flagging the existence of a broad line region in which clouds of photoionized gas reprocess part of the intense continuum emitted by the accretion flow. Due to the anticipated absorption of the γray radiation above a few tens of GeV through the interaction with soft (UV) photons (e.g., Liu & Bai 2006), one can put robust lower limits on the distance of the emitting region from the central black hole since, to avoid absorption, the emission is constrained to occur beyond the BLR radius R_{BLR}, of the order of ~ 0.1 pc. For a conical geometry of the jet (with semiaperture θ_{j} ≈ 0.1) this lower limit on the distance readily translates into a lower limit for the jet radius and (through lightcrossing time argument) to the minimum variability timescale, s (for which we use a Doppler factor^{1}), a value clearly not compatible with the most extreme minutetimescale events. One of the first FSRQ events showing rapid variability was PKS 1222+122 which, during a flare in June 2010, was detected by the MAGIC telescopes with a flux increasing with a doubling timescale of about 10 min (Aleksić et al. 2011). Since the discovery of PKS 1222+122 other FSRQ have been observed to vary on short timescales. Recently, 3C 279 was detected during a pointed Fermi/LAT observation to vary with a timescale of about 5 min (Ackermann et al. 2016).
In the framework of the classical shockinjet scenario for particle acceleration (e.g., Aller et al. 1985) it is clearly hard to explain very short flares. At first sight, the only possibility to explain the rapid variability in FSRQ seems to require the existence of compact (R ≲ 10^{15} cm) emission regions embedded in the jet at distances of the order of ~0.1 pc from the central engine (e.g., Ghisellini & Tavecchio 2008; Tavecchio et al. 2011). Generally, such compact regions have been identified in active substructures of the jet, such as plasmoids resulting from efficient reconnection of the magnetic field (Giannios 2013; Petropoulou et al. 2016) or in turbulent cells embedded in the relativistic flow (Narayan & Piran 2012; Marscher 2014). These scenarios have been widely discussed in the past (see also Aharonian et al. 2017). However, one can adopt a different view and explore the possibility that the emission occurs within the entire (or a large fraction of the) jet. This could be reconciled with the fast variability if the jet were to suffer a strong recollimation interacting with an external medium. Jet recollimation has been widely studied in the past (e.g., Komissarov & Falle 1998; Sokolov et al. 2004; Stawarz et al. 2006; Nalewajko & Sikora 2009). Bromberg & Levinson (2009) extended the calculations including the effect of strong radiative losses of the shocked plasma. The cooling, causing a loss of pressure in the compressed region, allows the jet to be “squeezed” by large factors. Bromberg & Levinson applied their scenario to the case of BL Lac objects and radiogalaxies but – as proposed in Tavecchio et al. (2011) – it can also be applied to FSRQ, for which the cooling of the plasma (dominated by the inverse Compton scattering of the ambient soft photons (e.g., Ghisellini et al. 1998) is expected to be even more severe than for other sources.
In this paper we intend to study the recollimation of a jet subject to important radiative losses extending the previous approximated semianalytical approach of Bromberg & Levinson (2009), by performing highresolution axisymmetric numerical simulations using the stateoftheart numerical code PLUTO (Mignone et al. 2007, 2012). The structure of the paper is as follows. In Sect. 2 we illustrate the model and the setup used, in Sect. 3 we show the results of simulations for stationary states and for the timedependent evolution of a density perturbation, and in Sect. 4 we conclude.
2. The model
We want to construct a numerical model of the confinement of a relativistic jet by the pressure of the surrounding medium. We assume axisymmetry and make use of the cylindrical coordinates r and z, where z is the coordinate along the jet axis and r is the radial coordinate. We inject a conical relativistic jet, with opening angle θ_{0}, in an overpressured confining region with a pressure profile (1)and a density profile (2)where z_{0} represents the height at which the jet starts to interact with the confining region.
The jet is injected in the confining region (z_{0} = 1) with a Lorentz factor γ_{j} and a radius r_{j} = θ_{0}z_{0}, and is subject to radiative losses by nonthermal processes. In this paper we do not consider the dynamical effects of a magnetic field, therefore the evolution equations determining the jet dynamics are the continuity equation (3)where ρ is the proper density and u^{μ} is the fourvelocity, and the energy momentum conservation (4)where is the stressenergy tensor of a perfect gas (5)where w and p denote, respectively, the proper enthalpy and pressure, and g^{μν} is the metric tensor for a flat space. S^{μ} denotes a source term associated with radiative losses, whose form is specified below. The system of Eqs. (3) and (4) is completed by providing an equation of state relating w, ρ, and p. Following Mignone et al. (2005), we adopt the following prescription: (6)which closely reproduces the thermodynamics of the Synge gas for a singlespecies, relativistic perfect fluid with a smooth transition from the adiabatic exponent Γ = 5/3 in the nonrelativistic limit to Γ = 4/3 in the ultra relativistic case (here and in the following we always put c = 1).
For radiative losses we assume a very simple form; they are taken to be proportional to the pressure and we additionally assume that it is only the gas above a certain threshold temperature T_{c} that radiates. The first assumption is widely adopted (e.g., Komissarov & Falle 1998; Bromberg & Levinson 2009). For the case under study here the cooling is dominated by the inverse Compton emission on a fixed target radiation field. The emissivity (the frequency integrated radiated power per unit volume) of the nonthermal particles measured in the plasma comoving frame, ϵ ≡ −S^{0}, can be expressed as: (7)where σ_{T} is the Thomson cross section, U_{ext} is the energy density of the external radiation field, γ is the electron Lorentz factor and n(γ) is the nonthermal electron energy distribution, for simplicity assumed to be a power law with slope 2, n(γ) = kγ^{2}, in the range γ_{min}<γ<γ_{max}, with γ_{min} ≪ γ_{max} (all quantities are expressed in the flow frame). Performing the integral one obtains ϵ = (4/3)σ_{T}cU_{ext}kγ_{max}. The energy density of the nonthermal relativistic electrons is: (8)where Λ = ln(γ_{max}/γ_{min}). Combining the two equations above, we can write: (9)where in the last step we use p_{e} = U_{e}/ 3 (valid for ultrarelativistic electrons) and we have defined an effective cooling time: (10)The pressure of the nonthermal particles is supposed to be a fraction ξ_{e} of the thermal one, that is, p_{e} = ξ_{e}p. We note that the relation ϵ ∝ p_{e} is valid if, as assumed here, the electron losses are determined by the inverse Compton emission on a fixed target radiation field. In case where the losses are dominated by synchrotron or synchrotronself Compton emission one should also consider the role of the magnetic field. The assumption on the critical temperature is intended to mimic the idea that high temperatures flag the presence of relativistic particles in the flow. A final point concerns the fact that (as also pointed out by Nalewajko & Sikora 2009) in a realistic treatment one should assume that relativistic nonthermal particles are accelerated close to the shock fronts and then advected in the other regions of the jet. Both the electron distribution and ξ_{e} should thus be treated as functions of the position in the flow. For simplicity this is not considered here. In any case, we expect that the global dynamical effects of cooling, such as those explored here, do not dramatically depend on these details. Similarly to Bromberg & Levinson (2009) the quantities characterizing the nonthermal population are assumed constant in the jet. Implicitly this choice assumes the existence of some mechanism (e.g., turbulence) continuously supplying the energy lost by particles through the emission.
With the assumptions above, and assuming that, in the rest frame of the gas, radiation is isotropic (for our applications we can safely neglect the anisotropy of the external Compton radiation field, Dermer (1995), Ghisellini & Tavecchio (2010)), we can express S^{μ} as (11)In the laboratory frame we can thus write (12)Equations (3) and (4) are solved numerically by using the adaptive mesh refinenment (AMR) version of the PLUTO code (Mignone et al. 2007, 2012), with a second order scheme and HLLC Riemann solver (Mignone & Bodo 2005). We perform twodimensional axisymmetric simulations using the cylindrical coordinates r and z on a grid that covers the domain 0 <r<r_{D}, 1 <z<z_{D}, where r_{D} = 1.5z_{0} and z_{D} = 9z_{0}. We make use of six levels of refinement, the base grid is made up of 64 × 384 points and we have an equivalent maximum resolution of 4096 × 24 576 points. The boundary conditions are reflective on the axis r = 0, outflow at the outer boundaries z_{D} = 9z_{0} and r_{D} = 1.5z_{0}. At the lower boundary, z_{0} = 1, and for r<r_{j} we have inflow conditions injecting the jet flow, while, for r>r_{j}, we have reflective conditions.
The jet is injected with a proper density ρ_{j} = 10^{6}ρ_{a}, where ρ_{a} is the ambient density at z = z_{0}, with a Lorentz factor γ_{j} = 10 and pressure is p_{j} = 7 × 10^{3}p_{a}. The corresponding energy flux is: (13)where w_{j}, according to Eq. (6), is (14)In the simulations we fix the parameters to values representative for the flaring state of PKS 1222+216 (e.g., Tavecchio et al. 2011). For the energy flux we use L_{j} = 10^{46} erg s^{1}. Using γ_{max} = 10^{5}, Λ ≃ 10 (quite insensitive on the exact value of γ_{min}) and U_{ext} = 3 × 10^{2} (suitable to reproduce the IR radiation field of the torus (e.g., Ghisellini & Tavecchio 2008) we obtain τ_{c} ≃ 3 × 10^{4} s. Since we assume that the jet is already beyond the BLR, we fix z_{0} to the BLR radius, z_{0} = 7 × 10^{17} cm. The cooling time is thus of the order of τ_{c} = 10^{3}z_{0}/c. We investigate the effects of radiative losses by performing a series of simulations with different values of ξ_{e} in the range 0.01−0.1.
3. Results
3.1. The equilibrium configuration
The characteristics of recollimation shocks resulting from the interaction of an underpressured relativistic jet with external confining material have been discussed in the past; in particular by Komissarov & Falle (1998), Bromberg & Levinson (2007, 2009), Nalewajko & Sikora (2009). The gross structure consists of a contact discontinuity separating the shocked jet layer and the ambient medium and a recollimation shock. Komissarov & Falle (1997) derived, under suitable approximations, simple analytical formulae for the geometrical properties of the shock. Assuming a power law profile for the pressure of the external gas, p_{ext}(z) = az^{− η}, the profile r(z) of the reconfinement shock follows the differential equation: (15)where δ = 1−η/ 2 and , where μ = 0.7. The solution of this equation, with the condition r = r_{0} at z = z_{0} is: (16)for η ≠ 2, and r(z) = zAln(z_{c}/z), for η = 2 (i.e., δ = 0) where z_{c} = z_{0}exp(r_{0}/z_{0}A). For the contact discontinuity surface a similar analytical expression can be derived (Bromberg & Levinson 2007), r_{c}(z) = r_{0}(z/z_{0})^{η/ 4}.
These analytical solutions do not include the possible cooling of the postshock plasma. Bromberg & Levinson (2007) and Bromberg & Levinson (2009) constructed a class of semianalytical models for the confinment of the jet by the external pressure including also the effect of the pressure loss due to radiation losses. As above, these solutions are subject to some approximations and limitations: in the shocked layer the flow parameters are assumed to depend only on the coordinate along the jet axis z (gradients along r have been included by Nalewajko & Sikora 2009) and, more critically, they do not get the structure of the solutions for values of z>z^{∗}, where z^{∗} is the point where the recollimation shock reaches the jet axis.
Here we use a different approach, by which we overcome the two limitations mentioned above, solving numerically the time dependent equations and following the jet evolution until it reaches a steady configuration. We start the simulations at z = z_{0} with a conical jet of opening angle θ_{0} established over the entire computational domain. As the simulation evolves in time, the jet is compressed by the higher pressure of the ambient medium, the recollimation shock is formed and the system evolves towards its steady structure, that is reached at about t = 300z_{0}/c. The main features of the steady solutions mentioned above can be observed in Fig. 1, where we show the Lorentz factor distribution for the case with the smaller value of ξ_{e} = 10^{2}. In the figure we see the region of unshocked jet material separated by the recollimation shock from the shocked jet layer. A contact discontinuity separates this shocked material from the external medium, which appears in black. The conical recollimation shock reaches the jet axis at z = z^{∗} ~ 3z_{0}, where it is reflected and gives rise to a strong jet deceleration. From the detailed view of the pressure distribution, shown in the left panel of Fig. 2, we can observe the formation of a complex shock structure, which leads to a pressure maximum for 3.6z_{0}<z< 3.8z_{0}. Having linked the emissivity of nonthermal particles to the value of the pressure (Eq. (11)), the region of maximum pressure corresponds also to a maximum of radiative losses. Therefore, in the region around the axis, the jet, having lost a fraction of its energy flux by radiation, will remain at a Lorentz factor quite a lot lower than its injection value (Fig. 1). This slower region on the axis is surrounded by a faster region, formed by material that has suffered less radiative losses.
Considering a jet propagating towards the observer with viewing angle θ_{v} = 1 /γ_{j}, the Doppler factor will be and, with the assumptions described above, the apparent (Doppler boosted) emissivity of each fluid element as measured by the observer is given by: (17)The location of the emissivity maximum depends both on the distribution of pressure and on the distribution of Lorentz factor within the flow. Since the Lorentz factor is lower on the axis and increases away from it, we expect to find the maximum of ϵ_{obs} away from the axis. In fact from the right panel of Fig. 2, where we display the distribution of the observed emissivity, we can see that most of the observed luminosity originates from an elongated region at some distance from the axis, where the product of the pressure and the Lorentz factor is maximized. We remark that, although the thickness of this region is very small (< 0.01 z_{0}), its length is of the order of ≈ 0.2 z_{0}.
Fig. 1
Distribution of the Lorentz factor for the case with ξ_{e} = 0.01, showing the main features of the steady solution. The red rectangle individuates the region shown in Fig. 2. We notice that the radial scale is strongly stretched. 
Fig. 2
Left panel: pressure distribution in the region marked by the red rectangle in Fig. 1 for the case with ξ_{e} = 0.01. The arrows indicate the main features already displayed in Fig. 1. Right panel: distribution of the observed emissivity (see Eq. (17) 
Fig. 3
Distribution of the Lorentz factors for the three cases with different values of ξ_{e}, that is ξ_{e} = 0.01 (left panel), ξ_{e} = 0.03 (middle panel), and ξ_{e} = 0.05 (right panel). 
The effects of increasing the radiative losses are illustrated by Fig. 3, where we show the distributions of the Lorentz factor for three cases with different values of ξ_{e}. From the Figure we can observe a progressive strong reduction in transverse size of the region of shocked jet material between the recollimation shock and the contact discontinuity, for increasing ξ_{e}. The “squeezing” of the jet also becomes more effective and the position of the minimum jet transverse radius comes closer to the point where the recollimation shock reaches the axis. The Lorentz factor on the axis is also more strongly decreased because of the increased energy losses. We note also that after a moderate reexpansion, following the development of the reflection shock, the jet is recollimated again.
Figure 4 shows a zoomed view of the distributions of pressure and observed emissivity for the case with ξ_{e} = 0.05. The pressure distribution shows an overall shape similar to the previous case, shown in Fig. 2, however it is much reduced in size and the maximum value is about one order of magnitude larger. Consequently, the emission also has a similar shape, with much reduced size and higher peak value.
In order to make this result more quantitative, in Fig. 5 we plot the radial width w of the emission region at the point where the value is half the maximum as a function of ξ_{e}. We observe a progressive decrease of the width for increasing values of the energy transferred to the nonthermal component. For the largest ξ_{e} = 0.05 the width is of the order of 10^{3}z_{0}, slightly smaller than the estimate of Bromberg & Levinson (2009). Therefore, even in nonextreme conditions, the focusing due to the recollimation is so effective that the brightest region of the jet is as small as ≈ 7 × 10^{14} cm, not far from the size estimated for the rapid variability observed in 1222+216 (Aleksić et al. 2011; Tavecchio et al. 2011).
Fig. 4
Left panel: pressure distribution in the region marked by the red rectangle in Fig. 1 for the case with ξ_{e} = 0.05. Right panel: distribution of the emissivity in the observer frame (see Eq. (17). 
For illustration, we also report in Fig. 6 the profiles (integrated in planes normal to the jet axis) along the jet axis of the observed emissivity (left panel) (18)where ϵ_{obs} is defined in Eq. (17) and of the jet energy flux (right panel) (19)for the three values of ξ_{e}. Both quantities are normalized to the jet power. The overall shape of the integrated emissivity L(z) (left panel) is the same for all values of the losses, displaying a broad maximum around z_{0} ~ 1.5 in correspondence with the onset of the recollimation shock and a second, narrower maximum marking the region of minimum transverse size. For the lowest value of ξ_{e} = 10^{2} (black line) this second peak is relatively broad and not very high. For larger values of ξ_{e} (red and green lines) the overall emissivity increases and the second maximum becomes quite narrow, flagging the compactness of the emission region. We note that, in all cases, the area under the two peaks is of the same order, indicating that the total received luminosity comprises similar contribution from the large (hence slowly varying) recollimation shock and the compact (rapidly variable) region. The profiles of the energy flux (right panel) show the increasing importance of energy losses for increasing values of ξ_{e}. After a smooth decrease, marking the progressive development of the recollimation shock, the curves (especially the two with ξ_{e} = 0.03 and 0.05) display a sudden jump, corresponding to the increased losses close to the region of smaller radius. After the expansion produced by the reflected shock, the energy flux remains approximately constant. The fraction of the jet energy flux lost through the emitted radiation goes from 15% for ξ_{e} = 0.01 to a fraction of about 30% for ξ_{e} = 0.05. These values are in agreement with the results of BL09.
Fig. 5
Transverse size of the emission region as a function of ξ^{e}. 
Fig. 6
Left panel: observed emissivity integrated over planes normal to the z axis as a function of the z coordinate. Right panel: energy flux integrated over planes normal to the z axis as a function of the z coordinate. 
It is also interesting to consider the efficiency ϵ_{diss} of the dissipation of the kinetic flux of the jet that we define, similarly to Nalewajko & Sikora (2009), as L_{kin} = ^{∫}ρc^{2}v_{z}γ_{j}(γ_{j}−1), where the integral is performed on planes normal to the z direction. Comparing the values of L_{kin} at injection and after the recollimation region, we derive ϵ_{diss} = 18% for ξ_{e} = 0.01 and ϵ_{diss} = 35% for ξ_{e} = 0.05. These values are larger than the average shock dissipation efficiency derived by Nalewajko & Sikora (2009), that found values lower than 10% for γ_{j}θ_{j} = 1. The reason for this higher dissipation efficiency should be looked for in the important radiation losses of the system, not included in the Nalewajko & Sikora treatment. In fact, radiation losses lead to a lower value of the gas pressure downstream from the recollimation shock. In order to maintain pressure balance with the external material, the system has thus to increase the fraction of jet luminosity dissipated into pressure. This is indeed naturally achieved through a narrower profile of the recollimation shock which increases the dissipation because of the larger component of the velocity normal to the shock.
3.2. Evolution of perturbations and variability
In the previous section we studied in detail the stationary equilibrium structure of the recollimation shock. For astrophysical applications, it is important to have a glimpse on the temporal behavior of the structure and its properties when the jet is subject to perturbations.
Fig. 7
Pressure distribution (logarithmic scale) for the simulation including a propagating density perturbation. The three panels show the pressure in the stationary (equilibrium) state, after a time t = 0.9z_{0}/c and, on the left, at t = 5z_{0}/c. The figure clearly shows the highpressure wake induced by the perturbation in the jet. In the last panel, the inward shift of the recollimation nozzle is clearly visible. 
Fig. 8
Profiles (integrated in planes normal to the jet axis) along the jet axis of the observed emissivity calculated for different times during the propagation of the perturbation. In both panels the red lines display the unperturbed profile. The left panel shows the profile a short time after the injection, while in the right panel we report the profiles at later times. 
In order to study the variability behavior, we injected a perturbation at the jet inlet and followed its evolution as propagated along the jet. We choose to change the jet density, keeping a constant temperature, with a perturbation of the form: (20)where ϵ is the relative amplitude of the perturbation, τ is the temporal width of the perturbation and t_{0} is the time at which the perturbation reaches its maximum value. From now on we take t_{0} as the origin of time. For τ we took the value 0.05z_{0}/c and we explored different values of ϵ. Our results show that only quite large values of ϵ give appreciable effects on the observed emissivity; for example ϵ = 1 gives only a 10% increase in the emissivity, therefore we discuss the results obtained for ϵ = 10. In Fig. 7 we show the logarithmic pressure distributions at three different times. In the left panel we have the equilibrium distribution. The middle panel is for t = 0.9z_{0}/c and we can see the propagating perturbation that leaves, on its sides, a high pressure wake that compresses and modifies the structure of the recollimation shock. Finally the right panel is for t = 5z_{0}/c, the perturbation is outside the region displayed in the Figure, however, its wake is still visible and the recollimation shock is still deformed; the point at which the recollimation shock reaches the axis has moved inward.
In Fig. 8 we show the profiles (integrated in planes normal to the jet axis) along the jet axis of the observed emissivity (Eq. (17)) for different times. The left panel shows the first phases of propagation (corresponding to t = 0.21,0.67and1.35z_{0}/c), with the red curve showing the unperturbed profile for comparison. The peak corresponding to the perturbation decreases its amplitude as it moves outward. This decrease is related to the jet expansion. Because of the higher pressure and the shock deformation in the wake of the perturbation (see Fig. 7), the inner jet maintains an emissivity that is approximately twice the unperturbed value even at larger times. In the right panel, the red curve always shows the unperturbed distribution for comparison, while the blue curve, corresponding to the time at which the perturbation reaches the recollimation region (t = 1.8z_{0}/c), shows a burst of emissivity, with an observed luminosity increasing by a factor of ≈ 5 with respect to the stationary case. Finally, the green curve in the right panel corresponds to a time when the perturbation has moved outside the displayed region, strongly decreasing its amplitude. In this case the peak in the emissivity at the recollimation point has moved inward, following the shift of the recollimation nozzle visible in the maps in Fig. 7.
These results, although obtained under rather idealized assumptions, outline a coherent scenario, qualitatively compatible with the behavior of powerful blazars during active states. In fact, the light curves of FSRQ during typical flare episodes show a longterm (weeks) increase of the flux (that in our simulations would be provided by the long duration active state of the perturbed inner jet), on top of which shorttimescale (days) variations are often detected. These rapid variations could be related to the excitation of the recollimation nozzle by the incoming perturbation. In real astrophysical conditions we could also imagine a series of perturbations with different properties injected into the flow that would account for the erratic behavior of the real light curves.
4. Discussion and conclusions
We have developed a series of axisymmetric simulations of a relativistic jet subject to recollimation caused by an external medium using the stateoftheart code PLUTO and including the radiative losses mediated by a population of nonthermal electrons in the jet. The set up of the simulations is tailored to model jets associated to FSRQ, in particular the source PKS 1222+216, displaying rapid variability at VHE. Our simulations overcome several of the approximations introduced in the semianalytical approach of Bromberg & Levinson (2009).
The radiative cooling is then assumed to be dominated by the inverse Compton scattering of external radiation (Sikora et al. 1994; Ghisellini & Tavecchio 2009). The estimated effective timescale for radiative losses, τ_{c}, is much lower then the dynamical timescale τ_{dyn} ~ z_{o}/c, ensuring the effective cooling of the shocked plasma. In order to allow a sizeable energy loss by the jet we have also assumed that the emitting electrons, while cooling, are continuously heated, efficiently converting the energy dissipated by the shock into radiation. We have shown that, even with a relatively minor contribution of the nonthermal particles to the total pressure of the shocked material, p_{e}/p = ξ_{e} = 0.05, the fraction of energy flux lost by the flow can be as large as 30%, resulting in a significant deceleration of the flow downstream the recollimation zone.
The significant pressure loss caused by the radiative cooling allows the jet to be focused and extremely “squeezed” by the external pressure. Consistently with the early findings by Bromberg & Levinson (2009), our simulations show that the minimum radius of the recollimation region can be as small as 10^{3}z_{0}, implying a light crossing time comparable to the timescale observed for PKS 1222+216. We have also produced maps for the observed emissivity of the flow, including the effects of the Doppler beaming. These show that the emission is concentrated in a very thin elongated region close to the jet axis. We need to caution that the release of the axisymmetry constraint may lead to instabilities of this narrow region, whose outcome can however be assessed only through threedimensional simulations.
We have also explored the behavior of the recollimation shock in response to density perturbations injected at the base of the jet. Although highly idealized, these simulations show a rather complex phenomenology, that could be, at least qualitatively, related to the typical flaring pattern of powerful FSRQ.
We would like to mention that our calculations, similarly to those of BL09, should be considered somewhat optimistic. In fact the radiative losses induced in the system are maximized by the specific form of the cooling term adopted in the simulations, simply proportional to the thermal gas pressure. As already noted in Sect. 2, such prescription implicitly assumes that the emitting electrons are continuously heated by some unspecified mechanism, allowing to radiate an amount of energy larger than that put into the nonthermal component just after the shock front. As discussed by BL09, possible heating mechanisms include nonlinear oscillations of the contact discontinuity surface or acceleration by shear flows. The investigation of these crucial aspects could be the subject of future specific simulations.
Acknowledgments
This work has been partly funded by a PRININAF 2014 grant, the PRIN INAF SKACTA grant and a PRIN MIUR 2015. We acknowledge the CINECA award under the ISCRA initiative, for the availability of highperformance computing resources.
References
 Ackermann, M., Anantua, R., Asano, K., et al. 2016, ApJ, 824, L20 [Google Scholar]
 Aharonian, F., Akhperjanian, A. G., BazerBachi, A. R., et al. 2007, ApJ, 664, L71 [Google Scholar]
 Aharonian, F. A., Barkov, M. V., & Khangulyan, D. 2017, ApJ, 841, 61 [NASA ADS] [CrossRef] [Google Scholar]
 Albert, J., Aliu, E., Anderhub, H., et al. 2007, ApJ, 669, 862 [NASA ADS] [CrossRef] [Google Scholar]
 Aleksić, J., Antonelli, L. A., Antoranz, P., et al. 2011, A&A, 530, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aleksić, J., Antonelli, L. A., Antoranz, P., et al. 2014, A&A, 563, A91 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aller, H. D., Aller, M. F., & Hughes, P. A. 1985, ApJ, 298, 296 [NASA ADS] [CrossRef] [Google Scholar]
 Arlen, T., Aune, T., Beilicke, M., et al. 2013, ApJ, 762, 92 [NASA ADS] [CrossRef] [Google Scholar]
 Begelman, M. C., Fabian, A. C., & Rees, M. J. 2008, MNRAS, 384, L19 [NASA ADS] [CrossRef] [Google Scholar]
 Bromberg, O., & Levinson, A. 2007, ApJ, 671, 678 [NASA ADS] [CrossRef] [Google Scholar]
 Bromberg, O., & Levinson, A. 2009, ApJ, 699, 1274 [NASA ADS] [CrossRef] [Google Scholar]
 Dermer, C. D. 1995, ApJ, 446, L63 [NASA ADS] [CrossRef] [Google Scholar]
 Ghisellini, G., & Tavecchio, F. 2008, MNRAS, 386, L28 [NASA ADS] [CrossRef] [Google Scholar]
 Ghisellini, G., & Tavecchio, F. 2009, MNRAS, 397, 985 [NASA ADS] [CrossRef] [Google Scholar]
 Ghisellini, G., & Tavecchio, F. 2010, MNRAS, 409, L79 [NASA ADS] [CrossRef] [Google Scholar]
 Ghisellini, G., Celotti, A., Fossati, G., Maraschi, L., & Comastri, A. 1998, MNRAS, 301, 451 [NASA ADS] [CrossRef] [Google Scholar]
 Giannios, D. 2013, MNRAS, 431, 355 [NASA ADS] [CrossRef] [Google Scholar]
 Komissarov, S. S., & Falle, S. A. E. G. 1998, MNRAS, 297, 1087 [NASA ADS] [CrossRef] [Google Scholar]
 Liu, H. T., & Bai, J. M. 2006, ApJ, 653, 1089 [NASA ADS] [CrossRef] [Google Scholar]
 Marscher, A. P. 2014, ApJ, 780, 87 [NASA ADS] [CrossRef] [Google Scholar]
 Mignone, A., & Bodo, G. 2005, MNRAS, 364, 126 [NASA ADS] [CrossRef] [Google Scholar]
 Mignone, A., Plewa, T., & Bodo, G. 2005, ApJS, 160, 199 [NASA ADS] [CrossRef] [Google Scholar]
 Mignone, A., Bodo, G., Massaglia, S., et al. 2007, ApJS, 170, 228 [Google Scholar]
 Mignone, A., Zanni, C., Tzeferacos, P., et al. 2012, ApJS, 198, 7 [Google Scholar]
 Nalewajko, K., & Sikora, M. 2009, MNRAS, 392, 1205 [NASA ADS] [CrossRef] [Google Scholar]
 Narayan, R., & Piran, T. 2012, MNRAS, 420, 604 [NASA ADS] [CrossRef] [Google Scholar]
 Petropoulou, M., Giannios, D., & Sironi, L. 2016, MNRAS, 462, 3325 [NASA ADS] [CrossRef] [Google Scholar]
 Romero, G. E., Boettcher, M., Markoff, S., & Tavecchio, F. 2017, Space Sci. Rev., 207, 5 [NASA ADS] [CrossRef] [Google Scholar]
 Sikora, M., Begelman, M. C., & Rees, M. J. 1994, ApJ, 421, 153 [NASA ADS] [CrossRef] [Google Scholar]
 Sokolov, A., Marscher, A. P., & McHardy, I. M. 2004, ApJ, 613, 725 [NASA ADS] [CrossRef] [Google Scholar]
 Stawarz, Ł., Aharonian, F., Kataoka, J., et al. 2006, MNRAS, 370, 981 [NASA ADS] [CrossRef] [Google Scholar]
 Tavecchio, F., BecerraGonzalez, J., Ghisellini, G., et al. 2011, A&A, 534, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Urry, C. M., & Padovani, P. 1995, PASP, 107, 803 [NASA ADS] [CrossRef] [Google Scholar]
 Vovk, I., & Babić, A. 2015, A&A, 578, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
All Figures
Fig. 1
Distribution of the Lorentz factor for the case with ξ_{e} = 0.01, showing the main features of the steady solution. The red rectangle individuates the region shown in Fig. 2. We notice that the radial scale is strongly stretched. 

In the text 
Fig. 2
Left panel: pressure distribution in the region marked by the red rectangle in Fig. 1 for the case with ξ_{e} = 0.01. The arrows indicate the main features already displayed in Fig. 1. Right panel: distribution of the observed emissivity (see Eq. (17) 

In the text 
Fig. 3
Distribution of the Lorentz factors for the three cases with different values of ξ_{e}, that is ξ_{e} = 0.01 (left panel), ξ_{e} = 0.03 (middle panel), and ξ_{e} = 0.05 (right panel). 

In the text 
Fig. 4
Left panel: pressure distribution in the region marked by the red rectangle in Fig. 1 for the case with ξ_{e} = 0.05. Right panel: distribution of the emissivity in the observer frame (see Eq. (17). 

In the text 
Fig. 5
Transverse size of the emission region as a function of ξ^{e}. 

In the text 
Fig. 6
Left panel: observed emissivity integrated over planes normal to the z axis as a function of the z coordinate. Right panel: energy flux integrated over planes normal to the z axis as a function of the z coordinate. 

In the text 
Fig. 7
Pressure distribution (logarithmic scale) for the simulation including a propagating density perturbation. The three panels show the pressure in the stationary (equilibrium) state, after a time t = 0.9z_{0}/c and, on the left, at t = 5z_{0}/c. The figure clearly shows the highpressure wake induced by the perturbation in the jet. In the last panel, the inward shift of the recollimation nozzle is clearly visible. 

In the text 
Fig. 8
Profiles (integrated in planes normal to the jet axis) along the jet axis of the observed emissivity calculated for different times during the propagation of the perturbation. In both panels the red lines display the unperturbed profile. The left panel shows the profile a short time after the injection, while in the right panel we report the profiles at later times. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.