Issue 
A&A
Volume 673, May 2023



Article Number  A159  
Number of page(s)  16  
Section  Extragalactic astronomy  
DOI  https://doi.org/10.1051/00046361/202142988  
Published online  24 May 2023 
Spectral analysis of a parsecscale jet in M 87: Observational constraint on the magnetic field strengths in the jet^{⋆,}^{⋆⋆}
^{1}
Department of Astronomy, Yonsei University, Yonseiro 50, Seodaemungu, Seoul 03722, Republic of Korea
email: hwro@yonsei.ac.kr
^{2}
Korea Astronomy & Space Science Institute, Daedeokdaero 776, Yuseonggu, Daejeon 34055, Republic of Korea
^{3}
National Astronomical Observatory of Japan, 2211 Osawa, Mitaka, Tokyo 1818588, Japan
^{4}
Kogakuin University of Technology & Engineering, Academic Support Center, 26651 Nakano, Hachioji, Tokyo 1920015, Japan
^{5}
University of Science and Technology, Gajeongro 217, Yuseonggu, Daejeon 34113, Republic of Korea
^{6}
Mizusawa VLBI Observatory, National Astronomical Observatory of Japan, 212 Hoshigaoka, Mizusawa, Oshu, Iwate 0230861, Japan
^{7}
Department of Astronomical Science, The Graduate University for Advanced Studies (SOKENDAI), 2211 Osawa, Mitaka, Tokyo 1818588, Japan
^{8}
Institute of Astronomy and Astrophysics, Academia Sinica, 11F of AstronomyMathematics Building, AS/NTU No. 1, Sec. 4, Roosevelt Rd, Taipei 10617, Taiwan, ROC
^{9}
Department of Physics and Astronomy, Seoul National University, 1 Gwanakro Gwanakgu, Seoul 08826, Republic of Korea
^{10}
Department of General Science & Education, National Institute of Technology, Hachinohe College, 161 Tamonoki, Uwanotai, Hachinohe City, Aomori 0391192, Japan
^{11}
Research Center for Intelligent Computing Platforms, Zhejiang Laboratory, Hangzhou 311100, PR China
^{12}
TsungDao Lee Institute, Shanghai Jiao Tong University, Shanghai 201210, PR China
^{13}
Department of Physics and Astronomy, Sejong University, 209 Neungdongro, Gwangjingu, Seoul 05006, Republic of Korea
^{14}
Institute for Cosmic Ray Research, The University of Tokyo, 515 Kashiwanoha, Kashiwa, Chiba 2778582, Japan
^{15}
Shanghai Astronomical Observatory, Chinese Academy of Sciences, 80 Nandan Road, Shanghai 200030, PR China
^{16}
SNU Astronomy Research Center (SNUARC), Seoul National University, 1 Gwanakro Gwanakgu, Seoul 08826, Republic of Korea
^{17}
Department of Physics, Faculty of Science, University of Malaya, 50603 Kuala Lumpur, Malaysia
^{18}
Department of Astronomy and Atmospheric Sciences, Kyungpook National University, Daegu 702701, Republic of Korea
^{19}
MaxPlanckInstitut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany
^{20}
Graduate School of Science, Osaka Metropolitan University, Osaka 5998531, Japan
^{21}
The Research Institute of Time Studies, Yamaguchi University, Yoshida 16771, Yamaguchicity, Yamaguchi 7538511, Japan
^{22}
Key Laboratory of Radio Astronomy, Chinese Academy of Sciences, Nanjing 210008, PR China
^{23}
Instituto de Astrofísica de Andalucía – CSIC, Glorieta de la Astronomía s/n, 18008 Granada, Spain
^{24}
Graduate School of Sciences and Technology for Innovation, Yamaguchi University, 16771 Yoshida, Yamaguchi, Yamaguchi 7538511, Japan
^{25}
Joint Institute for VLBI ERIC, 7991 PD Dwingeloo, The Netherlands
^{26}
Tokyo Electron Technology Solutions Limited, Iwate 0231101, Japan
^{27}
Massachusetts Institute of Technology Haystack Observatory, 99 Millstone Road, Westford, MA 01886, USA
^{28}
Black Hole Initiative at Harvard University, 20 Garden Street, Cambridge, MA 02138, USA
^{29}
Xinjiang Astronomical Observatory, Chinese Academy of Sciences, Urumqi 830011, PR China
^{30}
Toyo University, 52820 Hakusan, Bunkyoku, Tokyo 1128606, Japan
^{31}
Niigata University, 8050 Ikarashi 2nocho, Nishiku, Niigata 9502181, Japan
^{32}
National Astronomical Research Institute of Thailand (Public Organization), 260 Moo 4, T. Donkaew, A. Maerim, Chiangmai 50180, Thailand
^{33}
Department of Astronomy, The University of Tokyo, 731 Hongo, Bunkyo, Tokyo 1130033, Japan
^{34}
Department of Physics, UNIST, Ulsan 44919, Republic of Korea
^{35}
Basic Science Research Institute, Chungbuk National University, Chungdaero 1, SeowonGu, Cheongju, Chungbuk 28644, Republic of Korea
Received:
17
December
2021
Accepted:
21
November
2022
Context. Because of its proximity and the large size of its black hole, M 87 is one of the best targets for studying the launching mechanism of active galactic nucleus jets. Currently, magnetic fields are considered to be an essential factor in the launching and accelerating of the jet. However, current observational estimates of the magnetic field strength of the M 87 jet are limited to the innermost part of the jet (≲100 r_{s}) or to HST1 (∼10^{5} r_{s}). No attempt has yet been made to measure the magnetic field strength in between.
Aims. We aim to infer the magnetic field strength of the M 87 jet out to a distance of several thousand r_{s} by tracking the distancedependent changes in the synchrotron spectrum of the jet from highresolution very long baseline interferometry observations.
Methods. In order to obtain highquality spectral index maps, quasisimultaneous observations at 22 and 43 GHz were conducted using the KVN and VERA Array (KaVA) and the Very Long Baseline Array (VLBA). We compared the spectral index distributions obtained from the observations with a model and placed limits on the magnetic field strengths as a function of distance.
Results. The overall spectral morphology is broadly consistent over the course of these observations. The observed synchrotron spectrum rapidly steepens from α_{22 − 43 GHz} ∼ −0.7 at ∼2 mas to α_{22 − 43 GHz} ∼ −2.5 at ∼6 mas. In the KaVA observations, the spectral index remains unchanged until ∼10 mas, but this trend is unclear in the VLBA observations. A spectral index model in which nonthermal electron injections inside the jet decrease with distance can adequately reproduce the observed trend. This suggests the magnetic field strength of the jet at a distance of 2−10 mas (∼900 r_{s} − ∼4500 r_{s} in the deprojected distance) has a range of B = (0.3−1.0 G)(z/2mas)^{−0.73}. Extrapolating to the Event Horizon Telescope scale yields consistent results, suggesting that the majority of the magnetic flux of the jet near the black hole is preserved out to ∼4500 r_{s} without significant dissipation.
Key words: galaxies: active / galaxies: individual: M 87 / galaxies: jets / radio continuum: galaxies / relativistic processes / techniques: interferometric
Movie is available at https://www.aanda.org
A copy of the images is available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr (130.79.128.5) or via https://cdsarc.cds.unistra.fr/vizbin/cat/J/A+A/673/A159
© The Authors 2023
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
M 87 is known as the best target for investigating the active galactic nucleus (AGN) jet formation mechanism due to its proximity (distance = 16.7 Mpc; Blakeslee et al. 2009) and its large central supermassive black hole (SMBH; Macchetto et al. 1997; Gebhardt & Thomas 2009; Gebhardt et al. 2011; Walsh et al. 2013). Recently, the Event Horizon Telescope (EHT) successfully imaged the first ever black hole shadow (Event Horizon Telescope Collaboration 2019, 2021) and firmly determined the mass of the central SMBH in M 87 (M_{•} = 6.5 ± 0.7 × 10^{9} M_{⊙}; Event Horizon Telescope Collaboration 2019). At the distance of 16.7 Mpc, 1 milliarcsecond (mas) corresponds to ≈130 r_{s} ≈ 0.08 parsec, where r_{s} is the Schwarzschild radius of the SMBH. M 87 hosts a prominent radio jet that is detected from the radio to TeV γrays (e.g., Abramowski et al. 2012; EHT MWL Science Working Group 2021).
The current leading model for jet launching describes a jet that is driven by magnetic force and then accelerated by the relatively slow conversion of magnetic energy into kinetic energy (e.g., Blandford & Znajek 1977; McKinney 2006; Komissarov et al. 2007; Nakamura et al. 2018; Chatterjee et al. 2019). One way to test this scenario is to determine the radial profiles of the velocity field and the magnetic field along the jets. Over the last 15 years or so, the velocity field of the M 87 jet has been intensively explored via multiepoch very long baseline interferometry (VLBI) observations that have allowed us to directly probe the kinematics of the jet components (e.g., Kovalev et al. 2007; Ly et al. 2007; Asada et al. 2014; Mertens et al. 2016). Furthermore, in order to avoid errors in the identification of knot motions due to low cadence, highcadence observations were made and the global profile of the M 87 jet velocity field from 10^{2} r_{s} to 10^{7} r_{s} was revealed in great detail by Park et al. (2019a). Interestingly, the obtained velocity field data showed an apparent discrepancy with those predicted by general relativistic magnetohydrodynamics simulations (e.g., McKinney 2006; Nakamura et al. 2018), posing a new challenge to the jet formation mechanism. Future VLBI observations with higher angular resolution and higher cadence may allow us to address this issue.
Meanwhile, there has been recent progress in constraining the magnetic field strength of the M 87 jet (for a review, see Hawley et al. 2015). In the bright knot HST1, which is located in the ∼3.8 × 10^{5} r_{s} downstream of the SMBH, the magnetic field strength is constrained by the synchrotron cooling time in the Xray energy band (Harris et al. 2003, 2009). In the region within ∼100 r_{s}, there has been major progress in millimeter and submillimeter VLBI observations, which successfully measured the size and the flux of the jet base (the radio core), where the optical depth for synchrotron selfabsorption (SSA) becomes unity at the observed frequency (e.g., Doeleman et al. 2012; Hada et al. 2013, 2016). From this, the magnetic field strength of the radio core was estimated. The estimated strength suggests that the energy density of the magnetic field dominates the nonthermal electron energy density near the SMBH, which in turn suggests that magnetic fields are important for jet formation (Kino et al. 2014, 2015). The successful measurement of the frequencydependent position of the core (the socalled coreshift) provides another way of constraining the magnetic field strength of the jet (Hada et al. 2011; Zamaninasab et al. 2014; Zdziarski et al. 2015). Most recently, from images of the central black hole’s photon ring and its polarization distribution obtained from EHT observations (Event Horizon Telescope Collaboration 2019, 2021), the magnetic field strength of the surrounding plasma was constrained.
However, our understanding of the magnetic field profile along the jet is currently limited. This is because only the magnetic field strengths of notable individual components have been estimated so far. The magnetic field strength between the radio core and HST1 remains unknown. In particular, there has been no attempt to date to measure the magnetic field strength of the extended jet.
In this study we explore the magnetic field properties of the M 87 jet using spectral index maps. Although there have been many previous studies that have published spectral index maps of M 87 (Zavala & Taylor 2003; Dodson et al. 2006; Ly et al. 2007; Hada et al. 2011, 2012; Pushkarev & Kovalev 2012; Niinuma et al. 2014; Hovatta et al. 2014; Asada et al. 2016; Kravchenko et al. 2020), there has been no dedicated study trying to limit physical quantities using spectral index maps. In order to thoroughly constrain and characterize the radial profiles of physical quantities along the M 87 jet, we set up a dedicated observation program using the KVN and VERA Array (KaVA). After commissioning observations (Niinuma et al. 2014), M 87 has been observed as part of the KaVA AGN Large Program (hereafter, LP) since 2016.
This paper is organized as follows. In Sect. 2 we describe the multifrequency observations and data reduction of the KaVA LP and other archival data. In Sect. 3 we present the spectral index maps of M 87 between 22 and 43 GHz. In Sect. 4 we discuss a model of the spectral index distribution that can explain the observed spectral index profile of M 87 and constrain the magnetic field distribution of the jet. Furthermore, we compare our measurements to those from previous studies. Throughout the paper, we define the sign of the spectral index as S ∝ ν^{+α}.
2. Observations and data reduction
2.1. KaVA observations
M 87 observations using KaVA were performed nine times from February to June 2016. Each epoch consists of two sessions, at 22 and 43 GHz, which were observed within 1 to 2 days of each other. Each observation was allocated every 2−3 weeks. Among them, the data on June 1 were excluded due to antenna problems at Mizusawa and Ishigaki. As a result, eight out of the nine data sets are used for the spectral analysis. At each frequency, the total observation time at each epoch is approximately 7 h, and M 87’s onsource time is about 4 h and 30 min. More information about the observations is described in Park et al. (2019a). Data reduction was performed following the process described in Park et al. (2019a). The amplitude is calibrated via an a priori method using an opacitycorrected system temperature and the elevationdependent gain curve for each telescope. The amplitudes were scaled up by a factor of 1.3 to correct for amplitude losses that occurred during correlation (Lee et al. 2015; Hada et al. 2017). However, since the same factor is applied at both frequencies, the spectral index maps do not change after the correction. Imaging with CLEAN and selfcalibration was performed using the Difmap software package (Shepherd et al. 1994). When natural weighting is applied, the typical size of the full width at half maximum of the synthesized beam at 22 and 43 GHz is close to a circular shape with radii of 1.2 mas and 0.6 mas, respectively.
2.2. VLBA archival data
We also included archival quasisimultaneous data from the Very Long Baseline Array (VLBA). A total of five observations, with two in 2010 and three in 2014 were obtained. Details of the observations for the data in 2010 and 2014 are described in Hada et al. (2011) and Hada et al. (2016), respectively. We reperformed the data reduction for these data. Initial data calibrations (a priori amplitude correction, atmospheric opacity correction, fringefitting, and bandpass) were performed using the Astronomical Image Processing System (AIPS; Greisen et al. 2003). For the 2014 observations using the RDBE digital backend system, the autocorrelation amplitude was further corrected after bandpass correction, as suggested in VLBA Memo #37. The subsequent image reconstruction was performed in Difmap based on the usual CLEAN/selfcalibration procedure. The typical beam size is 0.8 × 0.4 mas at 22 GHz, and 0.45 × 0.22 mas at 43 GHz when natural weighing is applied.
3. Results
3.1. Spectral index maps
When observations are made at two frequencies in the same array, the lowfrequency data are distributed over a shorter distance in the (u, v)−plane than the highfrequency data (e.g., see Fig. 1 in Park et al. (2019a) for the (u, v)coverage of KaVA 22 and 43 GHz).
Creating spectral index maps using data with a nonidentical (u, v) range can cause the spectral index to become artificially steep in the extended region, where the sensitivity of the lower frequency is higher than that of higher frequency. Therefore, we excluded data on the long baseline of 43 GHz and the short baseline of 22 GHz, respectively. The resulting (u, v) range of the data is 33−170 Mλ for KaVA and 25−685 Mλ for VLBA. More details on the (u, v) range matching and its effects are described in Appendix A.
Then, we restored all the maps to the same circular beam of 1.2 mas × 1.2 mas (∼160 r_{s}), which is the comparable size of the synthesized beam of KaVA at 22 GHz. This beam size is 1.5 times larger than the original size of the VLBA’s synthesized beam. The selection of the size of the restored beam may affect the final spectral index map and further analysis. We explored this in Appendix B and found that the main conclusions are not changed.
During the selfcalibration process, the absolute coordinate position of the source is lost and the brightest component in the image (the radio core) is shifted to the center of the map. Therefore, an additional step is required to align images of different frequencies. For the M 87 jet, the positions of the radio core at different frequencies have been measured using the phasereferencing technique (Hada et al. 2011). According to their asymptotic relation, r_{RA} = Aν^{−α} + B (α = 0.94 ± 0.09, A = 1.40 ± 0.16 and B = −0.041 ± 0.012), the position difference between the 22 GHz core and the 43 GHz core is mas in right ascension. Assuming the jet position angle is −72° (Walker et al. 2018), then the difference in declination is mas. We took into account this shift when aligning the maps. Recently, it has been shown that the coreshift of AGN jets could change with time (Plavin et al. 2019). However, the phasereferencing observations of the M 87 jet showed little change in 22−43 GHz coreshifts from 2010 to 2019 (≲10 microarcsec; Hada et al. 2012, 2014; Jiang et al. 2021), and thus we use a constant coreshift in our analysis.^{1}
After aligning the images of the two frequencies, the spectral index map is created by calculating the spectral index of each pixel using the following relation:
where ν_{1} and ν_{2} are the two observed frequencies, and S_{ν1} and S_{ν2} are the flux densities at each frequency. The error of the spectral index at each pixel was estimated according to the discussion in Kim & Trippe (2014); The error of the intensity of each pixel (i, j) is considered to be the sum of the systematic error and the thermal random noise (i.e., σ_{Iν, ij} = δ_{ν}I_{ν, ij} + σ_{rmsν}). The factor of systematic amplitude error is empirically assumed as δ ∼ 10% (e.g., Hada et al. 2012; Niinuma et al. 2014; Hovatta et al. 2014; Cho et al. 2017). The thermal random noise (σ_{rmsν}) is obtained from the residual map after the CLEAN process for each individual epoch. Then, the error of the spectral index at each pixel is calculated as follows:
The spectral index error is typically σ_{α} ∼ 0.5 near the core and increases to σ_{α} ≳ 1 in the downstream jet. Because of the uncertainty in the coreshift, additional errors can occur during image alignment. We measured the amount of error due to image alignment for some epochs. As a result, σ_{α} ∼ 0.1 was found in the core region but it is almost negligible in the downstream jet (σ_{α} ≲ 0.02). Since the spectral index error is mainly resulting from the error of the intensity, the uncertainty of image alignment is not considered in this work.
The left and right panels in Fig. 1 are spectral index maps of the M 87 jet using 22 and 43 GHz images overlaid with 22 GHz contours observed by the VLBA and KaVA, respectively. All images have been rotated by −18° in order to align the jet central axis with the horizontal axis (e.g., Walker et al. 2018). The color indicates the spectral index between 22 and 43 GHz at each pixel. We blank pixels where the total intensity is less than 3σ_{rmsν} at one of the observed frequencies. In most maps, spectral indices can be obtained at distances up to 10 mas from the core, which is more than double than had been previously observed (Ly et al. 2007; Niinuma et al. 2014; Kravchenko et al. 2020). In particular, for some KaVA observations, spectral indices can also be obtained from structures located at ≈20 mas. However, we are not concentrating on this structure in this paper.
Fig. 1. Spectral index maps between 22 and 43 GHz obtained from VLBA and KaVA observations. All images have been rotated by −18°. The restoring beam size is 1.2 mas × 1.2 mas, drawn as a black circle in the bottomleft corner. Observing dates are shown to the left of each map. The contours represent the total intensity at 22 GHz. Contours start at 3σ_{rms}, increasing in steps of 2. 
It is notable that in general, the spectral index morphology appears to be quite similar over the course of the observations. While there are differences from epoch to epoch, the general structures remain remarkably consistent.
In the spectral index maps, the M 87 jet shows a flatter spectrum at the core, while the extended jet has a steeper spectrum. These spectral distributions are commonly found in many AGN jets (e.g., O’Sullivan & Gabuzda 2009; Pushkarev & Kovalev 2012; Hovatta et al. 2014). We determined the weighted average spectral indices of the core from a 0.5 × 0.5 mas box centered on the location of the peak flux density. We find α_{core, KaVA} = −0.25 ± 0.50 and α_{core, VLBA} = −0.20 ± 0.55, respectively. The spectral index values obtained from both facilities are identical within the error bars, giving an almost flat or slightly steep spectrum. The spectral index of the M 87 jet core at 22−43 GHz is more or less consistent with what has been previously reported (Hada et al. 2012; Niinuma et al. 2014; Kim et al. 2018a; Kravchenko et al. 2020), although it is not identical^{2}.
We note that there is a sudden spectral index increase at the edge of the jet for a few epochs (March 26, 2014, and May 7, 2014). However, these are possibly spurious features because they only appear intermittently in a few epochs and then disappear immediately. Such features are frequently found in many spectral index maps (e.g., O’Sullivan & Gabuzda 2009; Müller et al. 2011; Fromm et al. 2013; Hovatta et al. 2014; Boccardi 2015). One possible reason is the sparse sampling of the (u, v)coverage, especially on short ranges, leading to strong variations in the outermost structure. In the following analysis, we exclude these high spectral index features at the edge of the jet.
3.2. Spectral index distribution along the M 87 jet
The spectral distribution along the jet can be studied in several ways. One method is to extract the spectral index values along the ridge line of the jet (e.g., Fromm et al. 2013; Hovatta et al. 2014). However, it is not straightforward to apply this method in our study because the M 87 jet is known to consist of double ridges (e.g., Hada et al. 2013, 2016). Instead, we constructed the spectral index distribution of the jet by taking the weighted average of the spectral indices in the direction perpendicular to the jet as it moves down the jet pixel by pixel:
where α_{i} is the pixel value of the spectral index map, σ_{α, i} is the value obtained from the spectral index errors derived from Eq. (2), and n is the number of the pixels used for averaging. When calculating the weighted mean, we used pixel values located within one beam size from the jet axis (±1.2 mas). The reason is that it is wide enough to cover the spectral indices of the two ridges, while at the same time excluding edge structures with spurious features at some epochs (see Sect. 3.1). The error is obtained by multiplying the formal error of the weighted mean by to compensate for the biases due to correlations between pixel values (Park et al. 2019b).
Figure 2 shows the radial distributions of the spectral index as a function of deprojected distance from the SMBH in units of r_{s}. Distributions obtained from KaVA and the VLBA are indicated by red and blue lines, respectively. The positional difference between the map center and the SMBH is corrected using the results of Hada et al. (2011). We assumed a viewing angle of 17° (Walker et al. 2018) to convert the observed projected jet distance to a deprojected distance, where 10 mas corresponds to ∼4500 r_{s}. The spectrum of the innermost jet is relatively flat as was previously described in Sect. 3.1. Downstream of the jet base the spectrum becomes steeper. At a distance of about 2 mas (∼900 r_{s}), the spectral index becomes α ∼ −0.7, which is similar to the average spectral index of the M 87 jet between 22−43 GHz reported in previous studies using the same VLBA data as used in this study (Hada et al. 2012, 2016). This value is also similar to the typical spectral index of jets observed in many AGNs (e.g., Pushkarev & Kovalev 2012; Hovatta et al. 2014).
Fig. 2. Evolution of the spectral index between 22 and 43 GHz along the deprojected distance from the SMBH. Weighted mean values across the jet are used in this distribution (see Sect. 3.2). The projected distance is displayed on the upper axis in mas, and the corresponding deprojected distance is displayed on the lower axis in r_{s}. Red lines represent the spectral index distribution from 2016 KaVA observations and blue lines the spectral index distribution from VLBA observations from 2010 and 2014. The shaded area is the 1σ error of the spectral index. 
Two interesting features are found in Fig. 2 in the region farther than ∼2 mas. First, the spectral indices obtained from both VLBA and KaVA decrease rapidly with distance, reaching a steep spectrum of α ∼ −2.5. This high rate of decline in the spectral index with distance and the resulting very steep spectrum have not been reported in previous works on the M 87 jet. Second, there are hints that the spectral index stops declining at ∼6 mas in the KaVA data. Interestingly, a constant spectral index distribution after a certain distance is often found in many other AGN jets (e.g., Kovalev et al. 2008; Hovatta et al. 2014; Haga et al. 2015; Boccardi 2015; Lisakov et al. 2017; Pushkarev et al. 2019; Baczko et al. 2019; Park et al. 2021). However, it is not clear that this is also seen in the VLBA data. We cannot completely rule out the possibility that the spectral index may continue to decline in this region. There is a time difference of about ∼2 years between the KaVA and the VLBA observations. Therefore, the spectral distribution of the jet in this region could vary with the observation time.
4. Discussion
4.1. Comparison with previous works
In Sect. 3 we were able to obtain spectral index distributions of the M 87 jet up to 10 mas (∼4500 r_{s}) from the core for the first time thanks to the good image quality of KaVA and the VLBA at 22 and 43 GHz. An important characteristic of the distribution is that the spectral index decreases between 2−6 mas (∼900 r_{s} – ∼2500 r_{s}) from α ∼ −0.7 to α ∼ −2.5. In the KaVA data, it does not decrease after ∼6 mas.
The change of the spectral index can be obtained from the change in the energy distribution of nonthermal electrons N(γ), since the slope of the electron distribution (p) and the spectral index is related by α = (p + 1)/2 (Rybicki & Lightman 1979). The evolution of N(γ) as a function of time can be investigated by solving the transfer equation (or continuity equation) including nonthermal electron injection and energy losses.
One previous study discussing the synchrotron spectrum of parsecscale AGN jets has proposed two possibilities for spectral steepening using analytical solutions of the transfer equations for two nonthermal electron injection scenarios (Hovatta et al. 2014). (1) In the case of a constant and continuous injection of nonthermal electrons, the synchrotron spectrum has a break at γ_{break}, and the spectral index steepens by 0.5 before and after the break. (2) If there is no injection after the initial injection, the cutoff energy γ_{max} decreases with time. At the cutoff energy, the synchrotron spectrum decreases exponentially, and thus the spectral index decreases to an arbitrarily low value.
However, both scenarios cannot explain the properties observed in the spectral index distribution of the M 87 jet (Fig. 2). Firstly, the observed decrease in the spectral index is much greater than 0.5. Also, the decrease in the spectral index may not continue and may stop at a certain value, remaining at that value after that.
Therefore, in this study, we attempt to explain the properties of the observed spectral index distribution by numerically solving the transfer equation. Section 4.2 describes the transfer equation, and assumptions for physical parameters to solve the equation numerically, including magnetic field, jet radius, bulk jet velocity, and electron injection. In Sect. 4.3, the importance of electron injection is discussed in order to explain key features of the observed spectral index distribution. In Sect. 4.4 we try to constrain the magnetic field strength and the electron injection of the M 87 jet at the distance of 2−10 mas by comparing the observed spectral index distribution to the model. Then, in Sect. 4.5, we compare our estimate for the magnetic field strength with previous estimations.
4.2. Spectral index model
In the frame of a cross section of the jet comoving with the jet’s flow, the transfer equation is written as (e.g., Ginzburg & Syrovatskii 1964; Longair 2011; Blasi 2013)
Here, N(γ, τ) is the number density of the nonthermal electrons as a function of energy γ and time τ, v is the velocity of the system, Q is the nonthermal electron injection rate, and b is the energy loss rate. In the comoving frame, the bulk velocity along the jet direction is zero and only the expansion (or contraction) rate along the radial direction remains. In our discussion, we assume that adiabatic losses and synchrotron losses are the dominant energy loss processes. Then, b can be written as (e.g., Rybicki & Lightman 1979; Longair 2011)
where and is the coefficient of the adiabatic losses and the synchrotron losses, respectively. Here, R is the radius of the cross section, and B is the magnetic field strength in the comoving frame. We assumed that the pitch angle distribution is rapidly isotropized during synchrotron cooling^{3}.
In order to numerically solve Eqs. (4) and (5) along the jet, several physical quantities should be known as a function of distance: the bulk jet velocity profile (Γ(z), where Γ is the Lorentz factor of the bulk jet velocity), the jet radius profile (R(z)), the magnetic field strength profile (B(z)), and the nonthermal electron injection function (Q(z)). Among them, we take the bulk jet velocity and the jet radius as a function of jet distance from previous observations; Γ(z)∝z^{0.16} (Park et al. 2019a), R(z)∝z^{0.56} (Asada & Nakamura 2012; Hada et al. 2013, 2016; Nakamura et al. 2018). The magnetic field at this distance is assumed to be dominated by the toroidal component, in which case the magnetic field strength is inversely proportional to the radius of the jet as the magnetic flux is conserved (Marscher et al. 2010). Then, the magnetic field strength in the comoving frame is (Lyutikov et al. 2005; Zamaninasab et al. 2014)
where B_{i} is the strength of the magnetic field at z_{i}.
The energy spectrum of the nonthermal electron injection function, Q(γ), is generally assumed to be a power law with a slope p_{inj} limited to the minimum and maximum energies, γ_{inj, min} and γ_{inj, max}. In addition to this, we assume that the number of nonthermal electron injections varies with distance with a power of −q:
where Q_{0} is the nonthermal electron injection rate at z_{i} for electrons of γ = 1. For example, if q = 0, the same amount of nonthermal electrons as Q(z_{i}) is injected along the jet continuously. For q = ∞, there is no additional injection after the initial injection at z_{i}. The above two cases resemble the two scenarios of the spectral steepening suggested by Hovatta et al. (2014). If 0 < q < ∞, nonthermal electrons are injected continuously, but the injection volume decreases with distance. There are indeed theoretical models that suggest that electron acceleration along the jet is possible; one is the boundary shear acceleration between the jet and the surrounding matter (Ostrowski 1998), while the other is KelvinHelmholtz instability driven magnetic reconnection (e.g., Sironi et al. 2021). In Sect. 4.3, we discuss the effect of the electron injections on the spectral index distribution.
Given B_{i} and q, Eqs. (4) and (5) are numerically solved and N(γ, τ) is obtained, and the spectral index distribution α(z) between the two observed frequencies is obtained from N(γ, τ). Throughout the paper, the initial and final distance z_{i} and z_{f} are set as ∼900 r_{s} and ∼4500 r_{s}, which is the deprojected distance of 2 mas and 10 mas of the M 87 jet, respectively. The minimum and maximum energy of the injection function, γ_{inj, min} and γ_{inj, max}, is 10 and 10^{5}, and the energy slope of the injection function is p_{inj} = −2.4 (α_{inj} = −0.7), which is similar to the average value of the observed spectral index at 2 mas (Fig. 2). As Q_{0} does not affect the spectral index we arbitrarily selected Q_{0} as 1.00 × 10^{20}. Details of the model are described in Appendix D.
4.3. The effect of the nonthermal electron injections on the spectral index distribution
If we assume that the flattening of the spectral index seen in Fig. 2 is real, we can test the effect of the nonthermal electron injections on the spectral index distribution. If this assumption is wrong, we discuss the implications of this in Sect. 4.4. We made a number of spectral index distribution models by adopting different values of q. Here, we set the four different cases of the nonthermal electron injection function in the jet to be q = 0, 5, 10, and ∞, while the magnetic field distribution is assumed to be the same B_{i} with a value of 0.3 G.
Figure 3 summarizes the results of the test models. The left panel shows the evolution of electron energy distribution N(γ) as a function of the jet distance drawn as a vertical dashed line in the right panel. The color of the distribution represents the value of q used in the model. A blue and a red vertical line represent the electron energies emitting synchrotron radiation at 22 and 43 GHz in the observer’s frame. The slope of the nonthermal electron injection is p_{inj} = −2.4. The right panel of Fig. 3 shows the spectral index distribution between 22−43 GHz versus jet distance, α_{22 − 43 GHz}(z), obtained from N(γ, z). The horizontal dashed line represents the spectral index of the injection function (α_{inj} = −0.7).
Fig. 3. Results of the modeled electron energy distributions (N(γ)) and the spectral index distributions as a function of distance along the jet. This is done to test the effect of the nonthermal electron injection rate in the jet. Left: Modeled electron number density distributions with different values of q = 0, 5, 10, and ∞. The different colors represent different values of q in the model: q = 0 represents continuous new electron injections in the jet, q = ∞ no new electron injections in the jet, and q = 5 and 10 new electrons injected into the jet but the volume decreases with distance. In all cases, we set the same magnetic field distribution (B_{i} = 0.3 G). The slope of the nonthermal electron injection is p_{inj} = −2.4. The vertical red and blue lines represent the energies of the nonthermal electrons, which emit at 22 and 43 GHz. The location down the jet where the calculations were made is written in the top right corner. Right: Modeled spectral index distributions as a function of distance, which correspond to the slope of the electron number density distributions between the vertical lines in the left panel. The vertical dashed black line is the location down the jet where the calculations were made. The horizontal dashed line is the slope of the injection function (α_{inj} = −0.7). An animation of this figure is available. 
In the absence of electron injection (q = ∞), cutoff energy is shown in N(γ, z) where no higher energy electrons can be found. The cutoff energy moves to lower energy over distance. At z ∼ 2000 r_{s} (middle left panel), the cutoff energy is located in between γ(ν_{obs} = 43 GHz) and γ(ν_{obs} = 22 GHz). Therefore, α_{22 − 43 GHz}(z) diverges to −∞ around this distance. On the other hand, in the constant injection case (q = 0), the electrons still exist beyond the cutoff energy due to the newly injected one. Therefore, the spectral index changes only slightly with distance^{4}.
In the case of q = 5 and 10, the amount of newly injected electrons decreases with distance, thus the slope of N(γ, z) at energies higher than the cutoff energy is steeper than in the case of q = 0. As a result, the energy distribution becomes a broken power law, where the position of the break (γ_{break}) is the same as the cutoff energy in the q = ∞ case. Interestingly, α_{22 − 43 GHz}(z) for q = 5 and 10 show similar trends to the observed distributions, that is, α decreases until ∼5 mas and no longer decreasing or slowly decreases after that. This is because γ_{break} first passes through γ(ν_{obs} = 43 GHz) and then γ(ν_{obs} = 22 GHz), and accordingly, the slope between these two frequencies shifts from the first slope of the broken powerlaw to the second. The larger q, the steeper the second slope of N(γ), so the change in α will be greater. For example, if q = 10 then Δα ∼ 2, which is far greater than the expectation by Hovatta et al. (2014; Δα = 0.5), and similar to our observations.
From this experiment, we found that our model can reproduce the characteristic trend of the spectral index distribution seen in the KaVA observations for cases where the amount of nonthermal electron injections have a distance dependence. The combination of nonthermal electron injections and energy losses creates a broken powerlaw energy distribution in which the slope beyond γ_{break} is determined by the change in the amount of injection with distance (i.e., the parameter q of the electron injection function). The larger the q, the greater the slope beyond γ_{break} and the change in α. The continuous injection and noninjection cases proposed by Hovatta et al. (2014) correspond to two extreme cases in our model (q = 0 and ∞, respectively). In the former case, the change in α is the smallest, and in the latter, the break energy becomes the cutoff energy with no electrons above it, making the spectral index diverge.
4.4. Constraining the magnetic field strength and nonthermal electron injection along the M 87 jet
Now we constrain the physical quantities of the M 87 jet using our numerical model. In order to compare the modeled spectral index distributions with the observations quantitatively, we define an “allowed region” based on the observational spectral index distributions with three parameters: (1) The spectral index range at z_{i} (α_{i}), (2) the location of the end of the steepening region (z_{s, f}), and (3) the spectral index range after the steepening region (α_{f}). The gray shaded region in Fig. 4 is the allowed region of the spectral index of the M 87 jet drawn with parameters of −1.2 < α_{i} < −0.3, z_{s, f} = 5.5 mas, −3.0 < α_{f} < −1.5. More than 90% of the observed profiles lie within the allowed region.
Fig. 4. Constraining the magnetic field strength and nonthermal electron injection of the M 87 jet. Left: Modeled spectral index distributions constrained by the observations. The distributions must lie in the gray shaded allowed region (see Sect. 4.4 for details). Colored lines are the modeled spectral index distributions. Different colors represent different initial magnetic field strengths, B_{i}. The dashed and dotted red line is a modeled spectrum index distribution with q = 10 and q = 11, respectively, with B_{i} = 0.3 G. Right: Parameter space (B_{i}, q). The black region is the allowed values of B_{i} and q constrained by the observations. 
We then create many modeled spectral index distributions for different values of the parameters B_{i} and q. The parameters B_{i} were set from 0.1 G to 1.2 G in 0.1 G increments, and the parameter q was set from 5 to 13. We limited the parameters by excluding models that exist outside the allowed region. The colored lines in the left panel of Fig. 4 are the modeled spectral index distributions that lie completely inside the allowed region. The color indicates the initial magnetic field strength used in the model. In a strong magnetic field, the steepening region is located closer to the core, and the spectral index decreases less. As the magnetic field strength becomes weaker, the spectral index decreases over a longer distance, and the change in the spectral index are greater. For the same magnetic field strength, q essentially determines the final spectral index as discussed in Sect. 4.3. For example, for the same B_{i} = 0.3 G, the model with q = 11 (red dotted line) have steeper spectral index distribution than the q = 10 case (red dashed line). By constraining the modeled distributions to the allowed region, B_{i} is constrained to 0.3 G < B_{i} < 1.0 G and q is constrained to 7 < q < 11. The right panel of Fig. 4 shows the allowed parameter space of B_{i} and q.
We note that our method is highly sensitive to the choice of the allowed region. In Sect. 3.2 we found that a constant spectral index distribution from the VLBA data is less clear, so it is possible that the spectral index distribution could continue to decay farther down the jet (blue lines in Fig. 2). This would suggest two possibilities: (1) Spectral index flattening still occurs but farther downstream. In this case, the model indicates a weaker B_{i}. (2) The spectral index declines indefinitely. In this case, this would suggest that no electrons are injected beyond ∼6 mas. Considering that there are several years between the VLBA and KaVA observations, there is a possibility that the magnetic field strength and the electron injection rate could change with time. Since 2017, we have been monitoring the M 87 jet with the East Asia VLBI Network (EAVN; Cui et al. 2021) at 22 and 43 GHz quasisimultaneously with better sensitivity and (u, v)coverage than KaVA. These observations could help resolve these possibilities.
4.5. Radial magnetic field distribution of the M 87 jet
In Sect. 4.4 we successfully inferred the magnetic field strength from the modeled spectral index distributions in the M 87 jet from 2−10 mas (∼900 r_{s} − ∼4500 r_{s} in the deprojected distance). Therefore, the strength of the magnetic field evolves a function of distance as
where exponent −0.72 comes from Eq. (5).
Figure 5 shows the magnetic field strength of the M 87 jet as a function of the distance from the SMBH. The magnetic field distribution obtained by our models and the observations is drawn as a black area. The figure also summarizes the magnetic field strengths at various locations of the M 87 jet using observations over a wide range of wavelengths. The cross represents the theoretically expected magnetic field strength at the black hole surface based on BlandfordZnajek power estimation (Blandford et al. 2019)^{5}. A recent estimation by Kino et al. (2022) also indicates a similar magnetic field strength in this region. The circles represent magnetic field strengths in radio cores at various VLBI frequencies, obtained based on SSA theory (Kino et al. 2014, 2015; Hada et al. 2012, 2016; Reynolds et al. 1996; EHT MWL Science Working Group 2021)^{6}^{7}. Here, the distance from the core to the SMBH and its error are determined using the result of coreshift (Hada et al. 2011). The error of the magnetic field strength is calculated by applying 25% error to the radio core size and 5% error to the flux of the core (Hada et al. 2013). For the data shown as a ranged value (Kino et al. 2014, 2015), the error bars of the field strength is not displayed. The diamonds indicate field strengths obtained with VLBI observations using the coreshift method (Zamaninasab et al. 2014; Jiang et al. 2021)^{8}, light curve modeling (Acciari et al. 2009), brightness temperature (Kim et al. 2018b), and EHT observations (Event Horizon Telescope Collaboration 2019, 2021). The field strengths at HST1 estimated by applying the synchrotron cooling time to the variability of the Xray band are drawn as squares (Harris et al. 2003, 2009).
Fig. 5. Magnetic field distribution of the M 87 jet. The black area is the distribution of the magnetic field strength estimated in this study. The solid lines with the gray shaded region are the extrapolation of the upper and lower limits of the result in the upstream direction (B(z)∝z^{−0.72}). The dashed lines are the magnetic field distribution which assumes a steeper slope of the jet radius profile in the region closer to the SMBH (B(z)∝z^{−0.92}). The solid lines without a shaded region are the upper and lower bounds of the magnetic field distribution, assuming fast acceleration in the region closer within 2 mas (B(z)∝z^{−1.12}). The dotted line is the extrapolation of the maximum values of the result to the downstream region. Markers are estimated magnetic field strength from previous studies. The blue cross is the theoretically expected magnetic field strength in a BlandfordZnajek jet (Blandford et al. 2019). The circles are magnetic field strength values obtained in VLBI cores derived using SSA theory (Kino et al. 2014, 2015; Hada et al. 2012, 2016; Reynolds et al. 1996; EHT MWL Science Working Group 2021). The diamonds are values obtained from VLBI observations using various methods other than SSA theory (Acciari et al. 2009; Kim et al. 2018b; Zamaninasab et al. 2014; Event Horizon Telescope Collaboration 2019, 2021; Jiang et al. 2021). The squares are the magnetic field strengths at HST1 (Harris et al. 2003, 2009). The horizontal dotteddashed lines are the magnetic field strengths estimated by fitting the broadband SED of the jet using a singlezone SSC emission model (Abdo et al. 2009; MAGIC Collaboration 2020; EHT MWL Science Working Group 2021). 
The solid lines with gray shaded region are the extrapolation of the upper and lower bounds of the magnetic field distributions estimated in this study (B(z)∝z^{−0.72}; Eq. (8)). The dashed lines represent the distribution that reflects a steeper slope of the jet radius profile in the region close to the SMBH as R(z)∝z^{−0.76} (accordingly, B(z)∝z^{−0.92}; Hada et al. 2013). Interestingly, our estimation of B(z) on 900 − 4500 r_{s} scale smoothly connects to values based on synchrotron emission down to the EHT observing scale, keeping B(z)∝z^{−0.72}. The consistency of the slope of the magnetic field distributions with distance may imply that the majority of the magnetic flux of the jet near the SMBH is preserved without serious dissipation out to ∼4500 r_{s}.
Mertens et al. (2016) reported that the M 87 jet shows linear acceleration in the region closer than 2 mas, meaning that the bulk jet velocity profile in this region is proportional to the jet radius (Γ(z)∝R(z)∝z^{0.56}). The solid lines without shaded regions are the upper and the lower bounds of the magnetic field distribution estimated by applying it. In this case, due to the fast acceleration, the magnetic field strength at 2 mas should increase to 0.8 − 2.4 G in order to achieve a similar amount of spectral index steepening. Also, the slope of the magnetic field distribution becomes steeper (B(z)∝z^{−1.12}) in the region closer than 2 mas. As a result, a very strong magnetic field distribution is inferred that deviates from the magnetic field strengths estimated by observations and theory.
Extending our scaling relation downstream (dotted line) appears to be considerably higher than the magnetic field strength of HST1. HST1 is a complex structure in which recollimation occurs, and the magnetic field strength of the upstream jet before HST1 is expected to be smaller than that of HST1. Therefore, the magnetic field strength measured in HST1 should be considered as an upper limit. Given that our extrapolated maximum value to the location of HST1 is roughly an order of magnitude higher than the upper limits set by Xray variability studies in HST1, this suggests that assumptions such as the conservation of the magnetic flux along the jet may be wrong. For example, a recent largescale jet simulation showed that pinch instability can lead to magnetic flux dissipation along the jet (Chatterjee et al. 2019).
The horizontal dotteddashed lines are the magnetic field strengths estimated by fitting a singlezone synchrotron selfCompton (SSC) emission model to the broadband spectral energy distribution (SED) of the jet, specifically focusing on the very high energy (VHE) γrays (Abdo et al. 2009; MAGIC Collaboration 2020; EHT MWL Science Working Group 2021). Comparing these values directly with the magnetic field distributions in our study, the estimated location of the VHE emission is in the downstream jet between z ∼ 4500 r_{s} and ∼10^{5} r_{s}, which differs significantly from the most favored location of VHE flare activities (a few tens of r_{s} from the SMBH; see, e.g., Acciari et al. 2009; Abramowski et al. 2012; Hada et al. 2014; Akiyama et al. 2015). One possible explanation is that VHE γrays are emitted from the outer layer of the radioemitting region (e.g., Tavecchio & Ghisellini 2008). Indeed, a recent simultaneous broadband SED study found that the VHE emitting region is not overlapping with the radioemitting region, and requires a larger area (EHT MWL Science Working Group 2021). However, there are various physical processes that suggest the possibility of the VHE emission not only originating in the innermost jet base but also in the downstream jet (e.g., the HST1 knot; see Abramowski et al. 2012; Rieger & Levinson 2018, and references therein). The origin and location of the VHE emission in the M 87 jet are not yet conclusive.
5. Summary and conclusions
In this work we have investigated the spectral index distribution of the M 87 jet using KaVA LP observations from 2016 and archival VLBA observations from 2010 and 2014. This is the first systematic study of the spectral properties in the innermost regions of the M 87 jet at 22 and 43 GHz, and it allows us to investigate the spectral index out to ∼10 mas from the jet base. Our study can be summarized as follows:

