PROBING THE INNERMOST REGIONS OF AGN JETS AND THEIR MAGNETIC FIELDS WITH RADIOASTRON. I. IMAGING BL LACERTAE AT 21 μas RESOLUTION

RadioAstron observations of BL Lac at 22.2 GHz (K band) were performed on 2013 November 10–11 (from 21:30 to 13:00 UT). A total of 26 ground antennas were initially scheduled for the observations, but different technical problems at the sites limited the final number of correlated ground antennas to 15, namely Effelsberg (EF), Metsaehovi (MH), Onsala (ON), Svetloe (SV), Zelenchukskaya (ZC), Medicina (MC), Badary (BD), and VLBA antennas Brewster (BR), Hancock (HN), Kitt Peak (KP), Los Alamos (LA), North Liberty (NL), Owens Valley (OV), Pie Town (PT), and Mauna Kea (MK).

The data were recorded in two polarizations (left and right circularly polarized, LCP and RCP), with a total bandwidth of 32 MHz per polarization, split into two intermediate frequency (IF) bands of 16 MHz. The SRT data were recorded by the RadioAstron satellite tracking stations (Kardashev et al. 2013; Ford et al. 2014) in Puschino (21:30 – 06:10 UT) and Green Bank (07:30 – 13:00 UT), including some extended gaps required for cooling of the motor drive of the on board high-gain antenna of the Spektr-R satellite. These gaps were used for the ground-only observations of BL Lac at 15 and 43 GHz (see Section 2.2), as well as observations of several calibrator sources.

Correlation of the data was performed using the upgraded version of the DiFX correlator developed at the Max-Planck-Institut für Radiostronomie (MPIfR) in Bonn (Bruni et al. 2014), enabling accurate calibration of the instrumental polarization of the SRT (see also Lobanov et al. 2015 for a more detailed description of the imaging and correlation of RadioAstron observations).

The correlated data were reduced and imaged using a combination of the AIPS and Difmap (Shepherd 1997) software packages. The a priori amplitude calibration was applied using the measured system temperatures for the ground antennas and the SRT. Sensitivity parameters of the SRT (Kovalev et al. 2014) are measured regularly during maintenance sessions. Parallactic angle corrections were applied to the ground antennas to correct for the feed rotation.

2.1.1. Fringe Fitting

Fringe fitting of the data was performed by first manually solving for the instrumental phase offsets and single band delays using a short scan during the perigee of the SRT, when the shortest projected ground–space baseline distances (smaller than one Earth diameter) are obtained, providing the best fringe solutions for the space antenna. These solutions were applied before performing a global fringe search for the delays and rates of the ground array only.

Once the ground antennas were fully calibrated, they were coherently combined to improve the fringe detection sensitivity of the SRT (Kogan 1996). This baseline stacking was carried out setting DPARM(1) = 3 in AIPS's FRING task, also performing an exhaustive baseline search and combining both polarizations and IFs to improve the sensitivity. Progressively longer solution intervals were used for the fringe search, from one minute to maintain the coherence of the signal during the acceleration of the spacecraft in the perigee to four minutes to increase the sensitivity on the longer baselines to the SRT.

The phasing of a group of ground-based antennas allowed us to obtain reliable ground-space fringe detections up to projected baselines of 7.9 Earth diameters (${D}_{{\rm{E}}}$) in length, covering the duration of the experiment within which Puschino was used as the tracking station. No further fringes were obtained to the spacecraft once the tracking station changed to Green Bank, which is presumably due to a difference in clock setting between the two tracking stations. These were searched for by introducing trial clock offsets for the Green Bank tracking station, and performing new test correlations with a larger fringe-search window of up to 1024 channels and 0.1 s of integration time in width. However, no further fringes were detected to the SRT. We also note that 1.5 hr passed between the last Puschino scan and the first Green Bank scan, thus increasing the space baseline length and perhaps reducing the correlated flux density below the sensitivity threshold. The resulting fringe-fitted data visibility coverage of the Fourier domain (uv-coverage) is shown in Figure 1.

Zoom In Zoom Out Reset image size

After the fringe fitting the delay difference between the two polarizations was corrected using AIPS's task RLDLY, and a complex bandpass function was solved for the receiver.

