© ESO 2022

1. Introduction

Parsec-scale properties of jets in radio-loud active galactic nuclei (AGN) have been studied in much detail using the radio-astronomical technique of very long baseline interferometry (VLBI). In particular, the Monitoring Of Jets in Active galactic nuclei with VLBA Experiments (MOJAVE) program (Lister et al. 2009) has been highly successful in exploring the jet kinematics of sources at a frequency of 15 GHz (Lister et al. 2016). These data have been used to study the acceleration and collimation of parsec-scale jets (Homan et al. 2015), their opening angles (Pushkarev et al. 2017b), polarization properties (Lister et al. 2018; Pushkarev et al. 2017a), spectral distributions (Hovatta et al. 2014), and other applications.

By employing the model of a freely expanding relativistic jet (Blandford & Königl 1979; Königl 1981), the underlying physical properties such as the electron density, magnetic field strength, or power of the jets can be inferred from such observations. The strong observational biases due to relativistic bulk motion, such as aberation and Doppler boosting, can be corrected for by using statistical methods (Cohen et al. 2007). Beyond scaling models, VLBI observations can also be compared with detailed relativistic magnetohydrodynamical simulations (Porth et al. 2011; Fendt 2006), which describe the launching of jets from the black-hole accretion disk (Blandford & Payne 1982) and ergosphere (Blandford & Znajek 1977), as well as their acceleration and collimation. The simulations describe the fluid-dynamical evolution of the Poynting-flux dominated jets on subparsec scales all the way to the scales where they become visible on (multi) parsec scales as typical kinetic-energy-dominated jets to VLBI and connected radio interferometers. On VLBI scales, shock waves are commonly believed to accelerate particles to ultrarelativistic energies by the Fermi mechanism, causing moving knots in the jet flow and multiwavelength flares (e.g., Marscher et al. 2008).

In addition to brightness temperature (T_B) observations of the VLBI cores of AGN jets (Homan et al. 2021), the measurement of brightness-temperature gradients along their jets is a powerful independent observational tool. This approach can probe the internal and external physical conditions in AGN jets and has been employed in several single-source jet studies. For 3C 111, Kadler et al. (2008) measured brightness-temperature gradients and their evolution to study the particle and magnetic field density associated with individual jet components during flares. In the case of NGC 1052 (a twosided source), Kadler et al. (2004) and Baczko et al. (2019) studied the T_B gradients along both jets. Kadler et al. (2004) find a break of the T_B-gradients near the cores which they attribute to possible free-free absorption associated with a dusty torus. Baczko et al. (2019) further find that single power-law fits do not sufficiently describe the measured T_B behavior and adopt a two-zone model with a break in the brightness-temperature and diameter gradient of the jet. Studies of the brightness temperatures of VLBI jets and their gradients with distance from the core are still rare for larger samples.

In order to characterize jets close to their collimation zone, Pushkarev & Kovalev (2012) carried out dual-frequency observations at 2 GHz and 8 GHz, respectively, using a sample of 30 AGN jets. They found that at 8 GHz the jet diameter gradient shows flatter slopes than expected for freely expanding jets, while at 2 GHz the jet geometry is consistent with a conical geometry. They found brightness temperature gradients of T_B ∝ r^−2.2. Furthermore, Pushkarev & Kovalev (2012) found indications for the presence of a transition zone in blazar jets where the jet flow switches from a confined, accelerating jet (parabolic geometry) to a freely expanding jet (conical geometry). In the case of M87, Asada & Nakamura (2012) showed that such a geometry transition can be observed in the vicinity of the Bondi-sphere (Russell et al. 2015) where the external medium is expected to change its pressure gradient. More recent works also found these transitions in several different sources, NGC 6251 (Tseng et al. 2016), NGC 4261 (Nakahara et al. 2018), NGC 1052 (Nakahara et al. 2020; Baczko et al. 2022), NGC 315 (Boccardi et al. 2021; Park et al. 2021) and 1H0323+342 (Hada et al. 2018). In the cases of M87 (Mertens et al. 2016; Hada et al. 2018; Park et al. 2019), NGC 315 (Park et al. 2021) and 1H0323+342 (Hada et al. 2018) the jets also show accelerating jet components near the geometry-transition zones. Kovalev et al. (2020) systematically searched for these geometry transitions in 331 MOJAVE sources at 15 GHz and 1.4 GHz. They demonstrated the presence of geometry transitions from a parabolic to a conical geometry in 10 nearby sources on scales of several thousand Schwarzschild radii R_S. Potter & Cotter (2015) applied parameters, obtained from 42 simultaneous blazar spectral energy distributions (SED) to a 1d-relativistic fluid-dynamical jet model to predict a region where the jets experience a transition between a magnetically-dominated region (closer to the central engine), where the jet geometry is parabolic and a particle-dominated region of the jet (farther downstream) where the jet geometry is conical. Potter & Cotter (2015) predict a correlation between the radii of the jets in the transition zone and the jet power as well as the maximum bulk Lorentz factors. Consequently, they also predict a dichotomy of the blazar classes, BL Lacertae objects (BL Lacs) and flat-spectrum radio quasars (FSRQs), regarding the jet radii in the transition zone, where BL Lacs are expected to show smaller jet radii than FSRQs. Kovalev et al. (2020) argue that external pressure gradients most likely are not solely responsible for the geometry transitions within the VLBI jets and that the mechanism might involve a transition of the jets from a Poynting-flux to a particle-dominated configuration.

In order to investigate the T_B-gradients and morphological changes of parsec-scale jets at higher frequencies, we conducted a systematic dual-high-frequency (15 GHz and 43 GHz) analysis of a sample of blazars using data obtained with the Very Long Baseline Array (VLBA) data as part of the MOJAVE¹ (15 GHz) and Boston University² (BU, 43 GHz) jet monitoring programs. The following cosmological constants were adopted throughout the paper: h = 0.71, Ω_Λ = 0.73,  Ω_M = 0.27.

2. Sample

The objects and their observation epochs (i.e., characterizing the time interval during which the observation was carried out) for this study were selected by the following criteria:

• (a)

  The sources are in the MOJAVE data archive (Lister et al. 2018) and the BU sample (Jorstad & Marscher 2016).

• (b)

  Data exists for epochs between 2003 and 2013.

• (c)

  The imaged jet brightness distribution can be sufficiently represented and parametrized with at least three 2D-Gaussian jet components in at least 5 epochs³.

We chose the earliest epoch considered in our study to be at the onset of the new MOJAVE program (Lister & Homan 2005) at which the overall image fidelity and typical image noise levels improved systematically over earlier data sets from the progenitor 2cm-Survey program (Kellermann et al. 1998). The resulting sample consists of 15 flat-spectrum radio quasars (FSRQs), 11 BL Lacertae objecs (BL Lacs) and two radio galaxies (RGs), in total 28 objects. The object names, their classification, redshifts. and black hole masses are listed in Table 1.

3. Method

3.1. Parameter fitting

The optically thin synchrotron emission regions in the jet were fitted with circular Gaussians using DIFMAP⁴ (Shepherd 1997). Initially, the visibility map is opened and a starting model consisting of a circular Gaussian component in the center is fitted to best match the data. If the residual map still shows nonzero intensity after subtracting the established model, further Gaussians are added and fitted until the residual map shows only (white) noise. As an example, the corresponding jet model resulting from a single epoch of 4C +71.07 (0836+710) observations is shown in Fig. 1 for both frequencies. Our model fits are optimized to represent accurately the local surface brightness and jet width along the jet rather than individual moving knots within the jet. For this reason, our model fits differ from the ones used in standard MOJAVE kinematics studies. The diameter d, flux F = I × Ω_G (I is the intensity and Ω_G the solid angle of the fit Gaussian component), and the (angular) distance from the core r are measured for each component. Unresolved components are filtered according to

$$\begin{array}{c}\hfill r_{\lim }=b_{\mathrm{\Psi }}\sqrt{\frac{4\ln (2)}{\pi }\ln \left(\frac{S/N}{S/N-1}\right)}\end{array}$$(1)

Fig. 1. Fitted Gaussian-component representation of 4C +71.07 (0836+710) at 15 GHz (top) and at 43 GHz (bottom). The contours for both frequencies are shown at a 5σ level and are separated by a factor of 2. In the 15 GHz image the map noise level is ∼0.2 mJy. The contour levels are [−0.85,0.85,1.7,3.4,6.7, 13.7,27.5,54.9,109.8,219.6,430.2,878.5] mJy). At 43 GHz the map noise level is ∼0.8 mJy. The contour levels are [−3.8,3.8,7.5, 15.1,30.2,60.3,120.6,241.2,284.5,482.5,964.9, 1929.8, 3859.6] mJy.

which describes the minimal radius a resolved fitted component can achieve for a given signal-to-noise ratio (S/N) and beam size (b_Ψ) (cf. Kovalev et al. 2005). The uncertainties of the parameters of each component are approximated with the post-fit RMS σ_rms information (Fomalont 1999). The uncertainty of the flux density is further altered by also geometrically adding a systematic value of 5% of the flux density to account for absolute flux-scale uncertainties:

$$\begin{array}{cc}& {\sigma }_F=\sqrt{{\sigma }_{\mathrm{rms}}^2+{\sigma }_{\mathrm{rms}}^{2}S/N+{\left(0.05F\right)}^2}\hfill \end{array}$$(2)

$$\begin{array}{cc}& {\sigma }_d=d\frac{{\sigma }_F}{F}.\hfill \end{array}$$(3)

The diameters of the remaining components from each epoch are plotted together into a single graph against the distance from the radio core, at a 1-σ level (Fig. A.1 ff.). The radio core itself is the unresolved component closest to the (0/0)-position in the map. All distances are measured against the position of the respective core component. This ensures that the set of model components characterizing the jet is representative of the long-term time-averaged intensity of the evolving jet. By contrast, individual model components at a given epoch are representative of the local jet brightness and width at a given time. The components in these sources typically have considerable proper motions of up to several milliarcseconds along the jet per year (e.g., Lister et al. 2016), evolving in brightness and in width during this time. Our long-term multi-epoch data base allows us to parameterize the jet geometry and brightness temperature for each source in our sample at two separate high radio frequencies and at submilliarcsecond resolution. To this end, we fit the jet diameter with d ∝ r^l, that is as a power-law of distance r with index l. The brightness temperature is also fitted as a power law T_B ∝ r^s, characterized by the index s. This procedure is performed for all sources at both frequencies, 15 GHz and 43 GHz, along the respective entire jet length and in the overlap regions of both frequencies. To determine a physical overlap region, the core shift (Pushkarev et al. 2012) has to be taken into account since the turnover in the synchrotron spectrum from optically thick to optically thin emission depends on the frequency, that is, the radio core appears closer to the jet apex at higher frequencies. The core shift Δr_cs usually has to be calculated on a source-by-source basis between two frequencies. At lower frequencies (∼8 − 15 GHz), Pushkarev et al. (2012) studied the nuclear opacity of radio jets on a larger sample of 191 sources. For higher frequencies, O’Sullivan & Gabuzda (2009) analyzed the core shift between ∼15 and 43 GHz for three sources, 1418+546, 2007+777 and 2200+420. The core shifts are calculated assuming conical jet geometries and are measured to be between 0.04 ≲ Δr_cs/mas ≲ 0.1, which corresponds to 5 − 10% of the respective beam size of the used data. Typical beam sizes for the data used in this work are ∼0.5 mas for 15 GHz and ∼0.3 mas for 43 GHz, and thus by at least a factor of 3 larger than reported core shifts between these frequencies. The effect of such small core shifts in this study are negligible, however, in the case of 2200+420 we tested the impact of this core shift on the results which will be discussed in the corresponding paragraph in Sect. 6.