The overall spectral index morphology of the M 87 jet does not appear to have changed significantly over the course of these observations.

The observed spectral index distributions of the M 87 jet from the KaVA observations show a rapid steepening from ∼2 mas until ∼6 mas from α ∼ −0.7 to −2.5. Beyond ∼6 mas, the spectral steepening stops and remains at the same value until ∼10 mas. The trend of the spectral index distribution obtained from the VLBA data showed similar rapid spectral index decay, but the constant spectral index beyond ∼6 mas is less clear. This suggests that the spectral index distribution of the M 87 jet may change with time.

To understand the details of the spectral evolution in the M 87 jet, we modeled the synchrotron spectral index by solving the transfer equation and calculating the evolution of the energy distribution of the nonthermal electron densities in the jet. We find that if the amount of the nonthermal electron injections decreases with distance, the spectral index distribution reproduces the characteristics of the observed distributions.

We compare the model of the spectral index distributions with the allowed region, which is defined based on the observations. As a result, we find that the initial magnetic field strength at 2 mas of 0.3 G < B_{i} < 1.0 G and electron injection function with 7 < q < 11 reproduce the observed spectral index distribution appropriately. We conclude that the distribution of the magnetic field strengths as inferred from the model at 2−10 mas (∼900 r_{s} − ∼ 4500 r_{s} in the deprojected distance) in the M 87 jet is B = (0.3−1.0 G)(z/2 mas)^{−0.72}.