2.1.2. Polarization Calibration

Simultaneous ground-only observations of BL Lac at 15.4 and 43.1 GHz were obtained during gaps in the RadioAstron observations. Participating antennas were Effelsberg and VLBA antennas BR, HN, KP, LA, NL, OV, and PT. However, at 43 GHz, no fringes were detected on the intercontinental baselines with Effelsberg due to technical problems severely limiting the sensitivity and angular resolution, and in fact preventing the detection of polarization at this frequency. For this reason, we have used instead the 43 GHz data from the VLBA-BU-BLAZAR¹¹ monitoring program, performed on 2013 November 18, only one week after our observations.

Calibration of the 15 GHz data follows that outlined previously in Section 2.1, except for the particular steps related more directly to the space VLBI observations, such as the phasing of the ground array during the fringe fitting. Calibration of the absolute orientation of the EVPAs was also performed by comparison with Effelsberg single dish observations of BL Lac (${S}_{15\mathrm{GHz}}=5.57\pm 0.23$ Jy, $m=3.44\pm 0.3$%, EVPA = $3\buildrel{\circ}\over{.} 7\pm 1\buildrel{\circ}\over{.} 1$), with an estimated error of 5°–7°.

Figures 2 and 3 show the 15 and 43 GHz images of BL Lac obtained with the ground array. To characterize the emission structure, we have performed a fit of the complex visibilities by a set of components with circular Gaussian brightness distributions, listed in Table 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2.  Gaussian Model Fits for the 22 GHz RadioAstron and Ground-array Data at 43 and 15 GHz

Comp.	Flux	Distance	PA	Size	${T}_{{\rm{b}}}$	m	EVPA
	(mJy)	(mas)	(°)	(mas)	(K)	(%)	(°)
RadioAstron 22 GHz
Up	1578 ± 72	0.041 ± 0.003	−1 ± 4	0.050 ± 0.003	(1.56 ± 0.26) × 1012	4.0 ± 0.4	17 ± 1
Core	802 ± 43	...	...	<0.01	>2.0 × 1013	4.8 ± 0.4	16 ± 1
K1	1128 ± 37	0.164 ± 0.004	172 ± 1	0.067 ± 0.003	(5.82 ± 0.70) × 1011	5.7 ± 0.5	17 ± 1
K2	578 ± 12	0.320 ± 0.008	179 ± 1	0.122 ± 0.004	(9.79 ± 0.85) × 1010	8.6 ± 2.0	84 ± 5
K3	79 ± 12	0.79 ± 0.03	−161 ± 2	0.231 ± 0.012	(3.64 ± 0.93) × 109	...	...
K4	111 ± 13	1.26 ± 0.04	−170 ± 1	0.352 ± 0.018	(2.22 ± 0.49) × 109	...	...
K5	70 ± 13	1.66 ± 0.05	−173 ± 2	0.343 ± 0.017	(1.47 ± 0.42) × 109	...	...
K6	108 ± 12	2.39 ± 0.02	−176 ± 1	0.247 ± 0.012	(4.39 ± 0.91) × 109	20.8 ± 2.4	4 ± 4
K7	283 ± 14	2.57 ± 0.04	−171 ± 1	0.556 ± 0.056	(2.26 ± 0.57) × 109	26.2 ± 6.7	8 ± 3
K8	214 ± 15	3.37 ± 0.11	−168 ± 2	1.23 ± 0.12	(3.48 ± 0.92) × 108	...	...
VLBA-BU-BLAZAR Ground-array Data at 43 GHz
Core	1575 ± 72	...	...	0.026 ± 0.005	(1.48 ± 0.63) × 1012	1.7 ± 0.2	37 ± 6
Q1	1373 ± 65	0.091 ± 0.005	177 ± 8	0.048 ± 0.005	(3.92 ± 1.00) × 1011	2.0 ± 0.3	7 ± 4
Q2	786 ± 32	0.262 ± 0.005	173 ± 1	0.077 ± 0.007	(8.68 ± 1.93) × 1010	5.4 ± 0.1	18 ± 2
Q3	402 ± 18	0.378 ± 0.008	177 ± 1	0.105 ± 0.009	(2.40 ± 0.52) × 1010	3.1 ± 0.8	29 ± 4
Q4	50 ± 10	0.92 ± 0.07	−165 ± 5	0.239 ± 0.012	(5.79 ± 1.74) × 108	...	...
Q5	100 ± 11	1.51 ± 0.08	−170 ± 3	0.382 ± 0.019	(4.52 ± 0.95) × 108	...	...
Q6	209 ± 13	2.57 ± 0.09	−173 ± 2	0.562 ± 0.028	(4.34 ± 0.70) × 108	...	...
Ground-array Data at 15 GHz
Core	2705 ± 124	...	...	<0.06	>7.9 × 1012	2.0 ± 0.1	20 ± 1
U2	1193 ± 60	0.266 ± 0.019	−177 ± 5	<0.06	>3.5 × 1012	1.5 ± 0.1	30 ± 1
U3	266 ± 23	0.982 ± 0.025	−166 ± 2	0.523 ± 0.026	(5.03 ± 0.94) × 109	5.5 ± 1.7	67 ± 3
U4	84 ± 13	1.53 ± 0.02	−171 ± 1	0.235 ± 0.012	(7.87 ± 2.00) × 109	...	...
U5	462 ± 32	2.43 ± 0.02	−173 ± 1	0.503 ± 0.025	(9.44 ± 1.59) × 109	18.9 ± 4.0	5 ± 1
U6	332 ± 27	3.15 ± 0.04	−167 ± 1	1.03 ± 0.05	(1.62 ± 0.29) × 109	27.0 ± 18	16 ± 8
U7	391 ± 33	4.96 ± 0.12	−173 ± 2	2.84 ± 0.14	(2.52 ± 0.46) × 108	...	...