3.2. Scaling model for the jets

For the parameterization of the jet diameter and brightness temperature, we employ the relativististic jet model introduced by Blandford & Königl (1979) and Königl (1981). The basic assumptions of the model are:

• (a)

  The jet has a semi-opening angle Φ and an inclination angle ϑ to the line of sight. The diameter of the jet scales along the jet-axis r, d_jet ∝ r^l. The corresponding jet morphologies relate to a free expansion with constant advance speed (l = 1), accelerated motion (l < 1) or jet deceleration (l > 1).

• (b)

  Electrons are injected continuously into the jet and the electron energy density (n_e) decreases as a power law along the jet axis according to n_e ∝ rⁿ.

• (c)

  The magnetic field scales with a power law along the jet axis according to B ∝ r^b. A purely toroidal (azimuthal) magnetic field scales with the jet transverse radius R as B ∝ R⁻¹, resulting in a scaling along the jet axis of B ∝ r^−l. A purely poloidal magnetic field (parallel to the jet axis) would scale with the jet radius according to B ∝ R⁻² which translates to scaling along the jet-axis according to B ∝ r^−2 l.

The observed flux F_ν of a circular component with a diameter d_jet and radiating area $A\propto d_{\text{jet}}^2$ is proportional to

$$\begin{array}{c}\hfill F_{\nu }\propto A·d_{\mathrm{jet}}·j_{\nu }\end{array}$$(4)

where $j_{\nu }\propto n_e·B·{\left(\frac{\nu }{{\nu }_B}\right)}^{\alpha },({\nu }_B\propto B),$ is the emissivity of the synchrotron emission (cf. Krolik 1999, Sect. 9.2.1). The properties of the jet geometry, magnetic field, and the electron energy density (a)–(c), Eq. (4) implies that the flux F_ν of a model component has a power law dependency, scaling along the jet axis

$$\begin{array}{c}\hfill F_{\nu }\propto r^{3l+n+b(1-\alpha )}.\end{array}$$(5)

The brightness temperature T_B of a Gaussian-shaped component is given by

$$\begin{array}{c}\hfill T_{\mathrm{B}}=\frac{2\ln (2)}{\pi k}\frac{F_{\nu }{\lambda }^2(1+z)}{R^2}\end{array}$$(6)

where k denotes the Boltzmann constant, λ the observed wavelength, and r = d_jet/2 the radius of the component⁵. Equation (6) implies that T_B ∝ F_ν/R² (R² ∝ r^2l) and therefore

$$\begin{array}{c}\hfill T_{\mathrm{B}}\propto r^{l+n+b(1-\alpha )}.\end{array}$$(7)

For convenience, we define the index s that describes the scaling of T_B ∝ r^s by

$$\begin{array}{c}\hfill s=l+n+b(1-\alpha ).\end{array}$$(8)

Now, the complete set of scaling indices is given by the shape index l, the density index n, the magnetic field index b, and the brightness-temperature index s. Their interdependence depends on the degree of fulfillment of adiabatic expansion, particle number conservation, and magnetic-field equipartition.

3.2.1. Comparing the modelfit method with the stacked-image analysis

We compared the method of estimating the jet widths at a given position, described above with the stacked-image analysis as employed by Pushkarev et al. (2017b) and Kovalev et al. (2020). It is not a-priori clear that both methods should yield the same results, since it is known that stacked images map out a broad profile while at any given time traveling jet features do not fill out this entire stacked-profile range (Pushkarev et al. 2017b; Lister et al. 2021). On small scales, a varying jet nozzle can lead to strong variations of component position angles in the inner jet, especially in high-resolution data at high frequencies, which biases the apparent jet geometry in stacked images. In individual sources, further differences can arise when a source exhibits a strongly resolved wide jet emission region far downstream from the core. The modelfit method will typically yield one single large component to represent such an extended region in each individual epoch while image stacking and image-plane fitting of stacked images will benefit from enhanced signal-to-noise ratio and cover a large area. In result, a diameter gradient fit with the modelfit method will put much stronger weight on the (typically many) components from the compact jet and little weight on the single extended Gaussian component representing the diffuse emission. The extended jet regions in the image-stacking method are also weighted less because the method is sensitive to the signal to noise ratio (e.g., Kovalev et al. 2020; Casadio et al. 2021). However, since the jet is split into equidistant sections along the jet axis, the larger area of extended emission will yield more data points to be fit in these regions as the modelfit method. For such cases the modelfit method can derive flatter diameter gradients than the stacked image-plane method. We emphasize that a modelfit component in a single epoch at a given position in some cases will not fill the physical jet width. Therefore it is important to not derive jet parameters from single epoch fits. Therefore when studýing the gradients all epochs are considered at once. In consideration of the arguments made above, we consider the modelfit method in general as the most robust way to characterize the jet geometry.

A further difference between our analysis and the approach of Kovalev et al. (2020) is that we are not attempting to fit an additional parameter r₀ to consider the separation of the VLBI jet core from the true jet base. Median value of this parameter determined on large samples are of the order of 0.1 mas (e.g., Kravchenko et al., in prep.) but the generalized fitting of both r₀ and l does not converge in case of all of our data sets. We therefore choose to neglect this offset, understanding that this tends to lead to somewhat flattened apparent jet gradients. We tested the strength of this effect in case of 18 sources for which a generalized fitting was possible and found differences in |l| in the range of 1.4 × 10⁻³ to 0.27. In case of sources with well-determined and realistically small r₀ values, we found that both methods, model fitting and image stacking, yield comparable results and are overall in good agreement.

3.2.2. Changes in the geometrical shape

To measure possible changes in the geometrical shape of the jets, we subdivided each jet along the distance from the core into 10⁴ slices. The center of each slice represents a possible break point x_b in the jet. One of these break points is picked randomly, then two power-law fits are performed left and right of x_b, respectively. If the fit parameters are plausible, namely ${\chi }_{\text{red}}^2\sim 1$, the errors of the inner and outer jet geometry fits (l_in, l_out) are of the same order of magnitude and also if l_in < l_out, the solution is added to the sum of solutions. This is done 5000 times. To then find the best (l_in, l_out, x_b), the mean value of the respective parameter is calculated over all runs and weighted with the respective ${\chi }_{\text{red}}^2$ value. This set of parameters (l_in, l_out, x_b) is used as initial values for a smooth broken-power-law fit of the form

$$\begin{array}{c}\hfill f(x)=A{2}^{(l_{\mathrm{in}}-l_{\mathrm{out}})/\mathrm{\Delta }}{\left(\frac{x}{x_b}\right)}^{l_{\mathrm{in}}}{[1+{\left(\frac{x}{x_b}\right)}^{\mathrm{\Delta }}]}^{(l_{\mathrm{out}}-l_{\mathrm{in}})\mathrm{\Delta }}\end{array}$$(9)

(Nakahara et al. 2020; Astropy Collaboration 2013, 2018), where A = d_j(x_b) is the jet diameter at the break point x_b, l₁ is the diameter gradient index for the inner jet (closer to the radio core) and l₂ is the diameter gradient of the outer jet. This procedure turned out to be necessary because the fit of the smoothly broken power law is sensitive to properly guessed initial values. In order to not having to find them by visual inspection, we performed the search algorithm described above. There are different ways this function can be fit to data. Baczko et al. (2022) deploy the orthogonal distance regression while Casadio et al. (2021) use a gnuplot implementation of the Levenberg-Marquardt algorithm. In this work we use scipy (Virtanen et al. 2020) in order to fit the smoothly broken power law with a classical least square method. We used the trust-region reflective algorithm rather than the Levenberg-Marquardt algorithm as in Boccardi et al. (2021). The problem at hand is a bound problem in the sense that we let the fit algorithm vary the initial conditions. These initial conditions are obtained from the numerical splitting of the jet. A comparison between the least square fit results and the orthogonal distance regression can be found in Appendix B. We tested the parameter range of Δ between 10 (Baczko et al. 2022) and 100 (Boccardi et al. 2021). Within the errors the results are not distinguishable, therefore we set the parameter constant to 100.

4. Results

The results of the diameter and brightness-temperature gradient fits are shown in Figs. A.1–A.28 and summarized in Table 2. Results were obtained for the full jet length at each frequency and for the overlap regions probed at both frequencies.

A statistical comparison was carried out to find whether the scaling parameters show systematic differences between the FSRQ and BL Lac classes. The two RGs in the sample, 3C 111 and 3C 120, show a behavior similar to the FSRQs with a pronounced complex substructure consisting of bright radio knots. They will be discussed in more detail in Sect. 6. The distributions of the scaling parameters are shown in Figs. 2–5. All histograms are shown with a standard binning and a binning method based on Bayesian block analysis (Scargle et al. 2013), implemented in astropy (Astropy Collaboration 2013, 2018). The Bayesian binning method allows us to pinpoint significant differences between histograms at a glance. For all plots the indices (Q, B, G) denote FSRQs, BL Lacs and radio galaxies, respectively. The index “ol” denotes values that are calculated in the overlap region in contrast to the entire jet length. These were introduced in order to make the often densely packed plots more readable. In cases of uneven sample sizes, only the Kolmogorov-Smirnov (KS) test was performed. For samples of equal size, viz. the comparison of both frequencies, also the Spearman, Kendall and Pearson correlation coefficients were calculated. In cases where a difference between two classes are indicated by the KS-test, we checked the robustness of the difference by repeatedly performing the KS test while randomly excluding sources from each class (15 FSRQs, 11 BL Lacs) to find out which and how many sources are driving the possible differences. This was done by randomly excluding up to two sources from each sample. The results of the KS-test analysis is shown in Figs. 6 and 7.

Fig. 2. Distribution of l-values along the entire jet (top panel) and along the overlap-region (bottom). Values at 15 GHz are shown in blue and 43 GHz in red. The histograms are shown for equal bin width of Δl = 0.2 (dashed lines) and for Bayesian binning (solid lines; see text).