We compared the magnetic field strengths of the M 87 jet with those estimated in previous observational studies. The extrapolation of the magnetic field distribution is in good agreement with the estimates by VLBI in the inner jet region, indicating that the majority of the magnetic flux of the jet near the SMBH is preserved without serious dissipation up to ≲4500r_{s}. However, when extrapolating to the HST1 region, our maximum magnetic field strength is roughly an order of magnitude higher than the upper limits found from Xray variability studies.
In this work, the magnetic field strength of the M 87 jet at a distance of several thousands r_{s} has been estimated for the first time. This study showed that, if the bulk jet velocity profile and radius profile are known, the magnetic field strength and nonthermal electron injection of the extended jet can be inferred from the spectral index profile. Consequently, it provides the magnetic field profile of the jet, which is an important physical quantity to validate the launching model of the jet.
The coreshift can be estimated by performing twodimensional crosscorrelation on optically thin jet regions. However, Pushkarev et al. (2012) showed that this method has a large systematic alignment error in smooth, straight jets that exhibit significant spectral index gradients along the jet, such as the M 87 jet. Indeed, in all epochs, we have failed to estimate reliable values of coreshift of the M 87 jet using this method.
The spectrum at the peak of the M 87 jet reported by Hada et al. (2012) gives a slightly inverted value of α ∼ 0.1. The reason why different spectral indices were obtained despite their data being included in our study is that the size of the convolving beam used in our study is more than twice the size they used. The larger the beam size, the steeper the core spectral index is expected because the jet emission from the optically thin region is included as a structure blending effect. We confirmed that the result of Hada et al. (2012) is recovered when a convolving beam of the same size as used in their study.
See Appendix C for a discussion of when the pitch angle is not randomized.
Note that Δα in our calculations is much less than 0.5, which was suggested by Hovatta et al. (2014). However, Δα = 0.5 is achieved only if the magnetic field strength is constant over jet distance and no adiabatic losses exist. When the magnetic field strength decreases with distance and adiabatic losses exist, Δα is reduced.
This is somewhat larger than those estimated from nonthermal emission properties. This difference is probably due to the difference in the location of the jet of interest between radiationbased and energeticbased magnetic field estimation. The limbbrightened structure of the M 87 jet (e.g., Junor et al. 1999; Kovalev et al. 2007; Hada et al. 2016; Kim et al. 2018b; Walker et al. 2018) suggests that the majority of the radio emission is generated near the boundary of the jet, where the magnetization parameter holds (e.g., Nakamura et al. 2018; Chatterjee et al. 2019; Narayan et al. 2022). On the other hand, the estimation by Blandford et al. (2019) seems to focus on the high σ region.
The upper limit of the magnetic field strength of Kino et al. (2015) ∼124 G critically depends on their assumption of P_{jet} ∼ 5 × 10^{44} erg/s. If smaller a jet power is assumed (P_{jet} ∼ 1 × 10^{44} erg/s), then the maximum magnetic field strength is reduced to ∼65 G.
EHT MWL Science Working Group (2021) used two SED fitting models oriented at different wavelengths. One is oriented at radio wavelengths (model 1) and the other is at VHE emissions (model 2). Since model 1 describes the radio core flux using the conventional SSA model, we plot the value as a circle. We also note that the magnetic field strength of model 1 is a lower limit since they “maximally tuned” their parameters in order to maximize the γray flux.
Zamaninasab et al. (2014) calculated the magnetic field strength at 1 parsec based on the assumption of a conical jet, which is incorrect for the M 87 jet. Since this causes a larger error in the magnetic field strength when farther away from the core, we have shown the magnetic field strength at the 43 GHz core instead. Jiang et al. (2021) measured the magnetic field strength at the 88 GHz core.
Acknowledgments
We thank the anonymous referee for constructive comments to improve the quality of this paper. This work is based on observations made with the KaVA, which is operated by the Korea Astronomy and Space Science Institute and the National Astronomical Observatory of Japan. This research was supported by the Korea Astronomy and Space Science Institute under the R&D program supervised by the Ministry of Science and ICT. H.R., B.W.S., and A.C. acknowledge support from the KASIYonsei DRC program of the Korea Research Council of Fundamental Science and Technology (DRC122KASI). B.W.S. is grateful for the support by the National Research Foundation of Korea (NRF) funded by the Ministry of Science and ICT (MSIT) of Korea (NRF2020K1A3A1A78114060). A.C. acknowledges support from the National Research Foundation of Korea (NRF), grant Nos. 2022R1A2C100298211, and 2022R1A6A1A03053472. J.P. and I.C. acknowledge financial support from the Korean National Research Foundation (NRF) via Global Ph.D. Fellowship Grant 2014H1A2A1018695 and 2015H1A2A1033752, respectively. J.P. is supported by an EACOA Fellowship awarded by the East Asia Core Observatories Association, which consists of the Academia Sinica Institute of Astronomy and Astrophysics, the National Astronomical Observatory of Japan, the Center for Astronomical MegaScience, the Chinese Academy of Sciences, and the Korea Astronomy and Space Science Institute. I.C. acknowledges financial support by the Spanish Ministerio de Economía y Competitividad (grants PID2019108995GBC21). K.H. acknowledges support from the Mitsubishi Foundation grant No. 201911019. This work is partially supported by JSPS/MEXT KAKENHI (grants 18K03656, 18H03721, 19H01943, 18KK0090, 21H01137, 21H04488). S.S.S. is supported by JSPS KAKENHI grant No. 21K03628. K.L. and S.T. acknowledge financial support from the National Research Foundation (NRF) of Korea through grant 2022R1F1A1075115. J.O. was supported by Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (NRF2021R1A6A3A01086420). J.S.K. has been supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (2016R1A5A1013277 and 2020R1A2C1007219). J.Y.K. acknowledges support from the National Research Foundation (NRF) of Korea (grant no. 2022R1C1C1005255). S.S.L. was supported by the National Research Foundation of Korea (NRF) grant funded by the Korea government (MIST), (2020R1A2C2009003). Y.C. is funded by the China Postdoctoral Science Foundation (2022M712084). R.S.L. is supported by the Key Program of the National Natural Science Foundation of China (grant No. 11933007), the Key Research Program of Frontier Sciences, CAS (grant No. ZDBSLYSLH011), the Shanghai Pilot Program for Basic Research – Chinese Academy of Sciences, Shanghai Branch (JCYJSHFY2022013), and the Max Planck Partner Group of the MPG and the CAS. Z.S. acknowledges support from the Major Program of the National Natural Science Foundation of China (Grant No. 11590780, 11590784), and Key Research Program of Frontier Sciences, CAS (grant No. QYZDJSSWSLH057).
References
 Abdo, A. A., Ackermann, M., Ajello, M., et al. 2009, ApJ, 707, 55 [NASA ADS] [CrossRef] [Google Scholar]
 Abramowski, A., Acero, F., Aharonian, F., et al. 2012, ApJ, 746, 151 [NASA ADS] [CrossRef] [Google Scholar]
 Acciari, V. A., Aliu, E., Arlen, T., et al. 2009, Science, 325, 444 [Google Scholar]
 Akiyama, K., Lu, R.S., Fish, V. L., et al. 2015, ApJ, 807, 150 [NASA ADS] [CrossRef] [Google Scholar]
 Arfken, G. B., & Weber, H. J. 2005, Mathematical Methods for Physicists (Elsevier) [Google Scholar]
 Asada, K., & Nakamura, M. 2012, ApJ, 745, L28 [NASA ADS] [CrossRef] [Google Scholar]
 Asada, K., Nakamura, M., Doi, A., Nagai, H., & Inoue, M. 2014, ApJ, 781, L2 [Google Scholar]
 Asada, K., Nakamura, M., & Pu, H.Y. 2016, ApJ, 833, 56 [NASA ADS] [CrossRef] [Google Scholar]
 Baczko, A. K., Schulz, R., Kadler, M., et al. 2019, A&A, 623, A27 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Blakeslee, J. P., Jordán, A., Mei, S., et al. 2009, ApJ, 694, 556 [Google Scholar]
 Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433 [Google Scholar]
 Blandford, R., Meier, D., & Readhead, A. 2019, ARA&A, 57, 467 [NASA ADS] [CrossRef] [Google Scholar]
 Blasi, P. 2013, A&ARv, 21, 70 [Google Scholar]
 Boccardi, B. 2015, Ph.D. Thesis, Universität zu Köln, Germany [Google Scholar]
 Chatterjee, K., Liska, M., Tchekhovskoy, A., & Markoff, S. B. 2019, MNRAS, 490, 2200 [Google Scholar]
 Cheung, C. C., Harris, D. E., & Stawarz, Ł. 2007, ApJ, 663, L65 [Google Scholar]
 Cho, I., Jung, T., Zhao, G.Y., et al. 2017, PASJ, 69, 87 [Google Scholar]
 Cui, Y.Z., Hada, K., Kino, M., et al. 2021, Res. Astron. Astrophys., 21, 205 [Google Scholar]
 de Gasperin, F., Orrú, E., Murgia, M., et al. 2012, A&A, 547, A56 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Dodson, R., Edwards, P. G., & Hirabayashi, H. 2006, PASJ, 58, 243 [Google Scholar]
 Doeleman, S. S., Fish, V. L., Schenck, D. E., et al. 2012, Science, 338, 355 [NASA ADS] [CrossRef] [Google Scholar]
 EHT MWL Science Working Group (Algaba, J. C., et al.) 2021, ApJ, 911, L11 [NASA ADS] [CrossRef] [Google Scholar]
 Event Horizon Telescope Collaboration (Akiyama, K., et al.) 2019, ApJ, 875, L1 [Google Scholar]
 Event Horizon Telescope Collaboration (Akiyama, K., et al.) 2021, ApJ, 910, L12 [Google Scholar]
 Fromm, C. M., Ros, E., Perucho, M., et al. 2013, A&A, 557, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gebhardt, K., & Thomas, J. 2009, ApJ, 700, 1690 [NASA ADS] [CrossRef] [Google Scholar]
 Gebhardt, K., Adams, J., Richstone, D., et al. 2011, ApJ, 729, 119 [Google Scholar]
 Ghisellini, G. 2013, Radiative Processes in High Energy Astrophysics, 873 (Springer) [Google Scholar]
 Ginzburg, V. L., & Syrovatskii, S. I. 1964, The Origin of Cosmic Rays (New York, Macmillan) [Google Scholar]
 Greisen, E. W. 2003, in Information Handling in Astronomy  Historical Vistas, 285, 109 [NASA ADS] [CrossRef] [Google Scholar]
 Hada, K., Doi, A., Kino, M., et al. 2011, Nature, 477, 185 [NASA ADS] [CrossRef] [Google Scholar]
 Hada, K., Kino, M., Nagai, H., et al. 2012, ApJ, 760, 52 [NASA ADS] [CrossRef] [Google Scholar]
 Hada, K., Kino, M., Doi, A., et al. 2013, ApJ, 775, 70 [Google Scholar]
 Hada, K., Giroletti, M., Kino, M., et al. 2014, ApJ, 788, 165 [NASA ADS] [CrossRef] [Google Scholar]
 Hada, K., Kino, M., Doi, A., et al. 2016, ApJ, 817, 131 [Google Scholar]
 Hada, K., Park, J. H., Kino, M., et al. 2017, PASJ, 69, 71 [NASA ADS] [CrossRef] [Google Scholar]
 Haga, T., Doi, A., Murata, Y., et al. 2015, ApJ, 807, 15 [NASA ADS] [CrossRef] [Google Scholar]
 Hardcastle, M. J., & Krause, M. G. H. 2013, MNRAS, 430, 174 [NASA ADS] [CrossRef] [Google Scholar]
 Harris, D. E., Biretta, J. A., Junor, W., et al. 2003, ApJ, 586, L41 [NASA ADS] [CrossRef] [Google Scholar]
 Harris, D. E., Cheung, C. C., Stawarz, Ł., Biretta, J. A., & Perlman, E. S. 2009, ApJ, 699, 305 [NASA ADS] [CrossRef] [Google Scholar]
 Hawley, J. F., Fendt, C., Hardcastle, M., Nokhrina, E., & Tchekhovskoy, A. 2015, Space Sci. Rev., 191, 441 [NASA ADS] [CrossRef] [Google Scholar]
 Hovatta, T., Aller, M. F., Aller, H. D., et al. 2014, AJ, 147, 143 [Google Scholar]
 Jaffe, W. J., & Perola, G. C. 1973, A&A, 26, 423 [NASA ADS] [Google Scholar]
 Jiang, W., Shen, Z., MartíVidal, I., et al. 2021, ApJ, 922, L16 [CrossRef] [Google Scholar]
 Junor, W., Biretta, J. A., & Livio, M. 1999, Nature, 401, 891 [NASA ADS] [CrossRef] [Google Scholar]
 Kardashev, N. S. 1962, Sov. Ast., 6, 317 [NASA ADS] [Google Scholar]
 Kim, J.Y., & Trippe, S. 2014, J. Korean Astron. Soc., 47, 195 [Google Scholar]
 Kim, J.Y., Lee, S.S., Hodgson, J. A., et al. 2018a, A&A, 610, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kim, J. Y., Krichbaum, T. P., Lu, R. S., et al. 2018b, A&A, 616, A188 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kino, M., Takahara, F., Hada, K., & Doi, A. 2014, ApJ, 786, 5 [NASA ADS] [CrossRef] [Google Scholar]
 Kino, M., Takahara, F., Hada, K., et al. 2015, ApJ, 803, 30 [NASA ADS] [CrossRef] [Google Scholar]
 Kino, M., Takahashi, M., Kawashima, T., et al. 2022, ApJ, 939, 83 [NASA ADS] [CrossRef] [Google Scholar]
 Komissarov, S. S., Barkov, M. V., Vlahakis, N., & Königl, A. 2007, MNRAS, 380, 51 [Google Scholar]
 Kovalev, Y. Y., Lister, M. L., Homan, D. C., & Kellermann, K. I. 2007, ApJ, 668, L27 [NASA ADS] [CrossRef] [Google Scholar]
 Kovalev, Y. Y., Lobanov, A. P., Pushkarev, A. B., & Zensus, J. A. 2008, A&A, 483, 759 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kravchenko, E., Giroletti, M., Hada, K., et al. 2020, A&A, 637, L6 [EDP Sciences] [Google Scholar]
 Lee, S.S., Byun, D.Y., Oh, C. S., et al. 2015, J. Korean Astron. Soc., 48, 229 [NASA ADS] [CrossRef] [Google Scholar]
 Lisakov, M. M., Kovalev, Y. Y., Savolainen, T., Hovatta, T., & Kutkin, A. M. 2017, MNRAS, 468, 4478 [Google Scholar]
 Longair, M. S. 2011, High Energy Astrophysics (Cambridge, UK: Cambridge University Press) [Google Scholar]
 Ly, C., Walker, R. C., & Junor, W. 2007, ApJ, 660, 200 [NASA ADS] [CrossRef] [Google Scholar]
 Lyutikov, M., Pariev, V. I., & Gabuzda, D. C. 2005, MNRAS, 360, 869 [Google Scholar]
 Macchetto, F., Marconi, A., Axon, D. J., et al. 1997, ApJ, 489, 579 [NASA ADS] [CrossRef] [Google Scholar]
 MAGIC Collaboration (Acciari, V. A., et al.) 2020, MNRAS, 492, 5354 [NASA ADS] [CrossRef] [Google Scholar]
 Marscher, A. P. 2010, in Jets in Active Galactic Nuclei, ed. T. Belloni (Berlin, Heidelberg, Springer) 794, 173 [NASA ADS] [Google Scholar]
 McKinney, J. C. 2006, MNRAS, 368, 1561 [NASA ADS] [CrossRef] [Google Scholar]
 Mertens, F., Lobanov, A. P., Walker, R. C., & Hardee, P. E. 2016, A&A, 595, A54 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Müller, C., Kadler, M., Ojha, R., et al. 2011, A&A, 530, L11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Nakamura, M., Asada, K., Hada, K., et al. 2018, ApJ, 868, 146 [Google Scholar]
 Narayan, R., Chael, A., Chatterjee, K., Ricarte, A., & Curd, B. 2022, MNRAS, 511, 3795 [NASA ADS] [CrossRef] [Google Scholar]
 Niinuma, K., Lee, S.S., Kino, M., et al. 2014, PASJ, 66, 103 [Google Scholar]
 Ostrowski, M. 1998, A&A, 335, 134 [NASA ADS] [Google Scholar]
 O’Sullivan, S. P., & Gabuzda, D. C. 2009, MNRAS, 400, 26 [Google Scholar]
 Pacholczyk, A. G. 1970, Radio Aastrophysics. Nonthermal Processes in galactic and Extragalactic Sources (San Francisco: Freeman) [Google Scholar]
 Park, J., Hada, K., Kino, M., et al. 2019a, ApJ, 887, 147 [Google Scholar]
 Park, J., Hada, K., Kino, M., et al. 2019b, ApJ, 871, 257 [Google Scholar]
 Park, J., Hada, K., Nakamura, M., et al. 2021, ApJ, 909, 76 [Google Scholar]
 Perlman, E. S., & Wilson, A. S. 2005, ApJ, 627, 140 [NASA ADS] [CrossRef] [Google Scholar]
 Plavin, A. V., Kovalev, Y. Y., Pushkarev, A. B., & Lobanov, A. P. 2019, MNRAS, 485, 1822 [Google Scholar]
 Pushkarev, A. B., & Kovalev, Y. Y. 2012, A&A, 544, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pushkarev, A. B., Hovatta, T., Kovalev, Y. Y., et al. 2012, A&A, 545, A113 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pushkarev, A. B., Butuzova, M. S., Kovalev, Y. Y., & Hovatta, T. 2019, MNRAS, 482, 2336 [NASA ADS] [CrossRef] [Google Scholar]
 Reynolds, C. S., Fabian, A. C., Celotti, A., & Rees, M. J. 1996, MNRAS, 283, 873 [NASA ADS] [Google Scholar]
 Rieger, F., & Levinson, A. 2018, Galaxies, 6, 116 [NASA ADS] [CrossRef] [Google Scholar]
 Rybicki, G. B., & Lightman, A. P. 1979, Radiative Processes in Astrophysics (New York: Wiley) [Google Scholar]
 Shepherd, M. C., Pearson, T. J., & Taylor, G. B. 1994, BAAS, 26, 987 [NASA ADS] [Google Scholar]
 Sironi, L., Rowan, M. E., & Narayan, R. 2021, ApJ, 907, L44 [CrossRef] [Google Scholar]
 Tavecchio, F., & Ghisellini, G. 2008, MNRAS, 385, L98 [Google Scholar]
 Thompson, A. R., Moran, J. M., Swenson, G. W. J. 2017, Interferometryand Synthesis in Radio Astronomy, 3rd Edition (Springer) [CrossRef] [Google Scholar]
 Walker, R. C., Hardee, P. E., Davies, F. B., Ly, C., & Junor, W. 2018, ApJ, 855, 128 [Google Scholar]
 Walsh, J. L., Barth, A. J., Ho, L. C., & Sarzi, M. 2013, ApJ, 770, 86 [NASA ADS] [CrossRef] [Google Scholar]
 Zamaninasab, M., ClausenBrown, E., Savolainen, T., & Tchekhovskoy, A. 2014, Nature, 510, 126 [Google Scholar]
 Zavala, R. T., & Taylor, G. B. 2003, ApJ, 589, 126 [NASA ADS] [CrossRef] [Google Scholar]
 Zdziarski, A. A., Sikora, M., Pjanka, P., & Tchekhovskoy, A. 2015, MNRAS, 451, 927 [CrossRef] [Google Scholar]