Note. Tabulated data correspond to component's label, flux density, distance and position angle from the core, size, observed brightness temperature, degree of linear polarization, and electric vector position angle, uncorrected for Faraday rotation (see Section 4).

Download table as:  ASCIITypeset image

Previous high angular resolution monitoring programs of BL Lac systematically show the presence of two stationary features close to the core (Stirling et al. 2003; Jorstad et al. 2005; Mutel & Denn 2005). In particular, Jorstad et al. (2005) report, through an analysis of a sequence of 17 bimonthly VLBA observations at 43 GHz, two stationary components, labeled A1 and A2, located at a mean distance from the core (position angle) of 0.10 and 0.29 mas ($-160^\circ$ and $-159^\circ$), respectively. These can be associated with components C2 and C3 in Mutel & Denn (2005), respectively, obtained from independent 43 GHz VLBA observations.

More recently, Cohen et al. (2014, 2015) present an analysis of more than a decade of 15 GHz VLBA observations of BL Lac from the MOJAVE monitoring program. These observations confirm the existence of a stationary component, labeled C7 by these authors, located at a mean distance from the core of 0.26 mas and at a position angle of −1666, in agreement with component A2 reported by Jorstad et al. (2005).

This is also corroborated by our 15 and 43 GHz observations (see Table 2 and Figures 2 and 3), in which components U2 and Q2 would correspond to previously identified components A2 and C7, and component Q1 to A1.

Our measured position angles for the two stationary features (Q1 and Q2/U2) are slightly offset to the east by ∼20° with respect to the main values published by Jorstad et al. (2005) and Cohen et al. (2014). This may be associated with the jet precession reported by Stirling et al. (2003) and Mutel & Denn (2005), leading to a swing in the position angle of the innermost components. A similar variation in the position angle of component U2 is seen in MOJAVE observations by Cohen et al. (2014).

Model fitting of the innermost structure in our RadioAstron observations, listed in Table 2 and plotted in Figure 5, shows two close components near the core region, as well as two other components within the innermost 0.5 mas region. The fitted circular Gaussian components provide an accurate representation of the jet emission, yielding a residual map (uniform weighting) with an rms of 1.4 mJy/beam and minimum and maximum residuals of −45 mJy/beam and 53 mJy/beam, respectively. Allowing for elliptical Gaussian components provides very similar fitted values as those listed in Table 2 with no significant improvement in the residuals.