Fig. 3. Distributions of the geometry scaling parameter l (left) and its variance χ2 (right) shown with Bayesian (solid lines) and uniform binning (dashed lines) for FSRQs (red) BL Lacs (bue), and RGs (green). Upper four panels: distributions for the entire jet at both frequencies; Lower four panels: distributions for the overlap-region at both frequencies. Bin sizes were chosen as 0.2 for the l distributions and 1.0 for the variance distributions.

Fig. 4. Distribution of power-law slopes s of the brightness-temperature gradient along the entire jet (top panel) and along the overlap-region (bottom) at 15 GHz (blue) and 43 GHz (red). The histograms are shown for uniform bin width Δs = 0.5 (dashed lines) and for Bayesian binning (solid lines; see text).

Fig. 5. Bayesian binning (solid lines) and uniform binning (dashed lines) of the TB power-law slope distributions. The colors indicate the source classes (red: FSRQs, blue: BL Lacs, green: RGs). Top two rows: distributions along the entire jet for both frequencies; Bottom two lines: distributions along the overlap-region of both frequencies. Left column: distribution of s-values. For both frequencies, the bin size is 0.5. Right column: ${\chi }_{\text{red}}^2$ value distribution describing the brightness temperature. For both frequencies, a unit bin size was chosen.

Fig. 6. P-values for each pair of parameters involving the jet geometry compared with the two-sided KS test. The fields are color coded ranging from p < 5% (dark blue) to p ≥ 10% (yellow). The lower triangle, where all l values are compared and the upper triangle where all ${\chi }_{\text{red}}^2$-values are compared, in the plot have to be read separately. The indices (Q, B, G) denote FSRQs, BL Lacs and radio galaxies, respectively. The index “ol” denotes values that are calculated in the overlap region in contrast to the entire jet length. These were introduced in order to make the plots more readable.

Fig. 7. P-values for each pair of parameters involving the TB gradients compared with the two-sided KS test. The fields are color coded ranging from p < 5% (dark blue) to p ≥ 10% (yellow). The lower triangle, where all s values are compared and the upper triangle where all ${\chi }_{\text{red}}^2$-values are compared, in the plot have to be read separately.

4.1. Diameter gradients

For the diameter gradients the results of the statistical comparisons between the frequencies are listed in Table 3 and the results of the comparison of the source classes are shown in Fig. 6. We find typical values of l in the range between 0.5 and 1, i.e., collimated jet geometries. Uncertainties range between ∼0.08 in case of the apparently most uniform jet geometries to ∼0.3 for some complex cases where deviations from a simple r^l power-law behavior are most pronounced. The most common l value at 43 GHz is around 0.5 with typically somewhat larger values at 15 GHz.

Comparison of both frequencies. Figure 2 shows the l-value distribution at both frequencies along the entire jet detected at the respective frequency (top) and the overlap region (bottom). Along the entire jets there is a hint that the jets at 43 GHz are more collimated than the jets at 15 GHz (KS: p ∼ 2.6%). This tentative difference vanishes when comparing the jets in the overlap region (KS: p ∼ 94.4%).

Comparison of the source classes at 15 GHz. Figure 3 shows the l-value distribution for the 15 GHz jets. The measurements along the entire jet (first row) indicate a difference in jet geometry between the FSRQs and the BL Lac objects (KS: p ∼ 0.4%). The FSRQs show smaller l-values and thus stronger collimated jets than the BL Lac objects. The two radio galaxies 3C 111 and 3C 120 show collimated jets in the same regime as the FSRQs. As before in the case of the 15 GHz to 43 GHz comparison, however, when considering the overlap regions (third row), the differences between FSRQs and BL Lacs in terms of l become insignificant (KS: p ∼ 62.1%). When studying the distribution of the ${\chi }_{\text{red}}^2$ values for the jet-geometry power-law fits (Fig. 3, top row, right column) we find a pronounced difference between the FSRQs and the BL Lacs (KS: p ∼ 0.2%). A simple power law fit is a systematically worse description for the FSRQs than for the BL Lacs, regarding the ${\chi }_{\text{red}}^2$ values. This difference does not disappear when comparing the overlap regions in the jet (third row right column, KS: p ∼ 3.2%).

Comparison of the source classes at 43 GHz. At 43 GHz the l-values measured along the entire jet (Fig. 3, second row, left column) yield no significant distinction between FSRQs and BL Lacs (KS: p ∼ 11.9%) nor between FSRQs and RGs (KS: p ∼ 41.2%). Also the comparison of BL Lacs and RGs yield no hint for a difference between the source classes (KS: p ∼ 7.7%). In the overlap region, however, BL Lacs behave differently than FSRQs (KS: p ∼ 2.7%), showing a tendency for larger l-values. Regarding the ${\chi }_{\text{red}}^2$-distributions, no difference between FSRQs and BL Lacs can be constituted along the entire jet length and the overlap region, see the right column in Fig 2, second and fourth row.

4.2. Brightness-temperature gradients

For the power-law brightness-temperature gradient index s the results of the statistical comparisons between the frequencies are listed in Table 3 and the results of the comparison of the source classes are shown in Fig. 7. Most sources show T_b power law indices in the range of −3 ≲ s ≲ −1.5 with the width of the distribution somewhat increasing when only the overlap region at both frequencies is considered. Uncertainties of individual s-values range from 0.2 to 0.7, reflected in the large Bayesian-block widths. Most commonly, the slopes along the full jet widths lie in the range between −2 and −1.5.

Comparison of both frequencies. Figure 4 shows the s-value distribution of both frequencies evaluated along the entire detected jet at the respective frequency (top) and along the overlap region (bottom). The s-values are indistinguishable for both frequencies in the entire jet and overlap regions, respectively.

Comparison of the source classes at 15 GHz. Figure 5 (top row) shows the s-value distribution at 15 GHz along the full jet length and the corresponding overlap region, where the colors indicate the source classes. Regarding the distribution of the brightness temperature power law indices, the different source classes are indistinguishable, both along the entire jet length and the overlap region. The ${\chi }_{\text{red},{\text{s}}_{\text{15}}}^2$ distribution (Fig. 5, top row), however, indicates a difference between FSRQs and BL Lacs along the entire jet length. The BL Lac objects show a tendency for smaller ${\chi }_{\text{red},{\text{s}}_{\text{15}}}^2$ values than FSRQs, indicating that a single power law is a better description for the brightness-temperature gradient for BL Lacs than for FSRQs (KS: p ∼ 0.2%). With reference to the overlap region, this difference also becomes apparent (KS: p ∼ 0.1%).

Comparison of the source classes at 43 GHz. Figure 5 shows the s-value distribution at 43 GHz along the full jet length (second row) and the corresponding overlap region (fourth row). The s-value distributions are statistically indistinguishable both along the entire jet and the overlap region. As in the case of the jet-geometry power-law fits, the ${\chi }_{\text{red},{\text{s}}_{\text{43}}}^2$ value distributions (Fig. 5) indicate smaller values for BL Lac object than FSRQs along the entire jet lengths and the overlap regions (KS: p ∼ 1.4% in both cases).

5. Discussion

In this section, we discuss the general implications of the gradient fits at 15 GHz and 43 GHz that we have obtained in the previous section for our sample of 28 blazars and radio galaxies. A discussion of our findings regarding the individual sources will follow in the next section.

5.1. Brightness-temperature gradients

Pushkarev & Kovalev (2012) analyzed the T_B-gradients for 30 sources at 2 GHz and 8 GHz and found a mean value of ⟨s⟩∼ − 2.2. Furthermore Kravchenko et al (in prep., priv. comm.) find a similar value analyzing the entire MOJAVE data archive at 15 GHz. This is consistent with our smaller dual-frequency sample results (cf. Sect. 4.2). The T_B-gradients are determined by the jet geometry (see below), particle and magnetic-field gradients. We conclude that studies at different frequencies in general lead to consistent results for these physical quantities. Specifically, for both our frequencies, the s distributions are not clearly distinguishable along the entire jets and the overlap regions. Moreover, the different source classes of BL Lacs and FSRQs do not show any strong differences in terms of the average T_B-gradient trends. However, the ${\chi }_{\text{red}}^2$-distribution at both frequencies show that a single power law fit is more appropriate for BL Lacs than for FSRQs, indicating that the latter show a more complex T_B-distribution along the jet axis due to local surface brightness enhancements, leading to more scatter in our global fits. This is consistent with the findings of Jorstad et al. (2017) who find high T_B-values to be more common in knots of FSRQs as compared to BL Lacs. This gives support to the idea of systematically higher Doppler factors and stronger shocks in FSRQ jets and/or possibly more rapid expansion in the jets of BL Lacs.

5.2. Jet geometry

Measured along the entire jet lengths there is a systematic difference in collimation between the jets measured at 15 GHz and 43 GHz. The jets at 43 GHz appear more collimated (tend to show smaller l-values) than at 15 GHz. In the overlap region, this difference is not significant. It is therefore plausible, that the differences arise from the geometrical properties of the jet at larger scales, i.e. that the l-values change with r. Indeed, Kovalev et al. (2020) showed that in many cases, it is possible to measure a transition from parabolic to conical jet geometries.

At 15 GHz, the source classes show a notable difference when comparing them along the entire jet. FSRQs systematically show lower l-values than BL Lacs indicating a stronger jet collimation. This is also found by Kravchenko et al. (in prep., priv. comm.) on a larger 15 GHz AGN-jet sample. In our study, this difference is again not significantly found in the (smaller) overlap region of the images at both frequencies indicating again that the differences arise further down the jet and reinforcing the idea of a geometry transition on parsec scales. Analyzing the ${\chi }_{\text{red}}^2$ of the single power-law fits to the jet geometry also yields interesting results. The FSRQs tend to show larger ${\chi }_{\text{red}}^2$ values, both along the respective entire jet and in the overlap regions compared to BL Lacs. This indicates that the geometric structure in FSRQs is more complicated than the smooth power law observed in BL Lacs. The complexity could arise from regions where the full jets recollimate and widen again locally. Alternatively, some of the model components in FSRQ jets might indeed be representative of local brightness enhancements that are smaller than the full jet width at this distance. An example where this seems to be the case can be seen in Fig. 1 (top panel) at the relative coordinates (−6, −7).

At 43 GHz the difference between the source classes is not significantly detected along the entire jet nor in the overlap region. This is true for both, the l-values distribution and the ${\chi }_{\text{red}}^2$ – distribution. This could be simply a result of reduced image quality at the higher frequency. Alternatively, the high-frequency emission could be associated with plasma closer to the inner jet axis (the spine of the jet) while the lower-frequency emission might be dominated by outer layers (sheath). In such a stratified-jet scenario, the difference between FSRQ and BL Lac jet geometry that we find at 15 GHz can be interpreted as a process that is largely associated with the sheath while the spine remains unaffected.

5.3. Measurement constraints on jet-internal gradients