Appendix A: Effect of the (u, v)range matching on the spectral index distribution
As described in Sect. 3.1, it is important to match the (u, v) ranges of the two frequencies when producing spectral index maps. Considering the typical (u, v) ranges for KaVA and VLBA at 22 and 43 GHz given in Table A.1, we used data within 33−170 Mλ for KaVA and 25−685 Mλ for VLBA. This cuts out the data on average ∼ 21 % (∼ 46 %) for KaVA 22 GHz (43 GHz), and ∼ 5 % (∼ 8 %) for VLBA 22 GHz (43 GHz), respectively. As a result, the image noise increases on average by ∼ 7 % (∼ 33 %) for KaVA 22 GHz (43 GHz), and ∼ 0 % (3 %) for VLBA 22 and 43 GHz, respectively.
Typical (u, v) range and the angular size corresponding to the shortest baseline of KaVA and VLBA at 22 and 43 GHz.
The upper panel of Fig. A.1 compares the spectral index distributions before (red lines) and after (blue lines) matching the (u, v) ranges. This shows the distributions using the full (u, v) range have steeper spectra at longer distances. The difference between these two distributions – the artificial steepening along the jet distance (Δα_{arti}) – is shown in the lower panel. The gray lines are from individual epochs and the thick black line is their average. We also calculated the expected Δα_{arti} under the assumption that the loss of sensitivity for largescale structures at the shortest baseline is ∼sin(πθ/θ_{max})/(πθ/θ_{max}), (e.g., Thompson et al. 2017), where θ is the size of the structure and θ_{max} is the angular size corresponding to the shortest baseline. Here, we assume θ as jet widths and adopt it from Hada et al. (2013), and θ_{max} of KaVA and VLBA at 22 and 43 GHz are summarized in Table A.1. After calculating the expected Δα_{arti} for KaVA and VLBA, they were weighted and averaged by the number of observations to calculate the final value. This value is displayed in a red line in the lower panel.
Fig. A.1. Effect of the (u, v)range matching on the spectral index distribution. The upper panel shows spectral index distributions between 22 and 43 GHz before (red lines) and after (blue lines) matching the (u, v) range. Distributions using the full (u, v) range have steeper spectra at longer distances. Differences between the two distributions (i.e., artificial steeping due to inconsistent (u, v) ranges) are indicated by gray lines in the lower panel. The thick black line is the average difference, and the red line is the expected artificial steepening along the jet. See text for details. 
The averaged Δα_{arti} increases with distance, close to the expected value, and is ∼ − 0.2 at ∼ 10 mas. Individual cases show Δα_{arti} ∼ −0.6 at a distance of ∼ 7 mas in severe cases. Even taking this into account, the observed spectral index decreases more significantly (Δα ∼ −2). Therefore, we conclude that the decrease in the spectral index with the distance of the M 87 jet is intrinsic.
Appendix B: Effect of beam size on the spectral index distribution
Typical beam sizes of VLBA 22 and 43 GHz are 0.8 × 0.4 mas and 0.45 × 0.22 mas, respectively, which are much smaller than the typical beam sizes of KaVA. In Fig. B.1 the spectral index distribution between the VLBA 22 and 43 GHz maps was convolved into two different beams: 0.8 × 0.8 mas and 1.2 × 1.2 mas. The radius of each beam is equal to the major axis of the VLBA and KaVA beams at 22 GHz, respectively. Compared to the spectral index distribution using the KaVA beam, the distribution using the VLBA beam has a greater scatter with some local fluctuations. As the beam size increases, the local structures are averaged and the overall distribution is smoother. However, it can be seen that the global trend of the spectrum does not change with the size of the beam.
Fig. B.1. Spectral index distributions between VLBA 22 and 43 GHz maps convolved with beams of different sizes. The blue lines are convolved with a circular beam of 1.2 × 1.2 mas, where the radius is the same as the major axis of the KaVA 22 GHz beam. The red lines are convolved with 0.8 × 0.8 mas, which is equal to the major axis of the VLBA 22 GHz beam. The distribution with the VLBA beam has a larger scattering than the distribution with the KaVA beam. 
Appendix C: Alternative explanation of the spectrum steepening
In our model, we assumed a relationship of α = (p + 1)/2 between the slope of the electron distribution and the synchrotron spectrum (Rybicki & Lightman 1979). This assumption is maintained when the electron energy corresponding to the observation frequency is within a sufficiently wide range of γ in the optically thin region. The observation frequencies 22 and 43 GHz correspond to γ ∼ 10^{1.6 − 2.0} and are located within the electron energy range of the injection function (γ_{inj, min} = 10^{1}, γ_{inj, max} = 10^{5}). However, without electron injection, γ_{max} decreases over time and thus the above assumption may not hole. In this case, the synchrotron spectrum has a different shape depending on the synchrotron emission model we use.
The most widely used synchrotron emission models were proposed in Kardashev (1962) and Pacholczyk (1970; KP model) and Jaffe & Perola (1973; JP model). In the JP model, which is used in our spectral index model, the scattering due to turbulence in the magnetic field randomizes the electron pitch angles on shorter timescales relative to the radiative lifetime (e.g., Hardcastle & Krause 2013). Consequently, the cutoff energy is determined only by magnetic field strength, not by the pitch angle.
In the KP model, however, the pitch angle of individual electrons remains constant over the lifetime, and electrons with a pitch angle parallel to the magnetic field do not lose energy. Hence, the cutoff energy depends not only on the magnetic field strength but also on the pitch angle. Assuming an isotropic pitch angle, the synchrotron spectrum does not drop exponentially but has a break with the spectral index ranging from α_{inj} at low frequencies to 4/3α_{inj} − 1 at high frequencies (Kardashev 1962)^{9}.
Figure C.1 shows a schematic picture of the spectral index distribution using the KP and JP models in the absence of electron injection. Applying an initial slope of α_{inj} = −0.7, the KP type energy loss predicts the slope at energies higher than the break energy as 4/3α_{inj} − 1 ≈ −1.933, which is within the allowed region defined from observations. Therefore, this scenario may also explain the spectral index steepening.
Fig. C.1. Schematic diagram of the spectral index distributions without electron injection, assuming different synchrotron emission models. The blue line assumes the JP model, which is used in this study. The red line shows the spectral index distribution with the KP model. Assuming that the injection spectral index is α_{inj} = −0.7, the expected index after steepening is 4/3α_{inj} − 1 ≈ −1.933, which is within the allowable range obtained from observation. 
In summary, we suggest two scenarios to explain the observed spectral index distribution. (1) JPtype energy loss and particle injection with distance dependence. (2) KPtype energy loss and no electron injection. Both scenarios make an electron energy distribution a shape of a broken power law, and hence it is possible to explain the observed trend of the spectral index distributions. Further study is needed to discuss which scenario is more plausible on the M 87 jet. One of the expectations is that the brightness distribution between the JP model scenario and the KP model scenario is different since the JP model scenario requires continuous electron injection and the KP model scenario does not. Another distinguishing point between the JP and the KP scenarios is the fractional polarization. The JP model needs a fully random field plus an isotropic pitch angle. This will lead to very low fractional polarization. The KP model predicts linear polarization (depending on the magnetic field geometry). Therefore, fractional polarization could be a good indicator of the synchrotron loss model.
Appendix D: Details of the spectral index model
The general transfer equation of the electron number density distribution is written as follows (e.g., Ginzburg & Syrovatskii 1964; Longair 2011; Blasi 2013):
where D is the diffusion coefficient, v is the velocity of the system, b is the energy loss rate, and Q is the nonthermal electron injection rate. In our discussion, the frame of reference is a cross section of the jet comoving with the jet’s flow, then v ⋅ ∇N = 0. We further neglect the diffusion of the electron. Then, Eq. (D.1) becomes Eq. (4).
To solve the general transfer equation numerically, we change the partial derivative equation into a series of ordinary differential equations (e.g., Arfken & Weber 2005). First, we differentiate b(γ), (Eq. 5) with respect to γ,
Putting Eq. (D.2) into Eq. (D.1) and using the relation , we get
The lefthand side of Eq. (D.3) is the same as a total derivative of N(γ, τ) with respect to τ. Consequently, the transfer equation changes into a series of ordinary differential equations:
which are solvable numerically. The number density of electrons N(γ) changes with time with respect to Eq. (C.4.1), and at the same time, each electron in N(γ) loses its energy by following Eq. (C.4.2), (same as Eq. 5).
To solve the problem, we additionally need to know the relationship between the distance of the reference frame from the SMBH and the time of the frame (i.e., the proper time). This is because the physical quantities in the equation were obtained as a function of distance, while the transfer equation deals with the change with time. When τ is the proper time and z is the distance, τ and z are in the following relation:
where t_{obs} is the time in the observer’s frame, is the bulk jet speed in units of the speed of the light and Γ is the bulk Lorentz factor. Here, we use the equation of apparent jet speed and Doppler factor (e.g., Ghisellini 2013). The relation between τ and z is then calculated by integrating Eq. (C.5):
The distance range of our calculation is from z_{i} = 2 mas to z_{f} = 10 mas, which correspond to z_{i} ≈ 889 r_{s} and z_{f} ≈ 4446 r_{s} in deprojected units assuming a viewing angle of 17° (Walker et al. 2018). Applying the observed bulk jet velocity field of the M 87 jet (Park et al. 2019a) and setting τ_{i} = 0s, the travel time is τ_{f} ∼ 2.2 × 10^{8}s (∼ 6.9 years). If we apply the velocity field with faster acceleration in the region closer than 2 mas (Mertens et al. 2016), the travel time between 2 − 10 mas is ∼ 2.4 years.
Now we solve the Eqs. (C.4.1) and (C.4.2) to obtain the electron number density distribution as a function of time, N(γ, τ) ^{10}, which is further converted to N(γ, z) using Eq. (C.6). From N(γ, z), the slope of the number density distribution between γ(ν_{obs} = 22 GHz) and γ(ν_{obs} = 43 GHz) along the jet distance, p_{22 − 43GHz}(z), is obtained by using the relation of the observed frequency of synchrotron radiation and the energy of nonthermal electrons (Rybicki & Lightman 1979):
Finally, the spectral index distribution α_{22 − 43GHz}(z) is obtained by using the relation α = (p + 1)/2 (Rybicki & Lightman 1979).
Appendix E: Notes on the injection function Q(γ, z)
The electron injection function in our calculation is assumed to be a powerlaw distribution with a sharp cutoff at γ_{inj, min} = 10 and γ_{inj, max} = 10^{5}. The minimum and maximum range in the electron energy distributions are constrained by the observations. The M 87 jet has been detected in VLBI observations at frequencies as low as 1.6 GHz (e.g., Cheung et al. 2007), and as low as ∼20 MHz using LowFrequency Array (LOFAR; e.g., de Gasperin et al. 2012). From Eq. (C.7), the electron energy corresponding to synchrotron radiation at 20 MHz is γ(ν_{obs} = 20MHz)∼10^{0.1 − 0.4} when B_{i} = 0.3 − 1.0G and Γ = 2 is assumed. Likewise, the detection of Xray emission in the M 87 jet between the core and HST1 could constrain the maximum energy of the electron distribution (e.g., Perlman & Wilson 2005; EHT MWL Science Working Group 2021). From Eq. (C.7), the Xray synchrotron emission requires γ(ν_{obs} = 300PHz)∼10^{5.2 − 5.5}. Perlman & Wilson (2005) suggested that even higher energies of nonthermal electrons of γ ∼ 10^{6 − 8} is required in order to emit Xray synchrotron emission. We also found that electrons with energies greater than 10^{5} have little effect on the spectral index distribution at 22−43 GHz. Figure E.1 shows the spectral index distribution model using different γ_{inj, max} in the injection function. For all models, q = 10 and B_{i} = 0.3G were applied. As γ_{inj, max} increases, the spectral index distributions gradually converge, and when γ_{inj, max} is greater than 10^{4.5}, the spectral index distribution does not change. This is because the synchrotron cooling time of the electrons with energy γ ≳ 10^{5} is , which is much shorter than the travel time (∼2.2 × 10^{8}s). From this experiment, we found that, for γ_{inj, max} > 10^{5}, the spectral index distribution remains unchanged regardless of the shape of the highenergy tail of the injection function (e.g., a sharp cutoff or an exponentially decreasing). Therefore, in this study, we assume the injection function as a powerlaw energy distribution and the sharp cutoff at γ_{inj, max} = 10^{5}.
Fig. E.1. Comparison of model spectral index distributions using nonthermal electron injection functions with different γ_{inj, max}. Distributions with different colors represent different values of γ_{inj, max}. For all distributions, the magnetic field strength was assumed to be B_{i} = 0.3G. As γ_{inj, max} increases, the spectral index distributions converge and there is little change between the models with γ_{inj, max} = 10^{4.5} and γ_{inj, max} = 10^{6}. 
Movie
Movie 1 associated with Fig. 3 Access here
All Tables
Typical (u, v) range and the angular size corresponding to the shortest baseline of KaVA and VLBA at 22 and 43 GHz.
All Figures
Fig. 1. Spectral index maps between 22 and 43 GHz obtained from VLBA and KaVA observations. All images have been rotated by −18°. The restoring beam size is 1.2 mas × 1.2 mas, drawn as a black circle in the bottomleft corner. Observing dates are shown to the left of each map. The contours represent the total intensity at 22 GHz. Contours start at 3σ_{rms}, increasing in steps of 2. 