We can tentatively identify components K1 and K2 with the previously discussed stationary features Q1 and Q2 (U2). However, identification of the two-component structure in the core area requires a more detailed analysis of the evolution of the innermost structure of the jet.

For this, we have performed model fitting of the 43 GHz VLBA images from the VLBA-BU-BLAZAR program extending our analysis to cover the period between 2013 December and 2014 June, comprising a total of five more epochs with a cadence of roughly one month. The obtained model fits are listed in Table 3, and Figure 7 plots the component's distance from the core versus time. Stationary components Q1 and Q2 are detected at all epochs. Analysis of their flux densities show that Q1 becomes unusually bright (2.1 Jy) on 2013 December 16, followed by a similar increase in flux density of Q2 (1.1 Jy) on 2014 January 19. This could be associated with the passing of a new component, M1, through the standing features, leading to a brief increase in their flux densities (e.g., Gómez et al. 1997). Component M1 is identified in the last two epochs as a weak knot beyond 0.5 mas from the core. The estimated apparent speed for this new component is $7.9\pm 0.3\;c$ (1.76 ± 0.06 mas yr⁻¹), in agreement with values previously found in BL Lac (Jorstad et al. 2005), giving an ejection date of 2013.89 ± 0.05, or 2013 November 23 (±18 days), in coincidence with our RadioAstron observations within the errors.

Zoom In Zoom Out Reset image size

Considering the measured proper motion of the new component and the estimated time of crossing through the 43 GHz core, we can estimate that component M1 should be placed ∼50 μas upstream of the core in our RadioAstron image, or slightly smaller if we allow for some initial acceleration. Based on this, we consider that the core of the jet in our RadioAstron observations corresponds to the unresolved component labeled "Core," upstream of which component "Up" would correspond to component M1 identified also at 43 GHz. Evidence for emission upstream of the core is also found by looking at the 43 GHz polarization image (see Figure 3), which shows polarized emission to the north with a different orientation of the polarization vectors than the remaining core area.

We note that the appearance of new superluminal components upstream of the core in BL Lac has been previously reported by Marscher et al. (2008), associated with a multiwavelength outburst. Similarly, the radio, optical, and γ-ray light curves (see the VLBA-BU-BLAZAR web page¹² ) show a flare at the end of 2013, close in time to our RadioAstron observations.

Table 2 also lists the observed (i.e., not corrected for redshift or Doppler boosting) brightness temperatures, estimated from the model-fitted circular Gaussian components as ${T}_{{\rm{b}}}=1.36\times {10}^{9}S{\lambda }^{2}/{\theta }^{2}$ [K], where S (in Jansky) is the total flux density, θ (mas) is the size, and λ (cm) is the observing wavelength (e.g., Kovalev et al. 2005). Model-fitted data of the RadioAstron visibilities yield a brightness temperature of (1.56 ± 0.26 $)\times {10}^{12}$ K for the component upstream of the core, "Up," while for the unresolved core component, "Core," we obtain a lower limit of ${T}_{{\rm{b}}}\quad \gt$ 2.0 × 10¹³ K.

Zoom In Zoom Out Reset image size

The observed brightness temperature in jets is mostly affected by the transverse dimension of the flow and it may differ systematically from estimates obtained on the basis of representing the jet structure with Gaussian components (Lobanov 2015). In this case, constraints on the jet brightness temperature can also be found directly from visibility amplitudes and their errors, providing the minimum brightness temperature, ${T}_{{\rm{b}},\mathrm{min}}$, and an estimate of the formal maximum brightness temperature, ${T}_{{\rm{b}},\mathrm{max}}$ that can be obtained under the condition that the structural detail sampled by the given visibility is resolved. For further details, we refer the reader to Lobanov (2015).

These two estimates are compared in Figure 6 with the brightness temperatures calculated from the Gaussian components described in Table 2. One can see that the observed brightness temperature of the most compact structures in BL Lac, constrained by baselines longer than $5.3\;{\rm{G}}\lambda$, must indeed exceed $2\times {10}^{13}$ K and can reach as high as $\sim 3\times {10}^{14}$ K. As follows from Figure 1, these visibilities correspond to the structural scales of 30–40 μas oriented along position angles of 25°–30°. These values are indeed close to the width of the inner jet and the normal to its direction.