Figure 8 illustrates the allowed parameter space for the jet-internal gradients of the electron density n, the magnetic field strength b, and the spectral index α constrained by the measured values of the power-law slopes of the brightness temperature s and jet diameter l.

Fig. 8. Parameter space of n, b and α, spanned by the measured s and l values at both frequencies.

To investigate these dependencies in detail, we first reconsider some theoretical standard scenarios (cf. Kadler et al. 2004) in order to estimate plausible s-values, depending on l, n, b and α. We assume equipartition between the electron energy density and magnetic field strength density, thus n_e  ∝ B² and therefore the scaling parameter of the electron energy density n depends on the scaling of the magnetic field along the jet axis (n = 2b).

• Conical jet (l = 1):

  • torodial magnetic field (b = −1), equipartition (n = −2), α = −0.5 → s = −2.5.

  • polodial magnetic field (b = −2), equipartition (n = −4), α = −0.5 → s = −6.

In the case of a conical jet, equipartition, and a reasonable spectral index of α = −0.5 the s-values thus have a domain of reasonable values of −6 ≤ s ≤ −2.5.

• Parabolic jet (l = 0.5):

  • torodial magnetic field (b = −0.5), equipartition (n = −1), α = −0.5 → s = −1.25.

  • polodial magnetic field (b = −1), equipartition (n = −2), α = −0.5 → s = −3.

In the case of a parabolic jet, equipartition, and a reasonable spectral index of α = −0.5 the s-values have a domain of reasonable values of −3 ≤ s ≤ −1.25.

The measured s and l values for both frequencies seen in Fig. 8 show an anticorrelation (r_spearman ∼ −0.73, r_pearson ∼ −0.69) in the sense that larger l values lead to more negative s values. The r_pearson even indicates a simple linear connection of both parameters. This can also be constituted for the frequencies separately (r_{spearman, 15}, r_{pearson, 15}, r_{spearman, 43}, r_{pearson, 43})∼(−0.73, −0.69, −0.64, −0.69). Further splitting the values into source classes, ignoring the radio galaxies by number within the sample also keeps the anti-correlation measurable (r_{spearman, 15, Q}, r_{pearson, 15, Q}, r_{spearman, 43, Q}, r_{pearson, 43, Q})∼(−0.6,−0.62, −0.56, −0.66) and (r_{spearman, 15, B}, r_{pearson, 15, B}, r_{spearman, 43, B},r_{pearson, 43, B})∼(−0.56, −0.58, −0.86, −0.79). The anti-correlation of the s and l values suggests that more collimated jets imply flatter brightness-temperature gradients. It shows that the main impact of a resulting s-value comes from the jet geometry because other wise an uncorrelated data cloud would be expected. However the anti-correlation is not as clear, as the s values seem not to depend solely on l. The rank correlation indices show ranges of −0.8 ≲ r_s/p ≲ −0.5, which hints to a possible impact of the parameters n and b to the resulting s-values. This will be discussed source by source in the next section. The data are fit with a linear model (gray dashed lines) to measure the slope of the anti correlation which is ≃3.4. Equation (7) however suggests that this slope should be −1 if n, b and α where completely independent of the jet geometry. However as seen in the standard scenarios above and Sect. 3.2, at least the scaling of the magnetic field strength depends on the jet geometry. Depending on how n and b are de-coupled, the scaling of the electron energy density also depends on the jet geometry. The black dashed lines indicate the parameter space of the parameters (n, b, α)_min = ( − 2.5, −1.25, −1) and (n, b, α)_max = ( − 1, −0.5, 0). This can be concluded to be the parameter space of jet-internal l, b, and n gradients that is in agreement with the observational constraints for our sample analysis of VLBI jets.

5.4. Tests for redshift biases

The redshifts of the sources in our sample range from 0.0308 (Mrk 421) to 2.218 (4C +71.07). Moreover, the angular-to-linear scale conversion also strongly depends on the orientation of the jets to the line of sight. Thus, the angular scales probed with the VLBA in different sources can correspond to very different linear scales. To test whether our conclusions are affected by these effects and whether the overlap regions are comparable for the different sources, we first compute the apparent length of the jet from its maximum observed angular distance from the radio core. Second, we deproject the jet length from the plane of the sky using the inclination angle θ between the jet direction and the line of sight given in Pushkarev et al. (2009) to obtain its radial length in parsecs. For sources for which Pushkarev et al. (2009) does not list an inclination angle, we searched the literature in order to obtain an inclination angle, which is discussed on a source by source basis in Sect. 6. Assuming that the distance between the radio core and the jet base is small compared to the jet length, we find the overlap region typically at jet lengths 5 ≲ d_c/pc ≲ 100. Along the entire jet length and the overlap region, respectively, the brightness temperature and diameter gradients for each jet at 15 GHz were fitted. We tested the interesting pairs from Table 3 and Figs. 6, 7, regarding the jet collimation and complexity. The results from these tests are consistent with the findings from the projected mas-scale jets. Along the entire jet length, FSRQs are more collimated than BL Lac objects (KS: p ∼ 1.7%). Within the overlap region this difference is not significant (KS: p ∼ 28.7%). The interpretation is the same as before. The jets close to the base are similar while a possible difference in collimation can be measured further down stream. Regarding the ${\chi }_{\text{red}}^2$ distribution for the l₁₅ fits, FSRQs show a more complex structure along the entire jet length than BL Lac objects (KS: p ∼ 0.44%) which is again not found significant in the overlap region (KS: p ∼ 25.9%). Regarding the ${\chi }_{\text{red}}^2$ distribution from the s₁₅-fits, FSRQs consistently show larger ${\chi }_{\text{red}}^2$-values than BL Lacs, both along the entire jet length (KS: p ∼ 2.4%) and the overlap region (KS: p ∼ 0.04%). We conclude that using the lengths of the deprojected jets to describe the scaling behavior rather than their angular diameters does not alter our findings from the previous sections.

5.5. Possible biases form geometry transition sources

With the mean value and the variance of the l and s-value and the ${\chi }_{\text{red}}^2$ distributions, we simulated normal distributions with different class sizes. We tested, how much each class can be reduced in size in order for the KS-test to yield reliable results and therefore find the differences between the known distributions. We find, that each class should contain at least 8–10 values in order for the KS test to reliably find differences in known distributions. Knowing this, we excluded all sources that show a geometry transition (3 BL Lacs, 5 FSRQs and one radio galaxy). This reduces the sample sizes to 8 BL Lacs, 10 FSRQs and one radio galaxy. With the reduced sample, the difference in jet collimation at 15 GHz between FSRQs and BL Lacs become less significant along the entire jet length (KS: p ∼ 7%), however the complexity of the jet is still larger for FSRQs than for BL Lacs, for both the diameter and brightness temperature gradient fit, when comparing the ${\chi }_{\text{red}}^2$-values of the source classes for the different fits (KS: p-values in the range of (0.1 − 0.5)%). We want to emphasize, that this has to be taken with caution, since the sample sizes are reduced to an extent to where it is not entirely clear whether the KS-test is even able to yield interpretable results. Overall these findings however do not contradict the found differences between the source classes.

5.6. Geometry transition within the VLBI jets

A number of VLBI jets are known to exhibit transitions from a parabolic (l ∼ 0.5) to a conical (l∼ 1) geometry, for example for M87 (Asada & Nakamura 2012), NGC 6251 (Tseng et al. 2016), NGC 4261 (Nakahara et al. 2018), NGC 1052 (Nakahara et al. 2020), NGC 315 (Boccardi et al. 2021; Park et al. 2021) and 1H0323+342 (Hada et al. 2018). M87 (Mertens et al. 2016; Hada et al. 2018; Park et al. 2019), NGC 315 (Park et al. 2021) and 1H0323+342 (Hada et al. 2018) also feature accelerating jet components near the geometry transition zones. By comparing 2 GHz and 8 GHz data, Pushkarev & Kovalev (2012) found hints for a transition in the diameter scaling in jets. Kovalev et al. (2020) showed that geometry transitions from a parabolic to a conical geometry can be observed on scales of several thousand Schwarzschild radii R_S in individual AGN jets.

According to Potter & Cotter (2015), the four-momentum of a relativistic jet is composed of (I) the internal, kinetic and rest mass energy, (II) the magnetic energy and (III) the synchrotron, inverse Compton, and adiabatic losses. Depending on which term is larger than the magnetic term, the jet is said to be particle dominated or magnetically dominated, respectively. Transitions from parabolic to conical geometry as observed in the cases listed above may be explained by the jet switching from a magnetically dominated to a kinetically dominated flow. In the following, we use our systematic modeling of the 28 sources in our sample to put constraints on possible geometry breaks seen with this dual-frequency approach in AGN jets.

Figures 9–17, show the geometry breaks inferred from the fit procedure described above along with estimated positions of the Bondi radius r_b for each source. In order to obtain r_b, the mass of the SMBH and gas temperature of the matter accreting onto the SMBH has to be known: (Bondi 1952; Russell et al. 2015).

$$\begin{array}{c}\hfill \frac{r_b}{\mathrm{pc}}=31\times {\left(\frac{k_{b}T}{\mathrm{keV}}\right)}^{-1}\times \left(\frac{M}{{10}^9{M}_{\odot }}\right).\end{array}$$(10)

Fig. 9. Jet diameter vs. core distance for 3C 66A (0219+428). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.3 mas (dashed line). The dotted lines indicate the estimated location of the Bondi radius.

Fig. 10. Jet diameter vs. core distance for CTA 26 (0336-019). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.1 mas (dashed line). The dotted line indicate the estimated location of the Bondi radius.

Fig. 11. Jet diameter vs. core distance for 3C 111 (0415+379). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 5.4 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.

Fig. 12. Jet diameter vs. core distance for PKS 0528+134 (0528+134). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.1 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.

Fig. 13. Jet diameter vs. core distance for 4C +71.07 (0836+710). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.7 mas (dashed line). The dotted line indicate the estimated location of the Bondi radius.

Fig. 14. Jet diameter vs. core distance for Mrk 421 (1101+384). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 1.9 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.

Fig. 15. Jet diameter vs. core distance for 4C 29.45 (1156+295). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.5 mas (dashed line). The dotted line indicate the estimated location of the Bondi radius.

Fig. 16. Jet diameter vs. core distance for 3C 279 (1253-055). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 0.6 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.

Fig. 17. Jet diameter vs. core distance for BL Lac (2200+420). Top: single transition zone model (rigorous fit). The shape transitions at xb ∼ 1.8 mas from parabolic to conical. Bottom: two-transition model where the break points are ‘fit by eye’ parametrizes the complex structure.