In the text 
Fig. 2. Evolution of the spectral index between 22 and 43 GHz along the deprojected distance from the SMBH. Weighted mean values across the jet are used in this distribution (see Sect. 3.2). The projected distance is displayed on the upper axis in mas, and the corresponding deprojected distance is displayed on the lower axis in r_{s}. Red lines represent the spectral index distribution from 2016 KaVA observations and blue lines the spectral index distribution from VLBA observations from 2010 and 2014. The shaded area is the 1σ error of the spectral index. 

In the text 
Fig. 3. Results of the modeled electron energy distributions (N(γ)) and the spectral index distributions as a function of distance along the jet. This is done to test the effect of the nonthermal electron injection rate in the jet. Left: Modeled electron number density distributions with different values of q = 0, 5, 10, and ∞. The different colors represent different values of q in the model: q = 0 represents continuous new electron injections in the jet, q = ∞ no new electron injections in the jet, and q = 5 and 10 new electrons injected into the jet but the volume decreases with distance. In all cases, we set the same magnetic field distribution (B_{i} = 0.3 G). The slope of the nonthermal electron injection is p_{inj} = −2.4. The vertical red and blue lines represent the energies of the nonthermal electrons, which emit at 22 and 43 GHz. The location down the jet where the calculations were made is written in the top right corner. Right: Modeled spectral index distributions as a function of distance, which correspond to the slope of the electron number density distributions between the vertical lines in the left panel. The vertical dashed black line is the location down the jet where the calculations were made. The horizontal dashed line is the slope of the injection function (α_{inj} = −0.7). An animation of this figure is available. 

