High-Sensitivity 86GHz (3.5mm) VLBI Observations of M87: Deep Imaging of the Jet Base at a 10 Schwarzschild-Radius Resolution
Kazuhiro Hada, Motoki Kino, Akihiro Doi, Hiroshi Nagai, Mareki Honma, Kazunori Akiyama, Fumie Tazaki, Rocco Lico, Marcello Giroletti, Gabriele Giovannini, Monica Orienti, Yoshiaki Hagiwara
aa r X i v : . [ a s t r o - ph . H E ] D ec Draft Version
High-Sensitivity 86 GHz (3.5 mm) VLBI Observations of M87:Deep Imaging of the Jet Base at a 10 Schwarzschild-Radius Resolution
Kazuhiro Hada , , Motoki Kino , Akihiro Doi , Hiroshi Nagai , Mareki Honma , , KazunoriAkiyama , , Fumie Tazaki , Rocco Lico , , Marcello Giroletti , Gabriele Giovannini , , MonicaOrienti , and Yoshiaki Hagiwara , ABSTRACT
We report on results from new high-sensitivity, high-resolution 86 GHz (3.5 mil-limeter) observations of the jet base in the nearby radio galaxy M87, obtained by theVery Long Baseline Array in conjunction with the Green Bank Telescope. The re-sulting image has a dynamic range exceeding 1500 to 1, the highest ever achieved forthis jet at this frequency, resolving and imaging a detailed jet formation/collimationstructure down to ∼
10 Schwarzschild radii ( R s ). The obtained 86 GHz image clearlyconfirms some important jet features known at lower frequencies, i.e., a wide-openingangle jet base, a limb-brightened intensity profile, a parabola-shape collimation profileand a counter jet. The limb-brightened structure is already well developed at < . < R s , projected) from the core, where the corresponding apparent opening anglebecomes as wide as ∼ ◦ . The subsequent jet collimation near the black hole evolvesin a complicated manner; there is a “constricted” structure at tens R s from the core,where the jet cross section is locally shrinking. We suggest that an external pressuresupport from the inner part of radiatively-inefficient accretion flow may be dynamicallyimportant in shaping/confining the footprint of the magnetized jet. We also presentthe first VLBI 86 GHz polarimetric experiment for this source, where a highly polarized( ∼ Subject headings: galaxies: active — galaxies: individual (M87) — galaxies: jets —radio continuum: galaxies
1. Introduction
Accreting supermassive black holes at the center of active galaxies produce powerful relativisticjets that are observed as a collimated beam of plasma, often propagating beyond the host galaxies.Understanding the formation, collimation and propagation of relativistic jets is a longstanding con-cern in astrophysics (Blandford & Rees 1974; Blandford & Znajek 1977; Blandford & Payne 1982; 2 –Begelman et al. 1984), and the recent theoretical progress based on general relativistic magneto-hydrodynamical simulations has begun to elucidate the roles of the central black hole, surroundingaccretion flow and magnetic fields threading them, as well as the mutual interactions among thesecomponents, in generating and collimating a jet (e.g., McKinney 2006; Komissarov et al. 2007;McKinney & Blandford 2009; Tchekhovskoy et al. 2011; McKinney et al. 2012). To test the im-plications from such theories and then to better understand the jet formation, it is necessary topresent a detailed observation that can image the relevant scales.M87 is the first extragalactic jet discovered by Curtis nearly 100 years ago (Curtis 1918).This jet is exceptionally close to us ( D = 16 . γ -rays. These observational ad-vantages have allowed a broad range of studies associated with relativistic-jet physics, includingthe large-scale jet morphology (e.g., Owen et al. 1989), the nature of the optical jet emission andshocks (e.g., Biretta et al. 1999; Perlman et al. 2001) and the origin of the high-energy X-ray to γ -ray emission (e.g., Harris et al. 2006; Abramowski et al. 2012; Hada et al. 2014). Moreover, op-tical measurements of the nuclear stellar dynamics suggest the presence of a huge central blackhole of M BH = (6–6 . × M ⊙ (Gebhardt & Thomas 2009; Gebhardt et al. 2011), although gas-dynamical measurements derive a factor of two smaller M BH (Ford et al. 1994; Harms et al. 1994;Macchetto et al. 1997; Walsh et al. 2013). The combination of the proximity and the large blackhole yields a linear resolution down to 1 milliarcsecond (mas) = 0.08 pc = 140 Schwarzschild radii( R s ) (for D = 16 . M BH = 6 × M ⊙ ), which is typically 10 to 100 times finer thanthose accessible in distant quasars or blazars. Therefore, M87 offers a privileged opportunity forprobing the launch/formation scales of a relativistic jet with high-resolution Very-Long-Baseline-Interferometer (VLBI) observations.With the recent advent of the global short-millimeter VLBI project so-called the Event Horizon Mizusawa VLBI Observatory, National Astronomical Observatory of Japan, Osawa, Mitaka, Tokyo 181-8588,Japan; [email protected] INAF Istituto di Radioastronomia, via Gobetti 101, I-40129 Bologna, Italy Korea Astronomy and Space Science Institute (KASI), 776 Daedeokdae-ro, Yuseong-gu, Daejeon 305-348, Re-public of Korea Institute of Space and Astronautical Science, Japan Aerospace Exploration Agency, 3-1-1 Yoshinodai, Chuo,Sagamihara 252-5210, Japan National Astronomical Observatory of Japan, Osawa, Mitaka, Tokyo 181-8588, Japan Department of Astronomical Science, The Graduate University for Advanced Studies (SOKENDAI), 2-21-1 Os-awa, Mitaka, Tokyo 181-8588, Japan Department of Astronomy, Graduate School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo113-0033, Japan Dipartimento di Fisica e Astronomia, Universit`a di Bologna, via Ranzani 1, I-40127 Bologna, Italy Toyo University, 5-28-20 Hakusan, Bunkyo-ku, Tokyo 112-8606, Japan µ as, corresponding to 5 . R s (if M BH = 6 . × M ⊙ is adopted). More recently, interferometric closure phases for the corresponding structure have beenobtained at this wavelength (Akiyama et al. 2015), allowing a comparison between the observationand some horizon-scale theoretical models. Going to shorter wavelengths is also beneficial with re-spect to the synchrotron opacity, since the jet base becomes more transparent to the self-absorptioneffect (K¨onigl 1981). Indeed, the λ +0 . dependence of the M87 radio core position revealed by as-trometric measurements suggests the 1.3 mm core to be located within a few R s from the blackhole (Hada et al. 2011). Nevertheless, current VLBI experiments at such short wavelengths are stilltechnically challenging to synthesize interferometric images due to the limited number of availableantennas (thus only sparse uv -coverage) as well as the severe atmospheric disturbance (thus shortercoherence time). This prevents us from tracking the larger-scale propagation of the flow launchedfrom the central horizon-scale dimension, which is essential to fully understand the subsequentacceleration and collimation of the jet.So far, VLBI imaging studies of M87 have mostly been made at 7 mm (43 GHz), 1 cm (22–24 GHz), 2 cm (15 GHz) or longer. In these bands the M87 jet is bright enough and a sufficientnumber of VLBI stations are available for allowing an adequate uv -coverage. Previous high-qualityVLBI images of the M87 inner jet revealed some important features such as a wide-opening anglebase, a limb-brightened intensity profile and a counter jet (Junor et al. 1999; Ly et al. 2004, 2007;Kovalev et al. 2007; Hada et al. 2011), as well as a detailed movie near the jet base (Walker et al.2008). More recently, the inner jet was found to sustain a parabola-shape collimation profileover a wide range of distance from ∼
100 to ∼ R s from the nucleus (Asada & Nakamura 2012;Hada et al. 2013; Nakamura & Asada 2013). However in these bands, it is impossible to achieveangular resolution comparable to that of the short-millimeter VLBI unless relying on a space-VLBI satellite (Hirabayashi et al. 1998; Dodson et al. 2006; Kardashev et al. 2013). In addition,the higher optical depth at such long wavelengths precludes us from observing the close vicinity ofthe black hole. As a result, there still remains a large gap in our current understanding of this jetbetween the centimeter and the short-millimeter VLBI scales.In this context, an important “bridge” to connect this gap is observational study at 3.5 mm(86 GHz). At present, 3.5 mm is the shortest wavelength where one can reliably obtain synthesizedVLBI images, as represented by studies with the Global Millimeter VLBI Array (GMVA; e.g.,Giroletti et al. 2008; Lee et al. 2008; Agudo et al. 2007; Molina et al. 2014; Boccardi et al. 2015;Hodgson et al. 2015; Koyama et al. 2015). The angular resolution with 3.5 mm VLBI is typicallytwice better than that at 7 mm, and the transparency to a jet base is also higher. On the otherhand, compared to 1.3 mm, one can detect the extended (optically-thin) emission much furtherdown the jet due to the steep-spectral nature of the synchrotron radiation, which allows a bettermonitoring of the larger-scale jet propagation. Therefore, the use of 3.5 mm is currently an optimalchoice in terms of angular resolution, opacity and capability for imaging a jet. Nevertheless, there 4 –have been still less M87 observations in this band because 3.5 mm VLBI is generally less sensitivethan that at 7 mm due to the more rapid atmospheric fluctuations as well as the worse effectiveaperture efficiency of telescopes. Thus the brightness of the M87 jet base (typically hundredsmJy to < ∼ D = 16 . M BH = 6 . × M ⊙ for M87, correspondingto 1 mas = 0.08 pc = 140 R s . Spectral index α is defined as S ν ∝ ν + α . Also, any λ -related numbersare described in frequency unit in the rest of the paper.
2. Observations and data reduction2.1. 86 GHz data
In February 2014 we observed M87 at 86 GHz with VLBA in conjunction with GBT. FromVLBA eight out of the ten stations participated in this program since the other two (i.e., Hancockand Saint Croix) do not have the 86 GHz receiver. To increase the overall sensitivity, we madean 8-hour-long quasi-full-track observation twice on February 11 and 26. The observations weremade in dual (left/right-hand circular) polarization mode. The received signals were sampledwith 2-bit quantization and recorded at aggregate rate of 2 Gbps (a total bandwidth of 512 MHz)using the digital-downconverter-4 (DDC-4) wideband recording mode. The down-converted signalswere divided into two 128 MHz sub-bands in each polarization respectively. As an overall systemcalibrator (fringe check, bandpass, delay tracking, see below) of this program, we inserted 5-minutescans on the nearby bright source 3C 273 (10 degrees apart from M87 on the sky) every 30 minutes.3C 273 was observed also for the purpose of antenna pointing by adding another 6-minute-long scansevery 30–60 minutes. The second epoch has better weather conditions and system temperatures 5 –Table 1: VLBA observations of M87
UT Date ν Stations ∆ ν Beam size I peak I rms I peak /I rms (GHz) (MHz) (mas × mas, deg.) (cid:16) mJybeam (cid:17) (cid:16) mJybeam (cid:17) (a) (b) (c) (d) (e) (f) (g)2014 Feb 11............... 86.266 VLBA, GBT, − SC, − HN 512 0 . × . , −
12 (0 . × . , −
13) 549 0.63 8712014 Feb 26............... 86.266 VLBA, GBT, − SC, − HN 512 0 . × . , −
10 (0 . × . , −
11) 521 0.36 14472014 Feb (11+26) (h) .. 86.266 VLBA, GBT, − SC, − HN 512 0 . × . , −
10 (0 . × . , −
12) 547 0.29 18862014 Mar 8................. 43.230 VLBA, − MK, − FD 128 0 . × . ,
25 (0 . × . ,
27) 786 0.85 92423.830 VLBA, − MK, − FD 128 0 . × . ,
18 (0 . × . ,
21) 880 0.82 10732014 Mar 26............... 43.296 VLBA 256 0 . × . ,
12 (0 . × . ,
12) 795 0.50 159023.894 VLBA 256 0 . × . , . × . , −
1) 926 0.43 21532014 May 8................. 43.296 VLBA 256 0 . × . , − . × . , −
3) 735 0.56 131223.894 VLBA 256 0 . × . , − . × . , −
7) 869 0.64 1358
Notes: (a) central frequency; (b) participating stations. VLBA indicates all the ten VLBA sta-tions. GBT, SC, HN, MK and FD are the Green Bank Telescope, Saint Croix, Hancock, MaunaKea and Fort Davis, respectively. A minus sign before station name means the absence of thatstation; (c) total bandwidth; (d) synthesized beam with a naturally-weighting scheme. For refer-ence, a beam size with a uniformlly-weighting scheme is also shown in bracket; (e) peak intensity ofM87 images under naturally-weighting scheme; (f) off-source rms image noise level of M87 imagesunder naturally-weighting scheme; (g) dynamic range calculated with I peak and I rms ; (h) combinedvisibility data over the two epochs.over the array. The information on these data is summarized in Table 1.The initial data calibration was performed with the Astronomical Image Processing System(AIPS) developed at NRAO. We first corrected the visibility amplitude by applying the measuredsystem noise temperature and the elevation-gain curve of each antenna. Atmospheric opacitycorrections were made by solving for receiver temperature and zenith opacity for each antenna.We then calibrated the amplitude part of the bandpass characteristics for each station using theauto-correlation spectra of 3C 273.The calibration of the visibility phase was made following the amplitude calibration. To recover86 GHz fringes as much as possible, we performed the phase calibration in the following way. Wefirst corrected known phase variations due to parallactic angle effects. Next, the instrumental phaseand delay offsets for each antenna were derived using a scan of 3C 273, and subtracted from thewhole dataset assuming that they are constant with time (note that at 86 GHz no pulse-cal signalsare available to calibrate instrumental effects). After that, we ran a global fringe-fitting on 3C 273with its source structure model created by the procedure described below (a point source model wasassumed in the first round where the source model was not available), and derived time evolutionsof the residual delay, rate and phase for each IF separetely. We detected fringes for most of thescans at adequate SNRs. Since the derived residual delay component (that mainly comes fromunpredicted atmospheric fluctuations in the correlation stage) slowly varies on the sky, we can usethe 3C 273’s delay solutions as a good first-order approximation for those of M87’s residual delay.We thus interpolated the derived 3C 273’s delay solutions to the scans of M87. Implementing this 6 –Fig. 1.— Normalized histograms of signal-to-noise ratio for the detected fringes (by global fringefitting). The red histogram is for M87, while the blue one is for 3C 273. A solution interval of 15 secis used here.process, now we can go to fringe-fitting on M87 with a tight delay search window (as small as ±
10 nsec), which is quite effective in avoiding false signals. We performed a global fringe-fittingon M87 with a solution interval of 15 sec and an SNR threshold of 4.0. Since the phase and delayoffsets between IFs are already removed, we derived IF-averaged solutions for M87 to increase SNRby a factor of 1.4. Similarly to 3C 273, we ran two cycles of fringe-fitting iteratively; the first round(where no source model available) was performed with a point source model, and after creating acoarse image, the second round was executed with the source model. This indeed slightly increasedthe fringe detection rate. In every fringe-fitting the most sensitive station GBT was used as thereference antenna.Through this procedure, we recovered a number of fringes at sufficient SNRs. In Figure 1 weshow (normalized) histograms of SNR for the detected fringes on M87 and 3C 273. For 3C 273 themedian/mean SNRs were 38 and 38, respectively. For M87 the median/mean SNRs were 19 and20, respectively. Note that the SNR histogram for 3C 273 is relatively widely distributed. This is 7 –because the source structure of 3C 273 is highly complicated (see Appendix) and the correlated fluxdensity changes more drastically with baseline length, than in the case of M87.After the visibility data became coherent with time and frequency, the data were integratedover frequency (but the two IFs kept separated), and in time to 30 sec. These data were then usedfor creating images.We performed our imaging process in the following way. We first exported the averaged datato the Difmap software (Shepherd et al. 1994), and performed intensive flagging of the visibilitydata. At 86 GHz, antenna pointing is generally less accurate than that at lower frequencies, andthis causes an unwanted systematic decrease of visibility amplitude in some scans. These scansthen cause significant sidelobes in the initial stage of imaging process, which prevents us fromreconstructing a reliable CLEAN/deconvolution model. Thus, in addition to obvious outliers, wecarefully flagged such bad scans in an antenna-based manner.Following the initial exhaustive flagging, we then worked on iterative CLEAN/self-calibrationprocessing. We conducted this process using Difmap and AIPS in a hybrid manner. We used Difmaponly to perform CLEAN deconvolution (because CLEAN with Difmap is faster and more intuitivethan with AIPS). Then after obtaining a reasonable set of CLEAN components, the model and the(uncalibrated) visibility data were exported back to AIPS, and we performed self-calibration usingthe AIPS task CALIB. This is advantageous for better calibration because the self-calibration withCALIB can be handled more robustly and flexibly (e.g., in terms of weighting for each antennaand SNR cutoff setting) than that in Difmap. Self-calibration with CALIB was also necessaryfor proper polarization analysis (see below) since CALIB can solve the complex gain terms forLL/RR polarization separately, while Difmap cannot do that. Then the self-calibrated visibilitydata were again exported back to Difmap to create an improved CLEAN model. We repeatedthis CLEAN/self-calibration round until the reconstructed CLEAN model no longer improved sig-nificantly. In the first several CLEAN/self-calibration cycles, the phase-only self-calibration wasdone, and after the phase part was well corrected, the visibility amplitude was also self-calibratediteratively with solution intervals starting from several hours down to a few minutes.In Figure 2 we show a resulting uv -distance plot of the calibrated visibility amplitude of M87.Thanks to the excellent quality of the dataset, one can see that the visibilities are homogeneouslysampled over the entire uv -distance and thus the overall trend is well-defined. One can see thatthe whole visibility set primarily consists of two components. One is an extended componentwhich significantly contributes at baselines shorter than . λ , while the other is a compactcomponent that dominates beyond 0.2 G λ to the longest 2 G λ baseline. For the compact component,the correlated flux densities are monotonically decreasing with increasing uv -distance, and theamplitude is .
100 mJy at the longest baseline. Despite such low flux densities, the signals wererobustly detected since the longest baseline consists of GBT-MK pair. 8 –Fig. 2.— Visibility amplitude versus uv -distance plot of M87 for the VLBA+GBT 86 GHz obser-vations. The visibility data shown are post-self-calibrated ones. Although our 86 GHz observations were not ideally designed for accurate polarimetric study,we can still attempt a polarization analysis by taking advantage of the scans of 3C 273. Followingthe parallactic angle corrections described above, the cross-hand R-L phase and delay offsets werecalibrated using a scan of 3C 273. The feed polarization leakage for each antenna was corrected byusing the LPCAL method in AIPS (Lepp¨anen et al. 1995) with a total intensity model of 3C 273.The derived leakage values were different with station, polarization, IF and epoch, but we obtaineda value of 9 ± ± D -term”) scheme is validated. GBT was one of the best stations wherethe leakages were as small as ∼ ∼ ∼ ∼ σ in the core 9 –region (see Appendix). This is consistent with the above calculation.Regarding the electric vector polarization angle (EVPA), we cannot derive its absolute valuefrom our data alone, since we did not perform any additional EVPA-calibration observation. How-ever, we found a 43 GHz VLBA observation of 3C 273 that was carried out close in time with ours(on 2014 Feb 25) in the Boston University blazar monitoring program. We used this 43 GHz polar-ization image as a reference of our EVPA correction (see Figure 11 in Appendix). Here we assumethat the EVPA of the outermost polarized component (P3 at 1.5 mas from the core) is stable withboth time and frequency, and performed a nominal correction of EVPA by matching the observed86 GHz EVPA of P3 to that of the 43 GHz one. Note that this assumption may not be correct,since the previous concurrent 86/43 GHz VLBA polarimetric study of this source suggests a largerotation measure (RM) of ∼ × rad m − for inner jet components at ∼ ∼ × rad m − around P3 (see Appendix for more details).If this is the case, one would expect another systematic rotation of the 86 GHz EVPA about ∼ ◦ with respect to that at 43 GHz (i.e., ∆ χ RM ∼ ◦ ). As for the Boston 43 GHz polarization image,we adopt its EVPA uncertainty to be ∆ χ ∼ ◦ based on Jorstad et al. (2005). Therefore,we estimate that a potential total uncertainty of absolute EVPA in our 3C 273’s 86 GHz images is∆ χ ∼ ∆ χ + ∆ χ RM ∼ ± ◦ . Regarding M87, its EVPA uncertainty in 86 GHz imageswould be somewhat larger than this value, since the lower SNR of polarization signals (SNR ∼ χ therm (radian) ∼ σ p / P where σ p and P are rms noise level and polarized intensity inpolarization map (e.g., Roberts et al. 1994). With SNR= P/σ p ∼ χ therm , M87 ∼ ◦ .Assuming that ∆ χ and ∆ χ therm , M87 are statistically independent, we estimate a total errorbudget for M87 to be ∆ χ M87 ∼ ± ◦ . As supplementary datasets, we additionally made VLBA-alone observations of M87 at 24 and43 GHz close in time with the 86 GHz sessions. The observations were carried out on March 8, 26and May 8 2014, where both 24 and 43 GHz were quasi-simultaneously used by alternating eachreceiver quickly. On March 26 and May 8, all the VLBA stations were present, while on March 8the antennas at Mauna Kea and Fort Davis were absent. We received only RR polarization signalswith a total bandwidth of 128 MHz (on March 8) or 256 MHz (on March 26 and May 8). Amongthese sessions, the data on March 26 were the best in overall quality, while the data on March 8were relatively poor. The initial data calibration (apriori amplitude correction, fringe fitting andbandpass) was made in AIPS, and the subsequent image reconstruction was performed in Difmapbased on the usual CLEAN/self-calibration procedure. The basic information of these data is alsotabulated in Table 1. 10 –
3. Results3.1. New 86 GHz images
In Figure 3 we show a representative 86 GHz image of the M87 jet obtained by our VLBA+GBTobservations. For a better visualization, the image is produced by combining the visibility dataover the two epochs, and restored with a convolving beam of 0.25 mas × ◦ . A contour image with a naturally-weighting scheme is also displayed in the toppanel of Figure 4.Thanks to the significant improvement of sensitivity, a detailed jet structure was clearly imageddown to the weaker emission regions. The resulting image rms noise of the combined image was ∼ − . In this period the extended jet was substantially bright down to ∼ ∼ σ , and another ∼ σ level. The peak surface brightnessof the image was 500 mJy beam − at this resolution, corresponding to an image dynamic rangegreater than 1500 to 1 (the detailed value slightly varies as a function of the weighting schemeand convolving beam). This is the highest image dynamic range obtained so far at 86 GHz for thisjet, and is quite comparable to typical dynamic ranges in VLBA images at 43 GHz (e.g., Ly et al.2007). We describe a comparison of our 86 and 43 GHz images in the next subsection.Regarding the individual epochs, the second epoch was better in overall image quality thanthat at the first epoch. As listed in Table 1, the synthesized beam for the first epoch is moreelongated in the north-south direction than that for the second epoch. This is mainly becausethe Brewster station, which constitutes the longest baselines in the north-south direction, had thehigher system noise temperatures during the first session. Thus the less-weighted uv dataset onthis station creates a slightly larger fringe spacing along the north-south direction.Consistent with known lower-frequency images, most of the radio emission at 86 GHz is con-centrated on the compact radio core at the jet base. To quantify the structure of the core region, weperformed a single elliptical Gaussian modelfitting to the calibrated visibility data with the Difmaptask modelfit . The derived model parameters are summarized in Table 2. As a check we per-formed the same fitting to the three different datasets, i.e., the Feb/11 data, the Feb/26 data andthe combined one, but virtually the same result was obtained. Additionally, we performed anotherelliptical Gaussian model fitting on the image plane using the AIPS task JMFIT and examined thedeconvolved result, but this was also essentially the same within errors. The derived geometry ofthe core is close to a circular shape with a diameter of ∼ µ as, which is just in between the sizesobtained at 230 GHz (40 µ as; Doeleman et al. 2012) and 43 GHz (110–130 µ as; Hada et al. 2013).A size of 80 µ as is consistent with that obtained in our previous study based on an archival VLBA86 GHz dataset (Hada et al. 2013), but the result presented here is much more reliable. Adoptingthe parameters derived with the combined data, a brightness temperature of the 86 GHz core isestimated to be T B = 1 . × K. Although the observed epochs are different, this value is quite 11 –
M87
Fig. 3.— VLBA+GBT 86 GHz false-color total intensity image of the M87 jet. The image isproduced by combining the visibility data over the two epochs on 2014 February 11 and 26. Therestoring beam (0 . × .
08 mas in PA 0 ◦ ) is shown in the bottom-right corner of the image. Thepeak intensity is 500 mJy beam − and the off-source rms noise level is 0.28 mJy beam − , where theresulting dynamic range is greater than 1500 to 1. (A color version of this figure is available in theonline journal.) 12 –Fig. 4.— VLBA+GBT 86 GHz total intensity contour images of the M87 jet. The upper panelindicates a naturally-weighted image with a synthesized beam of 0 . × .
11 mas in PA − ◦ . Thelower panel shows a better-resolved image restored with a circular Gaussian beam whose FWHM isequivalent to that of the minor axis of the beam used in the top panel. In each image the restoredbeam is shown at the top-left corner of each panel. Contours in both images start from −
1, 1, 2 / ,2, 2 / , 4... times 0.86 mJy beam − . 13 –Table 2: Modelfit parameters for 86 GHz coreData θ maj θ min PA S core ( µ as) ( µ as) (deg.) (mJy)(a) (b) (c) (d)Feb 11 81 ± ±
10 56 ± ± ± ± ±
20 652 ± (e) ± ± ± ± modelfit in Difmap and JMFIT in AIPS. For (d), we adopt 10% uncertainty basedon the absolute typical amplitude calibration accuracy.similar to T B reported for the 230 GHz core ((1.2–1.4) × K; Akiyama et al. 2015). No significantvaliability was found in T B of the core between our two 86 GHz sessions.Downstream of the core, we clearly identified a limb-brightened jet profile as seen in Fig-ure 3 and Figure 4. While such a limb-brightened structure in M87 is repeatedly confirmed inprevious VLBA 43 GHz/15 GHz images (e.g., Junor et al. 1999; Ly et al. 2007; Kovalev et al. 2007;Hada et al. 2011, 2013), it was not so clear in previous VLBA 86 GHz images (Rioja & Dodson 2011;Nakamura & Asada 2013), although one of the early GMVA images presented by Krichbaum et al.(2006) suggested a hint of limb-brightening.For a better description of the near-core structure, in the lower panel of Figure 4 we displaythe same image as in the upper panel but restored with a circular Gaussian beam whose diameter isequal to the minor axis of the synthesized beam in the upper panel. The jet launching morphologywithin ∼ < . < R s , projected) from the core with a large opening angle. We investigate the transverse jetstructure in more detail in Section 3.3.At the eastern side of the core, we detected weak but significant emission (at a level of 6 σ inthe combined image) at ∼ −
1, 1, 2 / , 2, 2 / , 4... times 0.9 mJy beam − .In Figure 5 we also show a tapered image restored with a 0.4-mas-diameter beam to emphasizethe larger-scale emission more noticeably. Consistent with known lower-frequency images, the86 GHz jet is extending to the northwest direction on the large scale. In this period the southernlimb is brighter within ∼ ∼ ◦ , matched with that in lower-frequency images. Although our 86 and 43/24 GHz observations were not simultaneous but 10 days to 10 weeksapart, it is still useful to compare them to examine any structural consistency or variations. In fact,their comparable image dynamic ranges at a level greater than 1000 to 1 allow a reliable imagecomparison between 86 GHz and the lower frequencies for the first time.
In Figure 6 we show a sequence of our 86 and 43 GHz images, where all the images are re-stored with the same convolving beam of 0 . × .
11 mas in PA 0 ◦ , which is approximately anintermediate resolution between 86 and 43 GHz. The overall jet shape and characteristic structureis in good agreement with each other. It is known that the M87 jet is relatively smooth and lessknotty, but here we do see some noticeable features or patterns in the jet when imaged at the high 15 – N1 N2 N3 N4S1 S5S4S3S2CJ S1 S5S4S3S2CJ N1 N2 N3 N4CJ S1 S5S4S3S2N1 N2 N3 N4N1 N2 N3 N4S1 S5S4S3S2CJ
Fig. 6.— Multi-epoch images of the M87 jet. From the top, we show images observed on 2014February 11 at 86 GHz, 2014 February 26 at 86 GHz, 2014 March 26 at 43 GHz and 2014 May 8 at43 GHz, respectively. All the images are convolved with a common beam of 0 . × .
11 mas in PA 0 ◦ (shown at the bottm-right corner of each panel). Contours on each image are −
1, 1, 2 / , 2, 2 / , 4...times 1.0 mJy beam − (upper two panels) and 2.0 mJy beam − (lower two panels), respectively. Thecomponents with labels are the ones identified over the different epochs/frequencies consistently. 16 –resolution. As shown in Figure 6 we identified several components in both 86 and 43 GHz imagesconsistently, including a counter jet component (which is more prominent at 43 GHz). These fea-tures are marked as CJ (counter jet), N1–N4 (in the northern limb) and S1–S5 (in the southernlimb), respectively. Since these components are typically & & ∼ While these components were consistently identified over the observed period, they were indeedgradually moving in the jet. In Figure 7 we show the observed positions of these components withrespect to the core over the monitoring period. For each component the position was measuredby fitting an elliptical Gaussian model with the AIPS task JMFIT, and its position uncertaintywas estimated by the fitted size divided by the peak-to-noise ratio (Fomalont 1999). The measuredproper motion results are summarized in Table 3. Note that there may be an absolute positionoffset between the 86/43 GHz cores due to core-shift. However, such a shift is expected to beonly ∼ µ as between 86/43 GHz (assuming the ν − . dependence determined at low frequencies;Hada et al. 2011), which is substantially small compared to the proper motions observed in thepresent study.In the main jet a mean speed of the observed components is β app = 0 .
32, where most of themare moving at a similar speed in the range of 0.3–0.5. On the other hand, the counter jet CJ ismoving in the opposite (to the northeast) direction at a slightly slower speed of ∼ . c . Thesevalues observed both in the jet and counter jet are similar to those suggested by Ly et al. (2007),although at that time the measurement was based on only one pair of VLBA 43 GHz observationsseparated by more than 6 months. Note that S3 and S4 appear to be quite slow or quasi-stationarycompared to the rest of the features in the main jet. For these features we cannot completely ruleout the possibility of component misidentification among the different epochs. However, lookingat the overall evolution of the jet morphology within ∼ & c ) motionsas reported in Walker et al. (2008) and Acciari et al. (2009). 17 –Fig. 7.— Observed sky positions of each component relative to the core. The components in themain jet are all moving toward the western direction, while the counter jet component CJ is movingto the opposite direction (i.e., to the northeast).Table 3: Component motionFeatures r PA µ β app (mas) (deg.) (mas yr − ) ( v app /c )(a) (b) (c) (d) (e)CJ 0 . ± .
05 106 ± − . ± . − . ± . . ± .
01 322 ± . ± .
08 0 . ± . . ± .
04 317 ± . ± .
15 0 . ± . . ± .
16 316 ± . ± .
16 0 . ± . . ± .
08 304 ± . ± .
46 0 . ± . . ± .
02 277 ± . ± .
02 0 . ± . . ± .
05 275 ± . ± .
01 0 . ± . . ± .
05 272 ± . ± .
18 0 . ± . . ± .
13 274 ± . ± .
50 0 . ± . . ± .
11 274 ± . ± .
61 0 . ± . . ± .
11 306 ± . ± .
65 0 . ± . Comparing the closest pair of the data on Feb/26 (at 86 GHz) and Mar/8 (at 43 GHz), theobserved peak brightness in each image (when measured with a common 0 . × . − (at 86 GHz on Feb/26) and 546 mJy beam − (at 43 GHz onMar/8), respectively. This results in a non-simultaneous (but still close-in-time) 43/86 GHz spectralindex of the core (at this resolution) being flat ( α c , / ∼− . α j , / = − . ± . α j , / = − . ± .
3. These two values are consistent within the errors, although a possiblespectral steepening might exist at >
43 GHz due to the higher cooling efficiency.As for the counter jet, the observed spectrum between 43/86 GHz seems to be steeper thanthat of the main jet ( α cj , / ∼− . Given that the brightness profile of the main jet is very rich while the counter jet is weak andless characterized, there is a large uncertainty in determining the exact value of the jet-to-counter-jet brightness ratio (BR).Here we consider the following possibilities to estimate BR: (1) CJ is the counter part of S1:(2) CJ is the counter part of S2. The choice of (1) is because CJ and S1 are symmetrically locatedwith respect to the core, while the case (2) is also possible if we additionally consider the relativistic“arm-length ratio” that reflects a factor of ∼ ∼ Resolving a jet in the direction transverse to the jet axis is important for understanding theopening angle, collimation efficiency and possible velocity gradient across the jet. To perform this,one needs a high quality image at a sufficient angular resolution across the jet. The new 86 GHzimage presented here allows us to analyze a detailed transverse structure of the jet launch regionnear the black hole.In Figure 8 we show a close-up view of the innermost region of the M87 jet. For a betterdescription across the jet, the image is restored with a circular Gaussian beam whose FWHM isequivalent to that of the minor axis of the synthesized beam. This is the same contour image asthat in the lower panel of Figure 4, but here the image is rotated on the sky by − ◦ in order toalign the jet central axis to the horizontal axis.A strongly limb-brightened profile is evident beyond 0.5 mas from the core, as consistently seenin the 43 GHz images. Within 0.5 mas from the core, the limb-brightening is continuously visible.It is notable that the limb-brightened jet is already formed at 0.15 mas distance from the core.Moreover, one more remarkable feature near the core is that the limb-brightened jet is evolvingin a highly complicated manner; in particular, there is a “constricted” structure at ∼ W ( r ), whilethe lower panel plots the corresponding (apparent) opening angle profile φ ( r ) ≡ W ( r ) / r )assuming the size of the jet origin being infinitesimally small. As for W ( r ) here we defined it by 20 – Fig. 8.— Close-up view of the innermost region of the M87 jet observed with VLBA+GBT at86 GHz. The image is restored with a 0.11 mas circular Gaussian beam, which is equivalent to thatof the minor axis of the synthesized beam. The image is rotated on the sky by − ◦ in order toalign the jet central axis to the horizontal direction. For reference, on the origin of the coordinatewe superpose a green-colored circle having a diameter of 40 µ as, which is equivalent to the size ofthe 230 GHz core reported by Doeleman et al. (2012). The arrows point the place where the jetshape is locally shrinking. 21 –Fig. 9.— Radial profiles of the transverse structure of the M87 jet. The upper panel indicates ajet width profile W pp ( r ), which is obtained by measuring the peak-to-peak separation of the twolimbs at each r from the core. The grey dashed line represents ∝ r . dependence, which wasdetermined in our previous extensive jet width measurement between r ∼ r ∼
400 masfrom the core (Hada et al. 2013). The dashed line plotted here is not a fitted one to the data butis arbitrarily placed just for reference. The lower panel shows the corresponding apparent openingangle profile φ pp ( r ) ≡ W pp ( r ) / r ). The data points with grey color are previous resultsby Junor et al. (1999) and Biretta et al. (2002). In both panels the data with magenta color, bluecolor and green color refer to 86, 43 and 24 GHz, respectively. The arrow in each panel correspondsto the location of the arrow in Figure 8. 22 –the peak-to-peak separation of the two limbs perpendicular to the jet axis (thus for clarity, weredenote the present measurements as W pp ( r ) and φ pp ( r ) ≡ W pp ( r ) / r )) rather than theusual outer-edge-to-outer-edge separation that would be more appropriate for expressing an “entirejet width”. The reason why we use W pp here is because measurements of peak-to-peak separationare less affected by the applied convolving beam, while the measurements based on the outer edgesare more sensitive to the beam size. This is particularly relevant to the near-core region wherethe jet cross section is comparable to or even smaller than that of the synthesized beam. Ourpresent purpose is not to determine an absolute width of the jet but rather to measure a radialdependence of the jet evolution as close to the core as possible by using a specific streamline in thejet. In this respect, W pp is a proper way which permits to use a super-resolution image and thusto quantify the streamline closer to the core. We note that the radial dependence of W pp ( r ) canbe different from that of W ( r ). This should be in fact an interesting topic to be examined, butdetermining a detailed difference between them is beyond the scope of our present work. Such anadvanced analysis of the jet profile is indeed possible by taking in advantage of the sparse modelingtechnique mentioned above, and hence will be presented in the forthcoming paper.The following interesting features are found in Figure 9. Beyond 0.5–0.6 mas from the core, themeasured jet shape is well described by a parabolic collimation profile. There are several previousmeasurements for jet width and opening angle on this scale. For the jet width, Asada & Nakamura(2012) and Hada et al. (2013) determined a radial profile to be W ( r ) ∝ r . ± . over the distancefrom ∼ ∼
500 mas from the core, and the present W pp ( r ) beyond 0.5 mas is in good agreementwith this dependence. For the opening angle, it was previously measured by Junor et al. (1999)and Biretta et al. (2002), and their finding of φ ∼ ◦ opening angle at 0.5 mas from the core isconsistently seen in our result at the same distance . On the other hand, closer to the core whereour 86 GHz image can access, the jet profile becomes more complicated. From the distance of0.6 mas to 0.3–0.2 mas where the constriction exists, the observed opening angle remains roughlyconstant at φ pp ( r ) ∼ ◦ , indicating the jet having a conical geometry. Then, even closer to thecore within ∼ ∼ ◦ up to & ◦ in this region. This trend can be recognized inFigure 8 with the guide of polar coordinates. If the opening angle is defined with respect to W ( r )instead of W pp ( r ), the value should be even larger.Note that the opening angle may become smaller if the jet launch point has a finite size (crosssection). However, all the previous EHT observations of M87 at 230 GHz constrained a jet-launchsize to be remarkably small (40 µ as; Doeleman et al. 2012; Akiyama et al. 2015). For reference,we superpose the corresponding model of the 230 GHz core on the coordinate origin of Figure 8.As seen from this map, the modification of the opening angle is sufficiently small at our scale of In Junor et al. (1999) and Biretta et al. (2002) the opening angle is defined with respect to the full-width-quarter-maximum (FWQM) on jet intensity slice profiles. This gives a slightly larger opening angle than that in our methodat the same r .
23 –interest, so our assumption should be reasonablly valid. Therefore, this is the largest opening angleever observed in any astrophysical jets as well as in M87 itself.
In Figure 10 we show a result of our 86 GHz polarimetry analysis for the M87 jet. Herewe display the data on February 26, since on this epoch the data quality is better and also wecan perform the more reliable EVPA correction with the help of an external close-in-time VLBAobservation of the calibrator 3C 273 (see Appendix in more detail). The result on February 11 isessentially consistent with Figure 10, although the SNR is lower than that on February 26.Given that the M87 jet is only weakly polarized (or highly depolarized) on pc-to-subpc scalespresumably due to a dense foreground Faraday screen (Zavala & Taylor 2002), the present obser-vation is still challenging to reveal the whole polarimetric structure of this jet. Nevertheless, weindeed detected some significant polarized emission in a few parts of the jet at SNR ∼ − , and its fractional polarization is a level of 3–4%. On the other hand,we detected another polarized feature at 0.4 mas downstream in the jet. While this feature hasa similar polarized intensity (4.9 mJy beam − ) to that of the other one, the observed fractionalpolarization becomes as high as 20%. Since the previous VLBI polarimetric observations of this jet(which are usually made at 15 GHz or lower frequencies) reported a fractional polarization up to11.5% (Junor et al. 2001), this is the highest fractional polarization ever reported on pc-to-subpcscales of this jet. This polarized feature appears to be located at the boundary of S1 and the central“valley” of emission (so-called the “spine” part of the jet).Regarding EVPA, the near-core feature shows a mean EVPA direction along the jet, while thenear-S1 feature indicates a mean EVPA roughly perpendicular to the jet. If this is the case, thecorresponding magnetic-field-vector-polarization angle (MVPA) of the near-core/near-S1 featuresare perpendicular/parallel to the jet axis, respectively. However, we remind that there is still alarge uncertainty in our EVPA correction procedure (∆ χ M87 ∼± ◦ ; see Section 2.2 or Appendix).Moreover, there might be an additional external EVPA rotation if there is a significant amount ofthe foreground Faraday screen toward M87 (e.g., ∆ χ RM , M87 = RM λ . ∼ ◦ if RM ∼ rad m − ).This does not permit us to fix the intrinsic EVPA definitively.
4. Discussion4.1. Jet viewing angle and speed
The viewing angle θ of the M87 jet has been discussed for a long time and still remains as apuzzling issue for this jet. Owen et al. (1989) suggest that the jet is not too far out of the plane 24 –Fig. 10.— VLBA+GBT 86 GHz polarimetric result for M87 on 2014 February 26. The color map,the vectors and the contours indicate the observed polarized intensity, the observed EVPA and thetotal intensity distribution, respectively. The convolving beam is shown at the bottom-right corner.The polarized intensity (color map) is displayed from 4 σ rms noise level. The contours start from1, 2, 2 / , 4... times 2.1 mJy beam − . 25 –of the sky ( θ & ◦ ) based on the patterns of the helically-wrapped filaments seen in the kpc-scaleVLA jet. On pc-to-subpc scales Ly et al. (2007) suggest θ ∼ ◦ based on a proper motion andbrightness ratio measurement with VLBA at 43 GHz. On the other hand, a strong constraint isobtained from optical observations of the active knot HST-1 at ∼ ′′ from the nucleus, whereBiretta et al. (1999) found a superluminal motion up to 6 c , tightly requiring θ to be smaller than19 ◦ from our line of sight. Perlman et al. (2011) also suggest a similar range ( θ ∼ ◦ ) based onthe optical polarization properties of HST-1.In the present study, we estimate θ near the jet base based on the measured apparent motions inthe jet and the counter jet, such that µ j µ cj = β int cos θ − β int cos θ , in the assumption that the bidirectional jet isintrinsically symmetric. If we compare CJ with S1 or S2 (as explained in Section 3.2.4), the propermotion ratio results in µ j /µ cj =1.5–2.5. Searching for a common area with µ j = β int sin θ − β int cos θ , we thenobtain solutions of θ and β int to be θ = 29–45 ◦ and β int = 0 . .
50, respectively. With these rangesof θ and β int , we can also estimate an expected BR, and this results in a range of BR = 2 . . α = − . θ is rather similar to those suggested in Ly et al. (2007),whose measurement was also made at a similar distance from the core. In contrast, the derivedrange of θ is larger than that obtained from the optical HST-1 kinematics.We do not rule out the possibility of the smaller θ as suggested from the HST-1 observations,since there can be still an overlap ( θ > ◦ ) if we allow the maximum proper motion ratio ( µ j /µ cj =4 .
6) within the 1 σ errors. Also, a recent VLBA 43 GHz monitoring program of the inner jet byWalker et al. suggests a fast apparent motion (in the main jet) of ∼ (1–2) c (Walker et al. 2008),which might favor a small viewing angle. Thus the more accurate measurements of proper motionas well as brightness ratio are important in future studies.One issue, however, we should note is that the M87 jet is highly limb-brightened. A commonlyinvoked explanation for limb-brightening structure in relativistic jets is that there is a velocitygradient transverse to the jet, such that the flow speed becomes faster toward the jet centralaxis (e.g., Ghisellini et al. 2005; Nagai et al. 2014). According to this idea, the brighter part of thejet (i.e., the limb/sheath) has a larger Doppler factor δ ( δ ≡ [Γ(1 − β cos θ )] − where Γ is the bulkLorentz factor) to the observer, while the dim part of the jet (i.e., the central spine) has a smaller δ due to the lower beaming relative to the sheath. Since δ ( β ) reaches a maximum at β = cos θ , ifwe consider the case of θ = 10 ◦ –20 ◦ , the faster β yields the larger δ for most of β unless β spine isunrealistically higher than β sheath . This results in δ spine > δ sheath , indicating that the jet brightnessprofile would lead to a ridge-brightened structure. On the other hand, if we consider the case of θ = 29 ◦ –45 ◦ , δ reaches a maximum at β ∼ β , δ starts to decrease,which in principle can reproduce a limb-brightened intensity profile.Therefore, the observed limb-brightened structure of M87 may not be simply explained by atransverse velocity gradient alone if θ = 10 ◦ –20 ◦ , requiring some other process being at work (e.g.,Lobanov & Zensus 2001; Stawarz & Ostrowski 2002; Gopal-Krishna et al. 2007; Zakamska et al. 26 –2008; Clausen-Brown et al. 2011).An alternative hypothesis to accommodate the apparent discrepancy of θ between the innerjet and HST-1 is that the viewing angle of M87 is not constant all the way down the jet. Ourlong-term VLBI monitoring of HST-1 has recently revealed significant variations in the observedposition angles (from PA ∼ ◦ to PA ∼ ◦ ) of the substructures’ trajectories (Giroletti et al. 2012;Hada et al. 2015) . This implies a deprojected (intrinsic) change in direction to be ∼ ◦ (for a fixed θ = 15 ◦ ). Thus, it would not be surprising if the θ of HST-1 is also variable at this level, and thefastest ∼ c speed could be seen when its θ is maximally beamed to us. Such a local misalignmentof θ from the central jet axis can be realized if the HST-1 complex is traveling along a 3-dimensionalhelical trajectory with respect to the central jet axis.From the point of view of the high-energy emission, the jet base of M87 is proposed to be a likelysite of the very-high-enery (VHE) γ -ray production (e.g., Aharonian et al. 2006; Acciari et al. 2009;Abramowski et al. 2012; Hada et al. 2012, 2014). The detection of VHE γ -ray emission usuallyprefers a small viewing angle of the jet, which may also be opposed to the viewing angle derivedabove. However, it would be interesting to note the jet base of M87 has a very wide apparentopening angle up to φ app ∼ ◦ (Figure 8 and Figure 9). If we consider the case of θ = 30 ◦ , theintrinsic opening angle is estimated to be φ int ∼ φ app × sin θ = 50 ◦ . This means that the near sideof the sheath is almost pointing toward us at θ near ∼ ◦ close to the jet base (assuming an axially-symmetric jet). This value is quite similar to the typical viewing angle in blazars. Therefore, theobserved VHE emission from the jet base could be associated with a locally beamed substructurein the near side of the sheath (e.g., Lenain et al. 2008; Giannios et al. 2010).We also note that the bulk flow speed near the jet base may be somewhat faster than whatwe measure in the component speeds. In fact, β int derived from the proper motion ratio µ j /µ cj isapplicable to both bulk and pattern speeds, and one cannot distinguish between these two withsuch measurement alone. A bulk flow speed is more directly related to the brightness ratio BR.The expected BR of 8.7 described above (with β int = 0 . θ = 29 ◦ ) is actually only partiallyconsistent with the present observations (i.e., still inconsistent with the measured BR at 86 GHz),and this situation is essentially the same for the smaller angle of e.g., θ = 15 ◦ (BR = 11). Toreproduce the observed BR in our 43/86 GHz images fully consistently, we suggest that the bulkflow speed needs to be faster than ∼ . c at this location. One of the most intriguing features found from our 86 GHz observation is that the initial jetformation structure evolves in a quite complicated manner; there are multiple stages before the jet Changes in position angle of parsec-scale superluminal components are also seen in another nearby radio galaxy3C 120 (G´omez et al. 2000, 2001; Casadio et al. 2015)
27 –finally reaches the well-defined parabola in the outer scale. The jet is formed with φ app ∼ ◦ ,then rapidly collimated into φ app ∼ ◦ within ∼ ∼ R s ) fromthe core, and subsequently reaches the “constricted” point. From there, the jet expands roughlyconically at φ app ∼ ◦ until ∼ ∼ R s , projected) down the jet, and finally enters thelarge-scale parabola collimation zone. A possible structural change of the jet profile near the blackhole is originally suggested in Hada et al. (2013). While this kind of feature might be self-formedthrough some instabilities (such as the sausage-pinch instability in a magnetohydrodynamic flow;e.g., Begelman 1998), this may also be a signature of the interaction of the jet with the surroundingmedium. In fact, there are growing implications that a pressure support from an external mediumis necessary to build an efficient collimation of a jet (e.g., Nakamura et al. 2006; Komissarov et al.2007, 2009). Therefore, in what follows we discuss whether the formation of the M87 jet on thisscale can be subject to an external effect or not, based on a simple comparison of the pressurebalance between the jet ( p jet ) and the external medium ( p ext ).Here we consider the case that the internal pressure of the M87 jet is approximated by the sumof the leptons ( p ± ) and the magnetic fields ( p B ), and that the relative contribution from protonsis small (Reynolds et al. 1996). In Kino et al. (2014), we examined an allowed range of the energybalance between electrons and magnetic fields at the base of this jet based on the synchrotronself-absorption (SSA) theory, and showed that the radio core at 43 GHz can be highly magnetizedor at most in roughly equipartition (10 − . U e /U B . . unless the power index of the electronenergy distribution is too steep (i.e., unless q > n e ( E e ) ∝ E − qe ). The observed spectralindex for the optically-thin part of the jet in our 24/43/86 images ( α j ∼− . − . q ∼ q ≡ − α + 1 in the present definition) satisfies this condition. In this case, the totalpressure of the jet is predominantly due to the magnetic fields, or the particle pressure is at mostof the same order of magnitude of the magnetic one. Therefore, we can reasonably adopt that p jet ∼ p B = B / π at the 43 GHz core. As for the 86 GHz core, we can similarly estimate its B value based on the SSA formula (the equation (11) in Kino et al. 2014) in combination withthe modelfit parameters listed in Table 2, and this results in B core , ∼ B core , and B core , derived in Kino et al. (2014, 2015), and a magnetically-dominatedstate can be consistently satisfied. Thus adopting this B core , , we obtain p B, core , ∼ . − .For the jet downstream of the core, p B depends on the radial profile of the magnetic fields. In anycase, to support/confine the jet on these scales, the external medium needs to have a pressure thatcan be balanced with the suggested level of of p jet .According to the observed ν − . frequency dependence of the core shift (Hada et al. 2011),the radio core at 86 GHz is estimated to be located at ∼ R s from the black hole on the planeof the sky. This indicates a deprojected distance of the 86 GHz core to be between ∼ R s and The equipartition range cited here (that is originally from Kino et al. (2014)) was derived, for simplicity, in theassumption that the underlying magnetic-field configuration is isotropic (the equation (11) in Kino et al. (2014)).However, the same formula is also applicable for the case of an ordered magnetic-field geometry just by multiplyinga small modification factor of 1 / √
28 – ∼ R s for a range of θ = 15 ◦ or θ = 30 ◦ , respectively. On these scales a likely source of theexternal confinement medium may be the inner part of accretion flow or associated coronal re-gion (McKinney 2006; McKinney & Narayan 2007). Since the accretion rate onto the M87 nucleusis significantly sub-Eddington (Di Matteo et al. 2003), the accretion mode of the M87 black holeis thought to be an advection-dominated, hot accretion flow state (ADAF; e.g., Narayan & Yi1994). As ADAF is geometrically thick and roughly approximated by a spherically symmet-ric structure, such a configuration may be suitable in shaping/confining the initial stage of ajet. While the original ADAF assumed a radially-constant mass accretion rate ( ˙ M ), subsequenttheoretical studies favorably suggest that ˙ M decreases with decreasing radius due to convec-tion (Quataert & Gruzinov 2000) or outflows (Blandford & Begelman 1999). In fact, the recentpolarimetric study for the M87 nucleus at 230 GHz has derived ˙ M at r ∼ R s from the blackhole to be ˙ M < . × − M ⊙ yr − (Kuo et al. 2014), which is more than 100 times smallerthan that measured at the Bondi radius ( ˙ M Bondi ∼ M ⊙ yr − at r Bondi ∼
230 pc ∼ × R s ;Di Matteo et al. 2003). For such modified ADAF flows, Yuan et al. (2012) and Yuan & Narayan(2014) present an updated set of self-similar solutions by taking into account radially-variable ˙ M (i.e., ˙ M ( r ) ∝ r s where s ∼ p ADAF ( r ) ≈ . × α − m − ˙ m BH r − / s dyn cm − , where α visc , m BH and ˙ m BH are the dimen-sionless viscosity parameter, M BH /M ⊙ and ˙ M BH / ˙ M Edd ( ˙ M Edd ≡ L Edd / . c ), respectively . For˙ M BH onto the black hole, here we assume ˙ M BH ( r = R s ) ∼ − M ⊙ yr − together with s = 0 . M ( r ) at r = r Bondi and r = 21 R s consistently satisfy the above observations, and thus would be a reasonable combination of ˙ M BH and s . With these values, we finally obtain p ADAF ( r ) ∼ . α − r − . dyn cm − . Therefore, ifwe assume α visc to be of the order of 10 − (e.g., Yuan et al. 2012), the ADAF pressure results in p ADAF ∼ − at r = 10 R s . Interestingly, this is quite comparable to p jet estimated at the86 GHz core. Hence, we suggest that the pressure support from the inner part of the hot accretionflow may be dynamically important in shaping and confining the launch stage of the M87 jet.Given that the external pressure contribution is significant at the jet base, the observed con-stricted structure at a projected distance of ∼ R s from the core may reflect some importantphysical signature resulting from the jet-surrounding interaction. One possibility is that this fea-ture marks a reconfinement node of the flow (e.g., Daly & Marscher 1988; G´omez et al. 1997;Komissarov et al. 1997; Kohler et al. 2012; Matsumoto et al. 2012; Mizuno et al. 2015). Such aflow reconfinement can be realized if the radial profile of the external pressure downstream of thecore decreases more slowly than that of the jet pressure. Alternatively, this constriction might The self-similar solutions presented in Yuan et al. (2012), Yuan & Narayan (2014) and also in Narayan & Yi(1994) (for the original ADAF) are obtained using a height-integrated system of equations. Thus “ r ” in this casecorresponds to the cylindrical radius, not the spherical radius from the central black hole. However, the effect ofthe vertical integration is examined in detail by Narayan & Yi (1995), and they proved that the height-integratedsolutions are quite accurate approximations (within ∼
10% for the pressure) of the exact (sperically-averaged) solutionsin the limit of advection-dominated state. Hence, we can reasonably treat that r in the height-integrated solutionscorresponds to the radial distance from the black hole.
29 –be a signature of the sudden “breakout” from the central dense confining medium (presumablyADAF/hot corona), analogous to a jet in Gamma-Ray Bursts (e.g., Morsony et al. 2007). If thebreakout is the case, one may constrain a scale height (thickness) of the central confining mediumto be of the order of H ∼ × × (sin θ ) − = 140 R s (for θ = 30 ◦ ).We note that the above discussion still leaves a considerable uncertainty in each parameterspace and the actual force balance at the jet boundary may be more complicated (e.g., if the jetram pressure at the boundary is significant or if there exists some instability in the accretion flow).In particular, the deep 86 GHz images obtained here evoke the following simple question: why is thesuggested hot accretion flow still not detected in emission, despite the fact that an accretion scalewell below 100 R s near the black hole is already imaged? This issue should be explored by futurehigher-sensitivity imaging observations (e.g., including the phased-up ALMA; Fish et al. 2013).Also, the lack of additional 86 GHz images at different epochs with a similar sensitivity does notpermit us to conclude whether the observed constricted feature is a persistent structure or just atemporal one (although we note that the subsequent 43 GHz images in Figure 6 still show a hintof the corresponding feature at the same location). In any case, our simple order-of-magnitudeestimate discussed here implies a non-negligible contribution of the external medium to the initialevolution of the M87 jet. Further VLBI monitoring for the jet base at 86 GHz will enable us toaddress this issue more definitively.
Finally we briefly discuss the polarization structure of the M87 jet. While on kpc scales po-larimetric properties of this jet are intensively studied in radio and optical (e.g., Owen et al. 1989;Perlman et al. 1999; Chen et al. 2011; Perlman et al. 2011), on pc-to-subpc scales the polarizationstructure is still highly uncertain. This is because polarization signals from the M87 inner jet aregenerally quite low presumably due to a dense foreground medium on these scales, as suggested byZavala & Taylor (2002). They found an extreme RM distribution that varies from − × rad m − to > × rad m − in the jet at 20 mas (1.6 pc, projected) from the core. Recently, Kuo et al.(2014) performed the first RM measurement toward the nucleus in the 230 GHz band with the Sub-millimeter Array. Although the angular resolution in their study is limited to ∼ ′′ , they suggesteda value of RM near the central black hole to be within ( − . . × rad m − , by assumingthat the bulk of the 230 GHz emission originates in the jet base within ∼ R s from the black hole.The detection of polarized signals in our 86 GHz VLBI images provides some important insightsinto the close vicinity of the central black hole of M87. First, the detection of a polarized feature at ∼ ∼ R s , projected) from the black hole provide evidence that the RM associated with thisfeature is not larger than a certain value. If we assume that the depolarization has an external originand that the observed fractional polarization degree is given by m ( λ ) = m exp ( − σ λ ) (where m , m and σ RM are observed fractional polarization, intrinsic polarization and standard deviationof the RM fluctuations, respectively; Burn 1966), the maximum | σ RM | is constrained to be ∼ (5.8– 30 –17) × rad m − (adopting ∼ e − – e − damping of m ). This is consistent with | RM | < . × rad m − obtained at a similar scale by Kuo et al. (2014). As a future work more dedicatedVLBI polarimetric studies would be fruitful to reveal the more detailed spatial distribution of RMnear the jet base.The second point is that we detected a fractional polarization of ∼
20% in the region at ∼
5. Summary
We reported results obtained from a new high-sensitivity, high-resolution VLBA+GBT obser-vation of the M87 jet at 86 GHz. We summarize our main results as follows.1. We obtained a high-quality image of the jet-launch region of the M87 jet down to ∼ R s nearthe black hole. The resulting image dynamic range is greater than 1500 to 1, which is thehighest ever obtained for this source at 86 GHz. The high-sensitivity image clearly confirmedsome important well-known features of this jet such as a wide-opening angle jet launch, alimb-brightened intensity profile, a parabola-shape collimation profile and a counter jet. Thelimb-brightened structure is already well developed at < R s (projected) from the core, andthe corresponding apparent full-opening angle near the black hole becomes as broad as ∼ ◦ .This is the broadest opening angle ever seen in any astrophysical jets as well as in the M87jet itself. 31 –2. We discovered a complicated jet launch shape near the black hole ( r . R s ) in our 86 GHzimage, indicating multiple collimation stages before the jet finally reaches the well-definedparabola profile in the larger scale. In particular, there is a “constricted” structure at ∼ R s (projected) from the core, where the jet cross section is locally shrinking. We suggest thatan external pressure support/contribution from the inner part of accretion flow (presumablyan ADAF type hot accretion flow or associated corona) may be dynamically significant inshaping and confining the M87 jet on this scale.3. Complementing our 86 GHz data with close-in-time multi-epoch lower-frequency data, we de-tected proper motions in both the main jet and the counter jet, which were all subrelativistic.A mean speed of the main jet components were β app ∼ β app ∼ θ ∼ ◦ , although the more dedicatedproper motion studies are necessary.4. We reported on the first VLBI 86 GHz polarimetric result of the M87 jet. While it is stillchallenging to reveal the entire polarimetric property of this jet, we detected some polarizedfeatures near the jet base at this frequency. The detection of the polarization signals atthis frequency implies that the magnitude of the rotation measure toward these features arenot larger than ∼ (5–17) × rad m − , which is consistent with the result reported in therecent 230 GHz polarimetric study. Moreover, one of the polarized features has an observedfractional polarization up to ∼ A. The calibrator 3C 273
Here we describe our data analysis and imaging for the bright quasar 3C 273 obtained by the86 GHz VLBA+GBT program plus an additional 43 GHz archival dataset. Checking the total-intensity and polarimetric status of this source is important for validating the results for M87(particularly the polarimetric one).3C 273 was observed on 2014 February 11 and 26 with VLBA+GBT at 86 GHz as an overallcalibrator of our program. In the top-left panel of Figure 11 we show an 86 GHz total intensityimage of 3C 273 taken on February 26. The observed jet structure consists of the bright core withseveral discrete knots down the jet. We achieved a dynamic range of 670 to 1, allowing a firmdetection of the weaker emission down to ∼ ∼ × rad m − ) for the inner jet at ∼ χ = χ − χ . = 50 ◦ . This indicates that there is anabsolute difference in RM between P3 and P1 of at least | RM | = 2 . × rad m − . This level ofRM is very similar to that reported by Attridge et al. (2005). It seems that there is an additional(i.e., another ∼ ◦ ) rotation of EVPA at the near-core side of P1 in the 86 GHz image. This impliesthe further increase of RM toward the core, which is consistent with the non-detection of polarizedsignal in the corresponding region at 43 GHz.To estimate the absolute uncertainty of our EVPA correction procedure, it is necessary toknow the absolute difference of P3’s EVPA between 43 and 86 GHz. Unfortunately this is difficultto derive from the present data alone. If P3 has a similar level of RM as seen toward P1, the EVPAdistribution shown in the bottom-left panel of Figure 11 will have another ∼ ◦ rotation. However,this level of uncertainty should be regarded as an upper limit, and the actual uncertainty shouldbe smaller than this value, since P3 is located further downstream and shows the higher fractionalpolarization than that of P1, favoring the less amount of Faraday screen toward P3. On the largerscales (i.e., from a few to 10 mas from the core), several authors reporte RM values of hundreds to afew thousands rad m − (e.g., Asada et al. 2002; Zavala & Taylor 2005). Looking into the adjacentcomponent P2, there is actually a slight ( ∼ ◦ ) difference in EVPA between 43 and 86 GHz. Thissuggests an absolute difference in | RM | between P3 and P2 of 3 . × rad m − , which is just inbetween the values for P1 and the outer jet in the literature. Therefore, we regard that a levelof ∼ × rad m − would be a reasonable measure of | RM | at P3, and thus a likely uncertaintyof our 86 GHz EVPA correction relative to the 43 GHz data would be ∆ χ RM ∼ ◦ . In summary,adopting that the EVPA uncertainty of the Boston 43 GHz map is ∆ χ ∼ ◦ (Jorstad et al.2005), we estimate that a potential total uncertainty of EVPA in our 3C 273’s 86 GHz images is alevel of ∆ χ ∼ ∆ χ + ∆ χ RM ∼ ± ◦ . 34 – P1P2 P3 P1P2 P3 (EVPA reference)
Fig. 11.— Total intensity and polarization images of 3C 273 at 86 and 43 GHz. The top-left panelis a VLBA+GBT 86 GHz total intensity contour image observed on 2014 February 26. The beamsize and the peak intensity are 0 . × .
11 mas at PA − ◦ and 805 mJy beam − , respectively. Thetop-right panel is a 43 GHz VLBA total intensity contour image observed on 2014 February 25. Thebeam size and the peak intensity are 0 . × .
18 mas at PA − ◦ and 1100 mJy beam − , respectively.The bottom-left and the bottom-right panels are a close-up view of the inner jet at 86 and 43 GHz,respectively. In the bottom two panels, the images are restored with a common beam of 0 . × .
18 mas in PA − ◦ , corresponding to the synthesized beam of the 43 GHz data. The color map,the vectors and the contours indicate the observed polarized intensity, the observed EVPA and thetotal intensity distribution, respectively. P1, P2 and P3 indicate the polarized components whichare consistently identified in both images. For the upper panels contours start from -1, 1, 2, 2 / ,4, 2 / ... times 3 σ rms level of each image (1 σ = 1 . − and 3.3 mJy beam − at 86 and43 GHz, respectively), while for the bottom panels contours are 1, 2, 4, 8... times 3.6 mJy beam − and 9.9 mJy beam − at 86 and 43 GHz, respectively. 35 – REFERENCES
Abramowski, A., Acero, F., Aharonian, F., et al. 2012, ApJ, 746, 151Acciari, V. A., Aliu, E., Arlen, T., et al. 2009, Science, 325, 444Aharonian, F., Akhperjanian, A. G., Bazer-Bachi, A. R., et al. 2006, Science, 314, 1424Agudo, I., Bach, U., Krichbaum, T. P., et al. 2007, A&A, 476, L17Akiyama, K., Lu, R.-S., Fish, V. L., et al. 2015, arXiv:1505.03545Asada, K., Inoue, M., Uchida, Y., et al. 2002, PASJ, 54, L39Asada, N. & Nakamura, M. 2012, ApJ, 745, L28Attridge, J. M., Wardle, J. F. C., & Homan, D. C. 2005, ApJ, 633, 85Begelman, M. C., Blandford, R. D., & Rees, M. J. 1984, RvMP, 56, 255Begelman, M. C. 1998, ApJ, 493, 291Biretta, J. A., Sparks, W, B., & Macchetto, F. 1999, ApJ, 520, 621Biretta, J. A., Junor, W., & Livio, M. 2002, New. Astron. Rev., 46, 239Blakeslee, J. P., Jord´an, A., Mei, S., et al. 2009, ApJ, 694, 556Blandford, R. D., & Rees, M. J. 1974, MNRAS, 169, 395Blandford, R. D., & K¨onigl, A. 1979, ApJ, 232, 34Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433Blandford, R. D., & Payne, D. G. 1982, MNRAS, 199, 883Blandford, R. D., & Begelman, M. C. 1999, MNRAS, 303, L1Boccardi, B., Krichbaum, T., Bach, U., et al. 2015, arXiv:1504.01272Broderick, A. E, & Loeb, A. 2009, ApJ, 697, 1164Burn, B. J. 1966, MNRAS, 133, 67Casadio, C., G omez, J.-L., Grandi, P., et al. 2015, arXiv:1505.03871Chen, Y. J., Zhao, G.-Y., & Shen, Z.-Q. 2011, MNRAS, 416, L109Curtis, H. D. 1918, Publications of Lick Observatory, 13, 31Clausen-Brown, E., Lyutikov, M., & Kharb, P. 2011, MNRAS, 415, 2081 36 –Daly, R. A., & Marscher, A. P. 1988, ApJ, 334, 539Deller, A. T., Brisken, W. F., Phillips, C. J., et al. 2011, PASP, 123, 275Di Matteo, T., Allen, S. W., Fabien, A. C., & Wilson, A. S. 2003, ApJ, 582, 133Doeleman, S. S., Fish, V. L., Schenck, D. E., et al. 2012, Science, 338, 355Dodson, R., Edwards, P. G., & Hirabayashi, H. 2012, PASJ, 58, 243Fish, V., Alef, W., Anderson, J., et al. 2013, arXiv:1309.3519Fomalont, E. B. 1999, in ASP Conf. Ser. 180, Synthesis Imaging in Radio Astronomy II, A Collectionof Lectures from the Sixth NRAO/NMIMT Synthesis Imaging Summer School, ed. G. B.Taylor, C. L. Carilli, and R. A. Perley (San Francisco, CA: ASP), 301Ford, H. C., Harms, R. J., Tsvetanov, Z. I., et al. 1994, ApJ, 435, 27Gebhardt, K., & Thomas, J. 2009, ApJ, 700, 1690Gebhardt, K., Adams, J., Richstone, D., Lauer, T. R., Faber, S. M., 2011, ApJ, 729, 119Ghisellini, G., Tavecchio, F., & Chiaberge, M. 2005, A&A, 432, 401Giannios, D., Uzdensky, D. A., & Begelman, M. C. 2010, MNRAS, 402, 1649Giroletti, M., Giovannini, G., Cotton, W. D., et al. 2008, A&A, 488, 905Giroletti, M., Hada, K., Giovannini, G., et al. 2012, A&A, 538, L10G´omez, J.-L., Mart´ı, J. M., Marscher, A. P., Ib´a˜nez, J. M., & Alberdi, A. 1997, ApJ, 482, L33G´omez, J.-L., Marscher, A. P., Alberdi, A., Jorstad, S. G., & Garc´ıa-Mir´o, C., 2000, Science, 289,2317G´omez, J.-L., Marscher, A. P., Alberdi, A., Jorstad, S. G., & Agudo, I., 2001, ApJ, 561, L161Gopal-Krishna, Dhurde, S., Sircar, P., & Wiita, P. J. 2007, MNRAS, 377, 446Hada, K., Doi, A., Kino, M., et al. 2011, Nature, 477, 185 (H11)Hada, K., Kino, M., Nagai, H., et al. 2012, ApJ, 760, 52Hada, K., Kino, M., Doi, A., et al. 2013, ApJ, 775, 70Hada, K., Giroletti, M., Kino, M., et al. 2014, ApJ, 788, 165Hada, K., Giroletti, M., Giovannini, G., et al. 2015, arXiv:1504.01808Harms, R. J., Ford, H. C, Tsvetanov, Z. I., et al. 1994, ApJ, 435, 35 37 –Harris, D. E., Cheung C. C., Biretta, J. A., et al. 2006, ApJ, 640, 211Hirabayashi, H., Hirosawa, H., Kobayashi, H., et al. 1998, Science, 281, 1825Hodgson, J. A., Krichbaum, T. P., Marscher, A. P., et al. 2015, arXiv:1504.03136Honma, M., Akiyama, K., Uemura, M., & Ikeda, S. 2014, PASJ, 66, 95Hughes, D., J´auregui Correa, J.-C., Schloerb, F. P., et al. 2010, Proc. SPIE, 7733, 31Jorstad, S. G., Marscher, A. P., Lister, M. L., et al. 2005, 130, 1418Junor, W., Biretta, J. A., & Livio, M. 1999, Nature, 401, 891Junor, W., Biretta, J. A., & Wardle, J. F. C. 2001, in IAU Symp. 205, Galaxies and their Con-stituents at the Highest Angular Resolutions, ed. R. T. Schilizzi (San Francisco, CA: ASP),136Kardashev, N. S., Khartov, V. V., Abramov, V. V., et al. 2013, ARep, 57, 153Kino, M., Takahara, F., Hada, K., & Doi, A. 2014, ApJ, 786, 5Kino, M., Takahara, F., Hada, K., et al. 2015, ApJ, 803, 30Kohler, S., Begelman, M. C., & Beckwith, K. 2012, MNRAS, 422, 2282Komissarov, S. S., & Falle, S. A. E. G. 1997, MNRAS, 288, 833Komissarov, S. S., Barkov, M. V., Vlahakis, N., & K¨onigl, A. 2007, MNRAS, 380, 51Komissarov, S. S., Vlahakis, N., K¨onigl, A., & Barkov, M. V. 2009, MNRAS, 394, 1182K¨onigl, A. 1981, ApJ, 243, 700Kovalev, Y. Y., Lister, M. L., Homan, D. C, & Kellermann, K. I. 2007, ApJ, 668, L27Koyama, S., Kino, M., Giroletti, M., et al. 2015, arXiv:1505.04433Krichbaum, T. P., Graham, D. A., Bremer, M., et al. 2006, JPhCS, 54, 328Kuo, C. Y., Asada, K., Rao, R., et al. 2014, ApJL, 783, L33Laing, R. A. 1980, MNRAS, 193, 439Lee, S.-S., Lobanov, A. P., Krichbaum, T. P., et al. 2008, A&A, 136, 159Lenain, J.-P., Boisson, C., Sol, H., & Katarzy´nski, K. 2008, A&A, 478, 111Lepp¨anen, K. J., Zensus, J. A., & Diamond, P. J. 1995, AJ, 110, 2479 38 –Lobanov, A. P., & Zensus, J. A. 2001, Science, 294, 128Ly, C., Walker, R. C., & Wrobel, J. M. 2004, AJ, 127, 119Ly, C., Walker, R. C., & Junor, W. 2007, ApJ, 660, 200Macchetto, F., Marconi, A., Axon, D. J., et al. 1997, ApJ, 489, 579Mart´ı-Vidal, I., Krichbaum, T. P., Marscher, A., et al. 2012, A&A, 542, A107Matsumoto, J., Masada, Y., & Shibata, K. 2012, ApJ, 751, 140McKinney, J. C. 2006, MNRAS, 368, 1561McKinney, J. C., & Narayan, R. 2007, MNRAS, 375, 513McKinney, J. C., & Blandford, R. 2009, MNRAS, 394, L126McKinney, J. C., Tchekhovskoy, A., & Blandford, R. 2012, MNRAS, 423, 3083Mizuno, Y., Gomez, J.-L., Nishikawa, K., et al. 2015, arXiv:1505.00933Molina, S. N., Agudo, I., G´omez, J. L., et al. 2014, A&A, 566, 26Morsony, B. J., Lazzati, D., & Begelman, M. C. 2007, ApJ, 665, 569Nagai, H., Haga, T., Giovannini, G., et al. 2014, ApJ, 785, 53Nakamura, M., Li, H., & Li, S. 2006, ApJ, 652, 1059Nakamura, M., & Asada, K. 2013, ApJ, 775, 118Narayan R, & Yi I. 1994, ApJ, 428, L13Narayan R, & Yi I. 1995. ApJ, 444, 231Owen, F. N., Hardee, P. E., & Cornwell, T. J. 1989, ApJ, 340, 698Perlman, E. S., Biretta, J. A., Zhou, F., Sparks, W. B., & Macchetto, F. D. 1999, AJ, 117, 2185Perlman, E. S., Biretta, J. A., Sparks, W. B., Macchetto, F. D., & Leahy, J. P. 2001, ApJ, 551, 206Perlman, E. S., Adams, S. C., Cara, M., et al. 2011, ApJ, 743, 119Quataert, E., & Gruzinov, A. 2000, ApJ, 539, 809Reynolds, C. S., Fabian, A. C., Celotti, A., & Rees, M. J. 1996, MNRAS, 283, 873Rioja, M., & Dodson, R. 2011, AJ, 141, 114Roberts, D. H., Wardle, J. F. C., & Brown, L. F. 1994, ApJ, 427, 718 39 –Shepherd, M. C., Pearson, T. J., & Taylor, G. B. 1994, BAAS, 26, 987Stawarz, L., & Ostrowski, M. 2002, ApJ, 578, 763Tchekhovskoy, A., Narayan, R., & McKinney, J. C. 2011, MNRAS, 418, L79Walker, R. C., Ly, C., Junor, W., & Hardee, P. J. 2008, JPhCS, 131, 012053Walsh, J. L., Barth, A. J., Ho, L. C., & Sarzi, M. 2013, ApJ, 770, 86Yuan, F., Wu, M., & Bu, D. 2012, ApJ, 761, 129Yuan, F., & Narayan, R. 2014, ARA&A, 52, 529Zakamska, N. L., Begelman, M. C., & Blandford, R. D. 2008, ApJ, 679, 990Zavala, R. T., & Taylor, G. B. 2002, ApJ, 566, L9Zavala, R. T., & Taylor, G. B. 2005, ApJ, 626, L73