In the case of M87, the black-hole mass was determined by Gebhardt et al. (2011) (M = (6.6 ± 4.4)×10⁹ M_⊙) and Walsh et al. (2013) ($M={3.5}_{-0.7}^{+0.9}{10}^9{M}_{\odot }$). Later the Event Horizon Telescope Collaboration (2019) found the black-hole mass to be M = (6.5 ± 0.9) 10⁹ M_⊙ in accordance with Gebhardt et al. (2011). Russell et al. (2015) determined a gas temperature of T ∼ 0.9 keV in the center of M87 resulting in a Bondi radius of r_b ∼ 224 pc from X-ray observations. For blazars, the temperature of the coronal plasma near the SMBH must be estimated taking into account that the observed X-ray emission is dominated by the nonthermal emission from the jet. For a sample low-power radio galaxies, Balmaverde et al. (2008) measured gas temperatures between 0.4 ≲ E/keV ≲ 1.7. For high-power sources, typical gas temperatures near the SMBH are up to ∼7 keV, (e.g., Yaji et al. 2010; Gaspari et al. 2013). Since we did not measure the Bondi radius for the sample, discussed in this work, we allow the gas temperature to vary within the bounds 0.4 ≲ E/keV ≲ 7 for each source and compare the region of plausible Bondi radii with the geometry transition region. The respective jet-inclination (viewing) angle is taken into account by using the measured angles from Pushkarev et al. (2009) for the cases of CTA 26 (0336−019, θ ∼ 3.2 deg), PKS 0528+134 (0528+134, θ ∼ 1.7 deg), 4C +71.07 (0836+710, θ ∼ 3.2 deg), 4C +71.07 (1156+295, θ ∼ 2 deg) and BL Lacertae (2200+420, θ ∼ 7.5 deg). In the other cases, we took values of the inclination angles from the references cited in Sect. 6. The results of the broken geometry fits are listed in Table 4. In each case we compare the ${\chi }_{\text{red}}^2$ of smoothly broken power law fit with a single power law fit to both the 15 GHz and 43 GHz data. In all cases an improvement is achieved by using the smoothly broken power law fit to both frequencies. For each source we bracket the Bondi radius by the minimum/maximum black hole mass from Table 1, the minimum/maximum gas plausible temperature and the maximum/minimum inclination angle. For 3C 66A (0219+428), CTA 26 (0336−019), PKS 0528+134 (0528+134), 4C +71.07 (0836+710), Mrk 421 (1101+384), 4C +29.45 (1156+295), 3C 279 (1255−055) and BL Lac (2200+420) the transition zone lies beyond plausible Bondi radius ranges. In such cases Kovalev et al. (2020) argue that the geometry transition could mark the transition of a magnetically dominated to a particle-dominated jet. The transition zone of 3C 111 (0415+379) lies within the plausible Bondi radius ranges. In this case, the lowest gas temperatures has to be assumed in order for the Bond sphere to extent beyond the transition region. For CTA 26 (0336−019) the jet features a bright, stationary component in the brightness temperature gradient between d_c ∼ 1 − 2 mas during the ten year observation period, see Fig. A.2. This possibly indicates stationary shocks or jet curvature, where the jet in parts crosses the line of sight, creating the seemingly stationary features. The distribution of the break points shown in Fig. 18 (top) tentatively suggests relatively low values for BL Lacs as compared to a wider distribution ranging to about a factor of ten larger maximum values in the case of FSRQs. This is found to be significant by a two-sample KS test (p ∼ 3.7%). Furthermore, a KS-test with the H₀-hypothesis, that x_b, Q < x_b, B yields p ∼ 1.7% which suggests, that the break point of the geometry transition is further down the jet stream for FSRQs than for BL Lacs. This is also found on the entire MOJAVE data sample at 15 GHz by Kravchenko et al. (in prep., priv. comm.) A difference might indeed be expected according to Potter & Cotter (2015) based on the assumption that the geometry breaks result from the magnetic-to-kinetic transition. Figure 18 (bottom) shows the distribution of the jet radii at the transition zone (break point x_b). Again, we do not find a significant dichotomy between the two source classes (KS: p ∼ 14.3%). Again an alternative KS-test with H₀ : R_j(x_b, Q) < R_j(x_b, B) yields p ∼ 7%. This suggests that the jet radii at the break points of FSRQs and BL Lacs are likely drawn from the same distribution however FSRQs show a tendency for larger radii at the break point. It must be emphasized, however, that the sample contains only ten sources that do show a geometry break, of which only three are of the BL Lac type. Furthermore, we cautiously remark that different methods can lead to somewhat different results. For example, Boccardi et al. (2021) analyzed the jet of NGC315 in the visibility domain at different epochs and frequencies finding a transition from parabolic to conical jet geometry already at a distance of about half a parsec. Park et al. (2021) performing a multifrequency single-epoch study in the image plane, found the transition to occur at a distance about 10 times larger. In order to further shed light onto the issue of systematic jet geometry transitions, it would be interesting to study more sources from even higher (> 43 GHz) frequencies to lower frequencies (∼1 GHz) in order to study the collimation of a jet across a large scale and determine whether geometry breaks might happen on all scales or whether the transitions found are interesting physical points along the physical flow of the jet.

Fig. 18. FSRQs (in red), and BL Lacs (in blue). Top: distribution of the break points in linear scale and deprojected. Bottom: distribution of the jet radii at the respective geometry break point.

6. Comments on individual sources

0219+428 (3C 66A). The brightness-temperature along the jet of this source shows a local excess over the general power-law decrease between 2 ≲ d_core/mas ≲ 3, which is visible at both frequencies (Fig. A.1). This excess indicates a local deviation from the simple jet model (see Eq. (7)). The jet-diameter does decrease as a power law as described by the model but shows a dip at the same distance from the radio core where the T_B excess is observed. The dip is persistent over the epochs, indicated by the different colors, represented within the dip. This might hint at a stationary region where the jet recollimates and the local surface brightness temperature rises.

The s-values at both frequencies are indistinguishable within the errors, however the l-values differ along the respective entire jet length, while remaining indistinguishable in the overlap region. We find a transition of the geometry in this jet, which is shown in Fig. 9, featuring a parabolic jet shape up to d_c ∼ 2.3 mas that transitions to a conical shape jet downstream.

To estimate the angle to the line of sight, the variability Doppler factor can be used to estimate θ ∼ 9.8 deg as discussed in Hovatta et al. (2009). The scaled Bondi radius is shown in Fig. 9 projected for the derived angle, respective SMBH mass and the gas plausible gas temperature range (dotted lines). The break point x_b lies outside the Bondi-sphere, which is located close to the possible recollimation region.

0336−019 (CTA 26). The jet diameter along the jet axis in this source is characterized by strong scattering around the best-fit model (Fig. A.2). This behavior is typical for a FSRQ (see discussion above). When combining the information from both frequencies as described in Sect. 3.2.2, it is possible to constrain a geometry break point in this source (Fig. 10). The jet is best described through a parabolic shape out to d_c ∼ 2.1 mas, which opens up to a cone afterwards. The Bondi radius is smaller than the location of this break by two orders of magnitude and is therefore unlikely to be related in any way. The T_B profile shows an excess at d_c ∼ 1 − 2 mas which is persistent across the epochs (indicated by the different colored components in Fig. A.2), which might be associated with a standing shock at this location in the jet.

0415+379 (3C 111). This is one out of only two radio galaxies in our sample. Model-fitting results for this source are shown in Fig. A.3. The morphology and kinematics are rather complex and best described with a parabolic-jet model at both frequencies and relatively flat brightness temperature gradients −1.8 ≲ s ≲ −1.5. Kadler et al. (2004) fit a brightness temperature gradient of with a power law index of s ≃ 2.4 by fitting individual tracked components. The difference can be explained by the fitting of more epochs in this work and not excluding single components for the power law fit. A break in the jet geometry from parabolic to hyperbolic shape can be constrained at ∼5.4 mas from the radio core, see Fig. 11. We adopt the inclination angle of 19 deg from Kadler et al. (2008) to project the position of the Bondi radius, which may coincide with the transition region for the smallest plausible gas temperatures near the SMBH. Beuchert et al. (2018) fit the brightness temperature and the diameter along the jet axis for this source with two separate power laws in the region < 3 mas and > 3 mas from the radio core. They find a similar switching of the jet geometry consistently with our findings. Kovalev et al. (2020) also find a transition of geometry at x_b = (7.03 ± 0.5) mas. Boccardi et al. (2021) find a transition in this region. Taking the systematic errors at a 1-σ-limit into account, there is still a difference between the found break point in this study and the break point found in Kovalev et al. (2020) and Boccardi et al. (2021). We tested the impact of different established fitting methods (least square and orthogonal distance regression) in Appendix B and found that different methods yield the same results with the given data within the errors. We therefore attribute this difference between our result and the Kovalev et al. (2020) and Boccardi et al. (2021) results to the use of different data sets, especially the different epochs with different time coverage of the dynamically evolving sources and with different observational setups and image noise levels (cf. Sect. 2) and the different sets of model fits (cf. Sect. 3) parameterizing this highly complex source.

0430+052 (3C 120). The model fitting results for second radio galaxy in our sample are shown in Fig. A.4. We find a parabolic jet shape on all scales tested along the entire jet length. The s-values are flatter than expected from a freely expanding jet which is consistent with the narrowly confined jet geometry. Kovalev et al. (2020) and Boccardi et al. (2021) find a transition of the jet shape from parabolic to conical at ∼1.5 mas from the radio core. This region is also resolved by the data used in this work, however Kovalev et al. (2020) probe somewhat larger scales in their analysis of 15 GHz and 1.4 GHz data. Examining Fig. 1 in Kovalev et al. (2020) it becomes clear, that the geometry transition is seen primarily due to the use 1.4 GHz data. Overall our findings are consistent with Kovalev et al. (2020) in the sense that closer to the physical jet base they find a parabolic shape consistent with our results. Even when rigorously searching for a geometry transition, it is not possible to clearly constrain a significant break point on the scales covered by our higher-frequency data. This is due to the large scatter of the data in the region where Kovalev et al. (2020) and Boccardi et al. (2021) find the geometry transition.