In the text 
Fig. 4. Constraining the magnetic field strength and nonthermal electron injection of the M 87 jet. Left: Modeled spectral index distributions constrained by the observations. The distributions must lie in the gray shaded allowed region (see Sect. 4.4 for details). Colored lines are the modeled spectral index distributions. Different colors represent different initial magnetic field strengths, B_{i}. The dashed and dotted red line is a modeled spectrum index distribution with q = 10 and q = 11, respectively, with B_{i} = 0.3 G. Right: Parameter space (B_{i}, q). The black region is the allowed values of B_{i} and q constrained by the observations. 

In the text 
Fig. 5. Magnetic field distribution of the M 87 jet. The black area is the distribution of the magnetic field strength estimated in this study. The solid lines with the gray shaded region are the extrapolation of the upper and lower limits of the result in the upstream direction (B(z)∝z^{−0.72}). The dashed lines are the magnetic field distribution which assumes a steeper slope of the jet radius profile in the region closer to the SMBH (B(z)∝z^{−0.92}). The solid lines without a shaded region are the upper and lower bounds of the magnetic field distribution, assuming fast acceleration in the region closer within 2 mas (B(z)∝z^{−1.12}). The dotted line is the extrapolation of the maximum values of the result to the downstream region. Markers are estimated magnetic field strength from previous studies. The blue cross is the theoretically expected magnetic field strength in a BlandfordZnajek jet (Blandford et al. 2019). The circles are magnetic field strength values obtained in VLBI cores derived using SSA theory (Kino et al. 2014, 2015; Hada et al. 2012, 2016; Reynolds et al. 1996; EHT MWL Science Working Group 2021). The diamonds are values obtained from VLBI observations using various methods other than SSA theory (Acciari et al. 2009; Kim et al. 2018b; Zamaninasab et al. 2014; Event Horizon Telescope Collaboration 2019, 2021; Jiang et al. 2021). The squares are the magnetic field strengths at HST1 (Harris et al. 2003, 2009). The horizontal dotteddashed lines are the magnetic field strengths estimated by fitting the broadband SED of the jet using a singlezone SSC emission model (Abdo et al. 2009; MAGIC Collaboration 2020; EHT MWL Science Working Group 2021). 