The observed, ${T}_{{\rm{b}},\mathrm{obs}}$, and intrinsic, ${T}_{{\rm{b}},\mathrm{int}}$, brightness temperatures are related by

$$\begin{eqnarray*}&&{T}_{{\rm{b}},\mathrm{obs}}=\delta {(1+z)}^{-1}{T}_{{\rm{b}},\mathrm{int}}\end{eqnarray*}$$

where $\delta ={(1-{\beta }^{2})}^{1/2}{(1-\beta \mathrm{cos}\phi )}^{-1}$ is the Doppler factor, β is the jet bulk velocity in units of the speed of light, ϕ is the jet viewing angle, and z is the redshift of the source. Variability arguments (Jorstad et al. 2005; Hovatta et al. 2009) and kinematical analyses (Cohen et al. 2015) yield a remarkably consistent value of $\delta =7.2$, from which we estimate a lower limit of the intrinsic brightness temperature in the core component of our RadioAstron observations of ${T}_{{\rm{b}},\mathrm{int}}\ \gt 2.9\times {10}^{12}$ K.

It is commonly considered that inverse Compton losses limit the intrinsic brightness temperature for incoherent synchrotron sources, such as AGNs, to about 10¹² K (Kellermann & Pauliny-Toth 1969). In the case of a strong flare, the "Compton catastrophe" is calculated to take about one day to drive the brightness temperature below 10¹² K. Moreover, Readhead (1994) has argued that, for sources near the equipartition of energy between the magnetic field and radiating particles, a more accurate upper value for the intrinsic brightness temperature is about 10¹¹ K (see also Lähteenmäki et al. 1999; Hovatta et al. 2009), which is often called the equipartition brightness temperature.

Our estimated lower limit for the intrinsic brightness temperature of the core in the RadioAstron image of ${T}_{{\rm{b}},\mathrm{int}}\gt 2.9\times {10}^{12}$ K is, therefore, more than an order of magnitude larger than the equipartition brightness temperature limit established by Readhead (1994), and at least several times larger than the limit established by inverse Compton cooling. This suggests that the jet in BL Lac is not in equipartition, as may be expected in the case in which the source is flaring during the ejection of component M1 detected at 43 GHz, and raises the possibility that we are significantly underestimating its Doppler factor.

We also note that if our estimate of the maximum brightness temperature is closer to actual values, it would imply ${T}_{{\rm{b}},\mathrm{int}}\sim 5\times {10}^{13}$ K. This is difficult to reconcile with incoherent synchrotron emission models from relativistic electrons, requiring alternative models such as emission from relativistic protons (Jukes 1967), or coherent emission (Benford & Tzach 2000)—see also Kellermann (2002) and references therein.

The RM and RM-corrected EVPAs (${\chi }_{0}$) in the core area (delimited by the blue ellipse in Figure 9) exhibit a clear point symmetry around its centroid. To better analyze this structure, Figure 10 displays the probability distribution function of the two-dimensional histogram for the RM and ${\chi }_{0}$ as a function of the position angle of the pixels with respect to the centroid of the core, measured counterclockwise from north.

Zoom In Zoom Out Reset image size

By inspecting Figures 9 and 10(a), we note a gradient in RM with position angle from the centroid of the core, with positive RM values in the area upstream of the centroid, and negative downstream, in the direction of the jet. The largest RM values, of the order of 3000 rad m⁻², are found in the area northeast of the centroid (with a position angle of approximately−45°); smaller values are obtained as the position angle increases, reaching values of ∼1000 rad m⁻² northwest of the centroid. Downstream of the centroid the RM continues this trend, reaching values of approximately −900 rad m⁻² in the direction of the jet.

Similarly, Figures 9 and 10(b) show a progressive rotation in ${\chi }_{0}$ with position angle from the centroid of the core. Counterclockwise from east, ${\chi }_{0}$ rotate continuously from approximately−25° to $\sim 60^\circ$. On top of this, the two-dimensional histogram shows a concentration of ${\chi }_{0}$ values between approximately 30° and 50° in the area upstream of the centroid, while downstream these concentrate in values at approximately −25° and 0°.