0528+134 (PKS 0528+134). This FSRQ shows a prominent excess in brightness temperature between 0.3 ≲ d_core/mas ≲ 1 distance from the radio core at 43 GHz which is not seen in the 15 GHz data (Fig. A.5). Since this excess in the 43 GHz is seen in epochs with are close in time to each other (see color coding in Fig. A.5) this could be due to a radio flare in this source. This excess can be associated with a dip in the jet-diameter profile of the 43 GHz jet. This indicates a region where the jet locally contracts. The observations seem to probe the jet in the region where its geometry is still parabolic. The geometry fits at both frequencies are indistinguishable within their errors: l₄₃ = 0.45 ± 0.19,  l₁₅ = 0.79 ± 0.21. This also holds in the overlap region of both frequencies. The brightness temperature gradients however seem to be significantly flatter at 43 GHz than at 15 GHz, where the latter shows the canonical value, expected from a simple Blandford & Königl (1979)-jet (s₄₃ = −1.23 ± 0.43,  s₁₅ = −2.65 ± 0.46), when considering the respective entire jet lengths. The difference of the s-values can be attributed to a variation of one of the parameters n and/or b on the different angular scales probed. Due to the excess in T_b at 43 GHz and the subsequent steepening of the gradient’s fit, this difference vanishes in the overlap zone of both frequencies. At high frequencies above ∼30 GHz (Marscher 2016), the radio core is not necessarily the region where the synchrotron radiation becomes optically thin as in the Blandford & Königl (1979) jet model but might rather represent a standing recollimation shock which is backed by polarization measurements in the radio core of 1803+784 at 43 GHz (Cawthorne et al. 2013), which they link to conical recollimation shocks (Daly & Marscher 1988; Cawthorne & Cobb 1990; Cawthorne 2006). Jorstad et al. (2001) studied 42 γ-ray detected blazars at 22 GHz and 43 GHz of which 27 showed a sum of 45 nonmoving components besides the radio core which they also link to standing recollimation shocks if the component is close to the radio core and to jet bending accompanied by maximized Doppler boosting of the component is farther down to jet. We find a geometry break from parabolic to conical jet geometry in this source at x_b ∼ 2.1 mas. The radius of the Bondi sphere is closer to the radio core than x_b in this source. Combining 15 GHz and 1.4 GHz, Kovalev et al. (2020), do not find a geometry break in this source. We explain this by the additional use of 43 GHz data in this work which allow to resolve components closer to the physical jet base due to following reasons: With the single power law fit, Kovalev et al. (2020) find a diameter gradient power law index of ∼0.8 which is consistent with the value we find at 15 GHz (l₁₅ = 0.79 ± 0.21) along the entire jet length. In the overlap region between the 43 GHz and 15 GHz the jet shape is parabolic. This reveals that the jet features a geometry transition which can be constrained closer to the jet based, by combining the high frequency observations.

0716+714 (S5 0716+71). This source features a jet profile with 0.5 ≤ l ≤ 1 at both frequencies and a quasi-canonical s-value in the range of the simple Blandford & Königl (1979) jet model at both frequencies along the entire jet lengths and the overlap region (Fig. A.6). No strong excesses in the brightness-temperature or the jet-diameter profiles are found and no geometry break can be constrained in this source.

0735+178 (OI 158). The jet in this BL Lac object is slightly collimated at both frequencies, along the entire jet lengths and the overlap region. Combining 15 GHz and 1.4 GHz data Kovalev et al. (2020) also find this slight collimation in this source (l ∼ 0.93).

The s-values cannot be distinguished in the overlap region, however along the entire jet length the 15 GHz jet shows the canonical Blandford & Königl (1979)-expected value while the 43 GHz jet shows a flatter T_b-gradient (s₄₃ = −1.75 ± 0.36,  s₁₅ = −2.46 ± 0.37; see Fig. A.7). According to Eq. (7) this might be attributed to a variation in the magnetic field strength gradient between the frequencies..

0827+243 (OJ 248). At each of the two frequencies, this quasar jet features l and s-values which are indistinguishable within their errors, along the entire jet length and the overlap region (Fig. A.8). In comparison, the l₄₃-value is consistent with a parabolic jet structure while the l₁₅-value tends to be closer to a conical jet geometry. However, due to the large uncertainties and the overall consistency of the values with each other, a geometry transition zone cannot be constrained in this source. Combining 15 GHz and 1.4 GHz the jet geometry can be fit with a cone (Kovalev et al. 2020) which is consistent with the findings in this work in the 15 GHz jets.

0829+046 (OJ 049). This BL Lac object shows a slightly colliamted jet at 43 GHz (Fig. A.9). The l-values here are consistent with the 15 GHz results in the overlap region. Considering the entire jet length, the 15 GHz data shows a larger l value indicative of a conical jet on larger scales. However, a consistent geometry break cannot be constrained due to the large scatter of the data at 43 GHz. This scatter is largely caused by a T_b-excess accompanied by a recollimation of the jet which can be seen at 43 GHz between 0.2 ≲ d_c/mas ≲ 1, which is less pronounced at 15 GHz. This can be understood as a result of jet recollimation similar to the case of PKS 0528+134. Investigating the 43 GHz light curve of this source, provided on the website of thwe BU monitoring program⁶, reveals, that the source showed flaring activity within the time period, analyzed in this wortk between 2008 and 2009 and in the time period between 2012 and 2013. In these time periods the MOJAVE data base provides sparse data with about one epoch per year which might be the reason that this excess is not seen in the 15 GHz data.

0836+710 (4C +71.07). Model-fitting results for this distant (z = 2.218) quasar are shown in Fig. A.10. The jet is more collimated at 43 GHz than at 15 GHz if the full jet lengths in both data sets are considered. However, the l-value in the overlap region of the 15 GHz jet is smaller and consistent with the 43 GHz geometry, which indicates a more collimated jet on small scales. This indicates that there is a geometry transition from parabolic to conical within this jet. The s-values are consistent with the canonical value of the simple Blandford & Königl (1979) model except for the overlap region at 43 GHz. The fact that the s-values are consistent with the canonical value for a freely expanding jet while the jet itself features a parabolic shape at 43 GHz, suggests that either n or b are causing the flattening the s-values. Interestingly, Orienti et al. (2020) find limb-brightened polarized flux up to 1 mas from the core which indicates a complex jet structure and that the magnetic field gradient changes at this point.

A two-zone jet geometry model is shown in Fig. 13. The jet features a parabolic shape up to d_c ∼ 2.7 mas, transitioning to a conical shape further downstream. The projected Bondi radius along the jet axis is smaller than the location of the transitioning zone by almost two orders of magnitude. The change in the magnetic-field gradient, suggested by the results of Orienti et al. (2020) might indicate that the jet changes from magnetically dominated to particle dominated which in turn can explain the geometry transition (Kovalev et al. 2020). The complexity of this source is also well studied regarding growing Kelvin-Helmholtz instabilities within the jet at multiple scales (Perucho et al. 2012; Vega-García et al. 2019, 2020). They also show that the jet is likely particle dominated over large scales out to ∼100 mas in jet length.

0851+202 (OJ 287). This well known BL Lac object has been suggested as a likely candidate for a binary supermassive black hole system Sillanpaa et al. (1988), Lehto & Valtonen (1996). In this context, it is interesting to note that this source is among the objects with the absolute highest scatter in the jet-diameter fits of our whole sample (Fig. A.11). It features l and s-values which are indistinguishable within the uncertainties, along the entire jet length and the overlap region. On larger scales at 15 GHz, the average jet geometry appears quasi conical, which can be interpreted as a freely expanding jet, which is consistent with the findings at lower frequencies from Kovalev et al. (2020).

0954+658 (S4 0954+65). This BL Lac shows a jet geometry which in both frequencies is consistent with a conical jet along the entire jet length and the overlap region (Fig. A.12). Also the s-values are indistinguishable from each other, however they tend to be steeper than the canonical value of the simple Blandford & Königl (1979) model, especially at 43 GHz in the overlap region, and at 15 GHz. As the jet overall seems to feature a conical shape, this must be understood in the context of the parameters n and b. Indeed when allowing not only toroidal magnetic field components (Blandford & Königl 1979) but also poloidal magnetic field components (Königl 1981), even steeper T_b-gradients can be obtained.

1101+384 (Mrk 421). This is a nearby (z = 0.0308) high-peaked BL Lac object in which we are probing the smallest linear scales within our sample (0.61 pc mas⁻¹). The source features s-values which are not distinguishable within their uncertainties at a given frequency (Fig. A.13). However s₄₃ overall shows a clear tendency of being flatter than s₁₅. This is also resembled in that fact, that l₄₃ is consistent with a parabolic jet shape while l₁₅ is consistent with a conical jet shape. This hints to a possible geometry transition within this source which is indeed found by our algorithm (see Sect. 3.2.2). A two-zone jet geometry model is shown in Fig. 14 and features a parabolic jet shape out to d_c ∼ 1.9 mas where it transitions to a conical shape jet further downstream. On kpc-scales, the inclination angle of this source was estimated to be θ ∼ 19 deg, based on the jet to counter-jet flux density ratio (Giroletti et al. 2006). On pc-scales, inclination angles for this source were estimated based on the Doppler factors obtained from variability analysis (Hovatta et al. 2015) to be between 0.5 ≲ θ/deg ≲ 4, however they address a point made by Lico et al. (2012), that the maximum measured component velocity β_app, max = 0.218 ± 0.026 (Lister et al. 2019) is most likely not the flow velocity of the jet because the beaming properties predicted by the flux ratio of the jet and the counter jet require viewing angles close to zero. Lico et al. (2012) estimate the inclination angle to be 2 ≲ θ/deg ≲ 5 based on the argument that the jet shows limb brightening (Piner et al. 2010), an indication for a transverse velocity structure across the jet axis. For our analysis we adopt the angle range given by Lico et al. (2012). The resulting estimated location of the Bondi radius is well below the distance of the break.

1127−145 (PKS 1127−14). Model-fitting results for this quasar jet are shown in Fig. A.14. It features a parabolic jet shape at both frequencies. The s-values are flatter than expected for a freely expanding jet. The relatively large ${\chi }_{\text{red}}^2$ in the gradient fits, indicate a complex jet structure. Taking both frequencies together, a systematic geometry transition cannot be found in this source. The lower frequency data from Kovalev et al. (2020) indicate a steepening of the diameter gradient. Cmbining the 43 GHz, 15 GHz and 1.4 GHz data of this source might reveal a systematic geometry transition.

1156+295 (4C 29.45). This quasar features brightness temperature gradients which are consistent with a freely expanding jet (Fig. A.15). The diameter gradients, however, indicate a collimation of the jet at both frequencies. A possible jet geometry transition is indicated when comparing the overlap region in the 15 GHz jet with the entire jet length. The two-zone jet geometry model is shown in Fig. 15. The jet features a parabolic shape out to d_c ∼ 2.5 mas, transitioning to a conical shape further out. The projected Bondi radius along the jet axis is smaller than the location of the transitioning zone by almost two orders of magnitude. In this source no steady bright component can be seen in the T_b-gradient, which is consistent with the low scatter atypical for quasars (cf. Figs. 3 and 5). At lower frequencies Kovalev et al. (2020) do not find a geometry transition. The single power law fit they find is consistent with l_out ∼ 1.3 from our two zone model.

1219+285 (W Comae). This source features l and s-values which are indistinguishable within their uncertainties, along the entire jet length and the overlap region (Fig. A.16). No significant break is found in the jet geometry of this source. The jet geometry overall is consistent both with a conical and parabolic jet, due to the scattering of the data. The s-values are consistent with the canonical parameter values of the simple Blandford & Königl (1979) model. Allowing for the extensions of Königl (1981) and in combination with the conical jet shape, this can be interpreted as a freely expanding jet.