In the text 
Fig. A.1. Effect of the (u, v)range matching on the spectral index distribution. The upper panel shows spectral index distributions between 22 and 43 GHz before (red lines) and after (blue lines) matching the (u, v) range. Distributions using the full (u, v) range have steeper spectra at longer distances. Differences between the two distributions (i.e., artificial steeping due to inconsistent (u, v) ranges) are indicated by gray lines in the lower panel. The thick black line is the average difference, and the red line is the expected artificial steepening along the jet. See text for details. 

In the text 
Fig. B.1. Spectral index distributions between VLBA 22 and 43 GHz maps convolved with beams of different sizes. The blue lines are convolved with a circular beam of 1.2 × 1.2 mas, where the radius is the same as the major axis of the KaVA 22 GHz beam. The red lines are convolved with 0.8 × 0.8 mas, which is equal to the major axis of the VLBA 22 GHz beam. The distribution with the VLBA beam has a larger scattering than the distribution with the KaVA beam. 

In the text 
Fig. C.1. Schematic diagram of the spectral index distributions without electron injection, assuming different synchrotron emission models. The blue line assumes the JP model, which is used in this study. The red line shows the spectral index distribution with the KP model. Assuming that the injection spectral index is α_{inj} = −0.7, the expected index after steepening is 4/3α_{inj} − 1 ≈ −1.933, which is within the allowable range obtained from observation. 

In the text 
Fig. E.1. Comparison of model spectral index distributions using nonthermal electron injection functions with different γ_{inj, max}. Distributions with different colors represent different values of γ_{inj, max}. For all distributions, the magnetic field strength was assumed to be B_{i} = 0.3G. As γ_{inj, max} increases, the spectral index distributions converge and there is little change between the models with γ_{inj, max} = 10^{4.5} and γ_{inj, max} = 10^{6}. 

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