A similar dependence of the RM and ${\chi }_{0}$ with polar angle was found by Zamaninasab et al. (2013) in 3C 454.3, interpreted by these authors as the result of a helical magnetic field. Indeed, gradients in RM across the jet width are expected to arise in the case of helical magnetic fields (Laing 1981), as previously reported in a number of sources (e.g., Asada et al 2002; O'Sullivan & Gabuzda 2009; Hovatta et al. 2012).

Relativistic magnetohydrodynamic (RMHD) simulations have been used to study the RM and polarization distribution based on self-consistent models for jet formation and propagation in the presence of large-scale helical magnetic fields (Broderick & McKinney 2010; Porth et al. 2011). These simulations reproduce the expected gradients in RM across the jet width due to the toroidal component of the helical magnetic field and also provide detailed insights regarding the polarization structure throughout the jet, which depends strongly on the helical magnetic field pitch angle, jet viewing angle, Lorentz factor, and opacity.

As discussed in Zamaninasab et al. (2013), Broderick & McKinney (2010), and Porth et al. (2011), a large-scale helical magnetic field would lead to similar point symmetric structures around the centroid of the core of both RM and ${\chi }_{0}$ as found in Figures 9 and 10. This suggests that the core region in BL Lac is threaded by a large-scale helical magnetic field. A more detailed comparison between our observations and specific RMHD simulations using the estimated parameters for BL Lac is underway and will be published elsewhere.

The location of the stationary feature ∼0.26 mas from the core (see Section 3.1) marks a clear transition in the RM and polarization vectors between the core area and the remainder of the jet in BL Lac. Figure 9 shows a localized region of enhanced RM, reaching values of the order of $-2200\pm 300\;{\rm{rad}}$ m⁻², and a Faraday rotation corrected EVPA of $-40^\circ \pm 8^\circ$. ${\chi }_{0}$, therefore, becomes perpendicular to the local jet direction, suggesting a dominant component of the magnetic field that is aligned with the jet. Downstream of this location, ${\chi }_{0}$ rotates so that the magnetic field remains predominantly aligned with the local jet direction up to a distance of ∼1 mas from the core.

Further downstream polarization is detected again in a region at a distance of ∼3 mas from the core, corresponding to the location of components K6 and K7 (see Figure 5). This area has a mean RM and ${\chi }_{0}$ of −320 ± 140 rad m⁻² and 14° ± 6°, respectively. The polarization vectors are approximately parallel to the jet in this area, thereby suggesting that the magnetic field is predominantly perpendicular to the jet. This would be in agreement with a helical magnetic field in which the toroidal component dominates over the poloidal one, although other scenarios, like a plane perpendicular shock, cannot be ruled out.

Considering that the polarization properties in the core area and components K6 and K7 suggest that the jet in BL Lac is threaded by a helical magnetic field, the different polarization structure associated with the stationary feature at ∼0.26 mas from the core suggests that it may correspond to a rather particular jet feature. This would be in agreement with claims by Cohen et al. (2015), and references therein, in which these authors conclude that this stationary feature corresponds to a recollimation shock, downstream of which new superluminal components appear due to the propagation of Alfvén waves triggered by the swing in its position, in a similar way as exciting a wave on a whip by shaking the handle.

Previous multiwavelength observations of BL Lac suggest that the core may also correspond to a recollimation shock (Marscher et al. 2008). In that case, we can hypothesize that the upstream component found in our RadioAstron observation (see Figure 5 and Table 2) is located at the jet apex, so that the distance at which the first recollimation shock takes place, associated with the core, would be ∼40 μas. If so, the jet in BL Lac would contain a set of three recollimation shocks, spaced at progressively larger distances of approximately 40, 100, and 250 μas.

To investigate such a pattern of recollimation shocks as possibly observed in BL Lac, we have performed two-dimensional special RMHD simulations in cylindrical geometry using the RAISHIN code (Mizuno et al. 2006, 2011). The initial set-up follows Gómez et al. (1997) and Mizuno et al. (2015), in which a preexisting over-pressured flow is established across the simulation domain. We choose a rest-mass density ratio between the jet and ambient medium of $\eta ={\rho }_{{\rm{j}}}/{\rho }_{{\rm{a}}}=5\times {10}^{-3}$ with ${\rho }_{{\rm{a}}}=1.0{\rho }_{0}$, initial jet Lorentz factor ${\gamma }_{{\rm{j}}}=3$, and local Mach number ${M}_{{\rm{s}}}=1.69$. The gas pressure in the ambient medium decreases with axial distance from the jet base following ${p}_{{\rm{g}},{\rm{a}}}{(z)={p}_{{\rm{g}},0}/{[1+(z/{z}_{{\rm{c}}})}^{n}]}^{m/n}$, where ${z}_{{\rm{c}}}=60{R}_{{\rm{j}}}$ is the gas pressure scale height in the axial direction, ${R}_{{\rm{j}}}$ is the jet radius, n = 1.5, m = 2.3, and ${p}_{{\rm{g}},0}$ is in units of ${\rho }_{0}\;{c}^{2}$ (e.g., Gómez et al. 1995, 1997; Mimica et al. 2009). We assume that the jet is initially uniformly over-pressured with ${p}_{{\rm{g}},{\rm{j}}}=1.5{p}_{{\rm{g}},{\rm{a}}}$. We consider a force-free helical magnetic field with a weak magnetization B₀ = 0.05 (in units of $\sqrt{4\pi {\rho }_{0}\;{c}^{2}}$) and constant magnetic pitch ${P}_{0}={{RB}}_{z}/{R}_{{\rm{j}}}{B}_{\phi }=1/2$, where B_z and ${B}_{\phi }$ are the poloidal and toroidal magnetic field components, so that smaller P₀ refers to increased magnetic helicity.

Overexpansions and overcontractions caused by inertial overshooting past equilibrium of the jet lead to the formation of a pattern of standing oblique recollimation shocks and rarefactions, whose strength and spacing are governed by the external pressure gradient (see Figure 11). The jet is accelerated and conically expanded slightly by the propegation of rarefaction (e.g. Aloy & Rezzolla 2006; Mizuno et al. 2008) and the gas gradient force, reaching a maximum jet Lorentz factor of ∼8. The jet radius expands from the initial ${R}_{{\rm{j}}}=1$ at the jet base to ${R}_{{\rm{j}}}=1.5$ at $z=60{R}_{{\rm{j}}}$, yielding an opening angle of ${\theta }_{{\rm{j}}}=\mathrm{arctan}(0.5{R}_{{\rm{j}}}/60{R}_{{\rm{j}}})\sim 0\buildrel{\circ}\over{.} 48$. The recollimation shocks are located at progressively larger axial distances of $\sim 7{R}_{{\rm{j}}}$, $\sim 20{R}_{{\rm{j}}}$, and $\sim 42{R}_{{\rm{j}}}$. The relative distances between the second and third shocks with respect to the first one are ∼2.9 and ∼6, which roughly match those observed in BL Lac.

Zoom In Zoom Out Reset image size

The output from the RMHD simulation is used as input to compute the radio continuum synchrotron emission map shown in Figure 12. Following Gómez et al. (1995), the internal energy is distributed among the non-thermal electrons using a power law $N(E){dE}={N}_{0}{E}^{-p}{dE}$, with p = 2.4. The emission is computed for a viewing angle of 10^◦ and at an optically thin observing frequency, integrating the synchrotron transfer equations along the line of sight (e.g., Gómez et al. 1997). The pattern of recollimation shocks seen in Figure 11 leads to a set of knots in the total and linearly polarized intensity. The EVPAs in the knots are perpendicular to the jet direction, in agreement with the observations of component K2 in BL Lac. We note, however, that our RMHD simulations consider an already formed and collimated jet, and therefore do not provide an accurate account of the jet formation region, close to where the other two innermost recollimation shocks are located.

Zoom In Zoom Out Reset image size

Finally, although our simulations are in general agreement with the pattern of recollimation shocks observed in BL Lac, further numerical simulations, in progress, are required to better constrain the jet parameters of the model, such as the Lorentz factor, Mach number, gradient in external pressure, magnetic field helicity/strength, and viewing angle.