1222+216 (4C +21.35). This is a quasar that features T_b gradients which are consistent with the s-value expected in a canonical Blandford & Königl (1979) jet model, at both frequencies, as seen in Fig. A.17. In the overlap region, the data suggest a steeper T_B gradient, especially at 15 GHz. When applying the search, discussed in Sect. 5.6 it is not possible to constrain a break point between two different slopes for the geometry scaling. This might hint to a more complex, transverse jet structuring, where at 43 GHz in the overlap region, a more collimated substructure of the jet is measured than at 15 GHz.

1226+023 (3C 273). The T_b gradients in 3C 273 at both frequencies are flatter than expected in a freely expanding jet. Indeed, the diameter gradient at both frequencies is consistent with a parabolic jet shape (Fig. A.18), although the data show large scatter as it is typical for quasar jets. As it is apparent in Fig. 8, low l-values tend to be associated with flat s gradients. The complex structure of the 3C 273 jet has been analyzed in great detail by several authors. In particular, a double-helix structure out to kpc-scales was found (see Lobanov & Zensus 2001).

1253−055 (3C 279). Similarly to 3C 273, the s-values are consistent with a canonical Blandford & Königl (1979) jet model. However also in this source, the 43 GHz jet diameter gradient is consistent with a parabolic jet shape and is more collimated than the 15 GHz jet diameter gradient which is consistent with a conical jet, in both the overlap region and the respective entire jet lengths. Akiyama et al. (2018) find a transition from a parabolic to a conical jet shape however they also make use of VLBA data at 1.4 GHz and MERLIN data at 1.6 GHz, which probes the jet further downstream out to kpc scales. They constrain a transition zone which in Fig. A.19 would appear at ∼1 mas. Our higher-frequency VLBA data are more sensitive to smaller scales and we find a transition zone at x_b ∼ 0.6 mas. In the outer fit component, we still find an l-value smaller than unity and a possible second break at or above 1 mas is not probed by our data. The radius of the Bondi sphere for this source is well below the location of both the geometry break point reported by Akiyama et al. (2018), however also may coincide with the transition region, we found for the lowest gas temperatures. With the lower frequency VLBA data, Kovalev et al. (2020) did not constrain a geometry break. The single power law fit of the jet geometry is consistent with the l₁₅-value we find.

1308+326 (OP 313). This quasar features l and s-values which are indistinguishable within their uncertainties, along the entire jet length and the overlap region (Fig. A.20). The jet geometry overall appears conical and the s-values are consistent with the canonical parameter values from the Blandford & Königl (1979) model, which can be interpreted as a freely expanding jet. This is also consistent with Kovalev et al. (2020).

1633+382 (4C +38.41). This quasar jet features l and s-values which are indistinguishable within their uncertainties, along the entire jet length and the overlap region (Fig. A.21). In spite of strong scatter in the data, the jet geometry appears conical and the s-values are consistent with the canonical parameter values the Blandford & Königl (1979) model which can be interpreted as a freely expanding jet. Lower frequency data suggest an opening of the jet with a diameter gradient index of l ∼ 1.3 (Kovalev et al. 2020). Also for this source it would be interesting to combine the high and low frequency data in order to constrain a possible geometry transition.

1652+398 (Mrk 501). This is the second closest source in our sample and compared to most other sources is covered by only a relatively small number of observing epochs. It features diameter gradients at both frequencies which are consistent with a parabolic jet shape, as seen in Fig. A.22. This is suggestive of a jet that is observed on scales where it is still accelerating rather than freely expanding. This is also consistent with the T_b-gradients which are flatter than those expected for a freely expanding jet.

1730−130 (NRAO 530). This is a quasar jet that features l and s values which are consistent with a freely expanding jet, when analyzing the overlap region at both frequencies (Fig. A.23). Along the entire jet length, however, the 15 GHz jet seem to be more collimated and features a flatter T_b-gradient. This is dominated by the influence of a rather distant component between 20 ≲ d_c/mas ≲ 30 from the core which appears in epochs between ∼2006 and ∼2012. This component is relatively compact and features comparatively high T_b-values with respect to the distance from the core. This might be suggestive of a shocked region where the jet recollimates. Lu et al. (2011) find this component to be stationary. At lower frequenceie Kovalev et al. (2020) also do not constrain a geometry transition. The scaling of the jet diameter is consisntent with what is found in this study.

1749+096 (OT 081). Model-fit results for this BL Lac object are shown in Fig. A.24. Its diameter gradients are consistent with a conical jet shape while the brightness temperature gradients are steeper than expected within the simple Blandford & Königl (1979) jet model. This suggests that the parameters n and b are responsible for the steep T_B gradients. A poloidal magnetic field component could explain a steepening of the brightness temperature gradients while maintaining a conical jet shape.

2200+420 (BL Lacertae). The inner region of this prototypical BL Lac object features an overall jet shape which is consistent with a cone which is also confirmed by Casadio et al. (2021) by additionally taking GMVA data at 86 GHz into account. If fit with a simple power-law model (Fig. A.25), the s-values are consistent with the canonical value expected for a simple Blandford & Königl (1979) jet model. At ∼1.5 mas from the respective cores, we find a region where the jet locally shows a dip in its diameter and excess T_b values which is also consistent with the findings of Casadio et al. (2021). They associate this region with a recollimation shock which is connected to the changing external pressure gradient close to the Bondi radius (Kovalev et al. 2020). The large ${\chi }_{\text{red}}^2$ values suggest that the overall jet structure is complex and is not sufficiently represented by a single power law. The two-zone diameter-gradient model is shown in Fig. 17 (top) and features a parabolic jet shape out to d_c ∼ 2.0 mas transitioning to a hyperbolic shape jet further downstream which is consistent with the findings of Kovalev et al. (2020). Analyzing the jet in more detail it becomes clear that the jet of BL Lacertae is even more complex and still not fully described by a two-zone diameter-gradient model. To emphasize this we also show Fig. 17 (bottom). Here the break points are not found rigorously but are a ’fit by eye’. In this representation, the jet has a parabolic region up to roughly the point where the transition was fitted in the single-transition zone model, opening up to an extreme jet geometry (l ∼ 2.5), before transitioning again at ∼5 mas into a quasi-conical jet. The latter transition is roughly coinciding with the estimated location of the Bondi radius. This seems to contradict the findings of Casadio et al. (2021) where they do not find a geometry transition, additionally analyzing 86 GHz GMVA data. Casadio et al. (2021) argue that the region between 0.5 ≲ r/mas ≲ 1.5 should be excluded from a fit because in this region the jet shows a different expansion rate accompanied by a decrease in brightness which we also find in the brightness-temperature data (Fig. A.25). Not excluding this region may yield a parabolic geometry of this jet. This tentative contradiction can be resolved by the jet’s complex structure where clearly different regions with different behavior in the jet can be found. This would not affect the overall jet shape, measured over several decades in the distance from the respective radio cores, but would allow local geometry transitions of the jet.

We also tested the effect of the measured core shift of ∼0.1 mas, reported by O’Sullivan & Gabuzda (2009). The resulting jet geometry and brightness temperature gradients along the entire jet length are indistinguishable within the errors from the case where the core shift is not considered. Also the found geometry break position and the inner and outer diameter gradient do not differ.

2223−052 (3C 446). This quasar jet features s values which are indistinguishable regarding the two frequencies and are well described by the Blandford & Königl (1979) jet model. The diameter gradients are also indistinguishable, however there seems to be peculiarity within the 15 GHz data. While the 43 GHz jet shape is consistent with a cone, the 15 GHz jet shape is suggestive of a parabolic shape. This might be explained by the large scatter (as it is found to be characteristic for quasars) and the relatively low number of observing epochs for this source. A systematic jet shape transition however cannot be found in this source. The overall parabolic jet shape is consistent with the findings of Kovalev et al. (2020).

2230+114 (CTA 102). This quasar (see e.g., Fromm et al. 2011, 2013a,b) features s-values which are consistent with a parabolic jet in equipartition at both frequencies. The jet geometry can consistently be fit with a parabola at both frequencies, however along the entire jet axis at 15 GHz, the l-value becomes also consistent with a cone due to the larger uncertainties introduced by the larger scatter of the data beyond d_c ∼ 2 mas. A transition in this source cannot be constrained.

2251+158 (3C 454.3). Similarly to 1222+216 (4C +21.35) and 1253−055 (3C 279), the s-vlaues are consistent with a Königl (1981) jet model. However also in this source, the 43 GHz jet diameter gradient is consistent with a parabolic jet shape and is more collimated than the 15 GHz jet diameter gradient which is consistent with a conical jet, in both the overlap region and the respective entire jet lengths. Also here when rigorously searching for possible break points and an inner and outer slope, it cannot be constrained properly, due to the large scatter of the data points between 0.2 ≲ d_c/mas ≲ 3 in the 43 GHz jet. This jet also might be a candidate where at higher frequencies, more collimated, inner structures of a multilayered jet are observed. The overall parabolic jet shape is consistent with the findings of Kovalev et al. (2020).

7. Summary and conclusions

Using dual-frequency VLBI observational data for a sample of 15 FSRQs, 11 BL Lacs, and 2 radio galaxies at 15 GHz and 43 GHz, the geometry and brightness temperature gradients along the jet axes has been determined and interpreted. Our main observational findings are the following:

• The geometry of the VLBI jets in our sample can generally be described as power laws with indices l in the range between 0.5 and 1, i.e., collimated jet geometries.

• At 43 GHz, the jets tentatively show a stronger collimation than at 15 GHz. This is a significant difference if the full jet extent seen at both frequencies is considered, but not significant when the analysis is restricted to the overlap region where the jet is imaged at both frequencies (Table 3).

• FSRQs show stronger collimation than BL Lacs at both frequencies, but again this is not significant if restricted to the overlap region (Figs. 3 and 5).

• FSRQs show a more complex structure than BL Lac objects (with more local variations of the brightness temperature gradients leading to higher χ² values of the power law fits in Figs. 3 and 5). This is true for the jet geometry, indicating regions where the jets collimate and widen locally, as well as for the brightness temperature gradients indicating regions in the jet where the flux density is locally higher than in other regions.

• We found parabolic-to-conical geometry transitions in 5 FSRQs: CTA 26 (0336−019), PKS 0528+134, 4C +71.07 (0836+710), 4C +29.45 (1156+295), 3C 279 (1253−055). In addition, we found parabolic-to-conical geometry transitions in 3 BL Lacs: 3C 66A (0219+428), Mrk 421 (1101+384), and BL Lacertae (2200+420); and in a radio galaxy: 3C 111 (0415+379).

• The jet radius in the transition zone is smaller in BL Lac objects than in FSRQs and the geometry transitions occurs further down the jet for FSRQs which is likely due to systematically smaller SBMH masses in BL Lacs.

• The FSRQs CTA 26 (0336−019) and 4C +71.07 (0836+710) show a stationary component downstream the jet past the Bondi radius.

• BL Lacertae (2200+420) shows a jet structure of substantially greater complexity than a simple transition from a parabolic to a conical shape due to local effects like pressure gradient change and locally different expansion rates between 0.5 ≲ r/mas ≲ 1.5, (Casadio et al. 2021).

The high-frequency VLBI observations apparently resolve jets on physical scales where in many cases their shapes change from parabolic (collimated) to conical (free) morphologies. This observational finding is in agreement with relativistic magnetohydrodynamical jet models, e.g., Fendt (2006), Porth et al. (2011). These models show that the accretion disk in blazars with high accretion rates (FSRQs) flattens the magnetization profile at the jet base compared to jets from disks with low accretion rates (BL Lacs), either by advection of magnetic flux with the disk or by generating it with a disk dynamo. The radially flattened magnetization profile and the strong disk wind envelope are predicted to lead to (i) a collimation of the toroidal-field dominated jet extending to large distances from the core, (ii) morphological complexity as these jets are prone to produce axially propagating knots, and (iii) high terminal speeds. By contrast, jets in BL Lacs with steeper magnetization profiles at their base and weaker wind envelopes are expected more uniform, slower, and widen closer to the core (Potter & Cotter 2015), which is in good agreement with our observational results. The acceleration of the jets is driven by the conversion of Poynting flux to kinetic energy flux of particles entrained from the surrounding disk wind or injected by reconnection events or electrostatic gaps at the jet base (Wendel et al. 2021). The breakdown of the collimation (due to the kink instability) first affects the low-frequency sheath of the jet, in accord with the observational results reported here. On the kinetic level of description, the energy conversion is likely to be associated with turbulence, anomalous resistivity, and ultrarelativistic particle acceleration. These dissipative processes involve complex 3D-substructures in the jet plasma leading to fast flux variability at high energies. In this scenario, no direct correlation with the extent of the Bondi radius around the central supermassive black hole is predicted, which is in agreement with our results. Further improvements of this type of studies are possible: adding polarimetric analyses and increasing the sample size, the dynamical range, and the accessible bandwidth of VLBI observations, holds the promise to further advance our understanding of the physical processes involved in the blazar jet phenomenon.

Acknowledgments

We are grateful to Till Steinbring, who did important preliminary work on which this study partly has been based. We appreciate valuable scientific discussions with the MOJAVE team, especially Matthew L. Lister, Yuri Y. Kovalev, and Evgeniya V. Kravchenko, which helped us to improve this paper. We thank the referee C. Casadio for reviewing the manuscript which helped to improve the quality of this paper significantly. This research has made use of data from the MOJAVE database that is maintained by the MOJAVE team (Lister et al. 2009, 2018). This study makes use of VLBA data from the VLBA-BU Blazar Monitoring Program (BEAM-ME and VLBA-BU-BLAZAR; http://www.bu.edu/blazars/BEAM-ME.html), funded by NASA through the Fermi Guest Investigator Program. The VLBA is an instrument of the National Radio Astronomy Observatory. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated by Associated Universities, Inc.

References

Appendix A: Single power law fits

Fig. A.1. 0219+428 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.2. 0336-019 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.3. 0415+379 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.4. 0430+052 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.5. 0528+134 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.6. 0716+714 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.7. 0735+178 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.8. 0827+243 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.9. 0827+243 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.10. 0836+710 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.11. 0851+202 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.12. 0954+658 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.13. 1101+384 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.14. 1127-145 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.15. 1156+295 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.16. 1219+285 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.17. 1222+216 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.18. 1226+023 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.19. 1253-055 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.20. 1308+326 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.21. 1633+382 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.22. 1652+398 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.23. 1730-130 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.24. 1749+096 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.25. 2200+420 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.26. 2223-052 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.27. 2230+114 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Fig. A.28. 2251+158 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.

Appendix B: Comparing different fit methods

In this study we made use of the orthogonal distance regression (ODR) as used in Baczko et al. (2016), Baczko et al. (2019) and Baczko et al. (2022) in order to fit the smoothly broken power law in the geometry transition sources. For the discussion in this paper, we use the least square implementation in scipy with the trust-region reflective algorithm in order to fit the smoothly broken power law which allows us to set bounds on each parameter which we get from the numerical method when searching for a possible break in geometry in each jet. Figure B.1 shows the fit parameters acquired from the ODR method and the least square method. Both well established methods yield results that are undistinguishable within the error bars.

Fig. B.1. Comparison of the orthogonal distance regression and the least square fitting method. Overall the fitting results for the different parameters are in agreement within the 1-σ error bars. Top left: Jet radius at the break point; Top right: break point; Bottom left: inner diameter gradient; Bottom right: outer diameter gradient.

All Tables

All Figures

Fig. 1. Fitted Gaussian-component representation of 4C +71.07 (0836+710) at 15 GHz (top) and at 43 GHz (bottom). The contours for both frequencies are shown at a 5σ level and are separated by a factor of 2. In the 15 GHz image the map noise level is ∼0.2 mJy. The contour levels are [−0.85,0.85,1.7,3.4,6.7, 13.7,27.5,54.9,109.8,219.6,430.2,878.5] mJy). At 43 GHz the map noise level is ∼0.8 mJy. The contour levels are [−3.8,3.8,7.5, 15.1,30.2,60.3,120.6,241.2,284.5,482.5,964.9, 1929.8, 3859.6] mJy.

In the text

	Fig. 2. Distribution of l-values along the entire jet (top panel) and along the overlap-region (bottom). Values at 15 GHz are shown in blue and 43 GHz in red. The histograms are shown for equal bin width of Δl = 0.2 (dashed lines) and for Bayesian binning (solid lines; see text).
In the text

Fig. 3. Distributions of the geometry scaling parameter l (left) and its variance χ2 (right) shown with Bayesian (solid lines) and uniform binning (dashed lines) for FSRQs (red) BL Lacs (bue), and RGs (green). Upper four panels: distributions for the entire jet at both frequencies; Lower four panels: distributions for the overlap-region at both frequencies. Bin sizes were chosen as 0.2 for the l distributions and 1.0 for the variance distributions.

In the text

	Fig. 4. Distribution of power-law slopes s of the brightness-temperature gradient along the entire jet (top panel) and along the overlap-region (bottom) at 15 GHz (blue) and 43 GHz (red). The histograms are shown for uniform bin width Δs = 0.5 (dashed lines) and for Bayesian binning (solid lines; see text).
In the text

Fig. 5. Bayesian binning (solid lines) and uniform binning (dashed lines) of the TB power-law slope distributions. The colors indicate the source classes (red: FSRQs, blue: BL Lacs, green: RGs). Top two rows: distributions along the entire jet for both frequencies; Bottom two lines: distributions along the overlap-region of both frequencies. Left column: distribution of s-values. For both frequencies, the bin size is 0.5. Right column: ${\chi }_{\text{red}}^2$ value distribution describing the brightness temperature. For both frequencies, a unit bin size was chosen.

In the text

Fig. 6. P-values for each pair of parameters involving the jet geometry compared with the two-sided KS test. The fields are color coded ranging from p < 5% (dark blue) to p ≥ 10% (yellow). The lower triangle, where all l values are compared and the upper triangle where all ${\chi }_{\text{red}}^2$-values are compared, in the plot have to be read separately. The indices (Q, B, G) denote FSRQs, BL Lacs and radio galaxies, respectively. The index “ol” denotes values that are calculated in the overlap region in contrast to the entire jet length. These were introduced in order to make the plots more readable.

In the text

	Fig. 7. P-values for each pair of parameters involving the TB gradients compared with the two-sided KS test. The fields are color coded ranging from p < 5% (dark blue) to p ≥ 10% (yellow). The lower triangle, where all s values are compared and the upper triangle where all ${\chi }_{\text{red}}^2$-values are compared, in the plot have to be read separately.
In the text

	Fig. 8. Parameter space of n, b and α, spanned by the measured s and l values at both frequencies.
In the text

	Fig. 9. Jet diameter vs. core distance for 3C 66A (0219+428). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.3 mas (dashed line). The dotted lines indicate the estimated location of the Bondi radius.
In the text

	Fig. 10. Jet diameter vs. core distance for CTA 26 (0336-019). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.1 mas (dashed line). The dotted line indicate the estimated location of the Bondi radius.
In the text

	Fig. 11. Jet diameter vs. core distance for 3C 111 (0415+379). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 5.4 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.
In the text

	Fig. 12. Jet diameter vs. core distance for PKS 0528+134 (0528+134). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.1 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.
In the text

	Fig. 13. Jet diameter vs. core distance for 4C +71.07 (0836+710). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.7 mas (dashed line). The dotted line indicate the estimated location of the Bondi radius.
In the text

	Fig. 14. Jet diameter vs. core distance for Mrk 421 (1101+384). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 1.9 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.
In the text

	Fig. 15. Jet diameter vs. core distance for 4C 29.45 (1156+295). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 2.5 mas (dashed line). The dotted line indicate the estimated location of the Bondi radius.
In the text

	Fig. 16. Jet diameter vs. core distance for 3C 279 (1253-055). The jet shows a transition from a parabolic shape to a conical shape at xb ∼ 0.6 mas (dashed line). The dotted line indicates the estimated location of the Bondi radius.
In the text

	Fig. 17. Jet diameter vs. core distance for BL Lac (2200+420). Top: single transition zone model (rigorous fit). The shape transitions at xb ∼ 1.8 mas from parabolic to conical. Bottom: two-transition model where the break points are ‘fit by eye’ parametrizes the complex structure.
In the text

	Fig. 18. FSRQs (in red), and BL Lacs (in blue). Top: distribution of the break points in linear scale and deprojected. Bottom: distribution of the jet radii at the respective geometry break point.
In the text

	Fig. A.1. 0219+428 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.2. 0336-019 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.3. 0415+379 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.4. 0430+052 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.5. 0528+134 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.6. 0716+714 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.7. 0735+178 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.8. 0827+243 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.9. 0827+243 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.10. 0836+710 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.11. 0851+202 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.12. 0954+658 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.13. 1101+384 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.14. 1127-145 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.15. 1156+295 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.16. 1219+285 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.17. 1222+216 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.18. 1226+023 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.19. 1253-055 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.20. 1308+326 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.21. 1633+382 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.22. 1652+398 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.23. 1730-130 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.24. 1749+096 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.25. 2200+420 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.26. 2223-052 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.27. 2230+114 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. A.28. 2251+158 jet diameter dj as a function of the distance from the core dc and the brightness temperature TB as function of dc. The color coding for both frequencies indicates the epoch, where lighter colors indicate more recent epochs.
In the text

	Fig. B.1. Comparison of the orthogonal distance regression and the least square fitting method. Overall the fitting results for the different parameters are in agreement within the 1-σ error bars. Top left: Jet radius at the break point; Top right: break point; Bottom left: inner diameter gradient; Bottom right: outer diameter gradient.
In the text