Multiband RadioAstron space VLBI imaging of the jet in quasar S5 0836+710

Detailed studies of relativistic jets in active galactic nuclei (AGN) require high-fidelity imaging at the highest possible resolution. This can be achieved using very long baseline interferometry (VLBI) at radio frequencies, combining worldwide (global) VLBI arrays of radio telescopes with a space-borne antenna on board a satellite. We present multiwavelength images made of the radio emission in the powerful quasar S5 0836+710, obtained using a global VLBI array and the antenna Spektr-R of the RadioAstron mission of the Russian Space Agency, with the goal of studying the internal structure and physics of the relativistic jet in this object. The RadioAstron observations at wavelengths of 18cm, 6cm, and 1.3cm are part of the Key Science Program for imaging radio emission in strong AGN. The internal structure of the jet is studied by analyzing transverse intensity profiles and modeling the structural patterns developing in the flow. The RadioAstron images reveal a wealth of structural detail in the jet of S5 0836+710 on angular scales ranging from 0.02mas to 200mas. Brightness temperatures in excess of $10^{13}$\,K are measured in the jet, requiring Doppler factors of $\ge 100$ for reconciling them with the inverse Compton limit. Several oscillatory patterns are identified in the ridge line of the jet and can be explained in terms of the Kelvin-Helmholtz (KH) instability. The oscillatory patterns are interpreted as the surface and body wavelengths of the helical mode of the KH instability. The interpretation provides estimates of the jet Mach number and of the ratio of the jet to the ambient density, which are found to be $M_\mathrm{j}\approx 12$ and $\eta\approx 0.33$. The ratio of the jet to the ambient density should be conservatively considered an upper limit because its estimate relies on approximations.


Introduction
Very long baseline interferometry (VLBI) that combines groundbased radio telescopes with an orbiting antenna to engage in space VLBI observations is an established tool for reaching an unprecedentedly high angular resolution of astronomical measurements (see Arsentev et al. 1982;Levy et al. 1986;Hirabayashi et al. 2000 for overviews of earlier space VLBI programs). The ongoing international space VLBI mission Ra-dioAstron (Kardashev et al. 2013) led by the Astro Space Center (ASC, Moscow, Russia) reaches a resolution of ∼ 10 µas at a frequency of 22 GHz. The space segment of RadioAstron employs a 10-meter antenna (Space Radio Telescope, SRT) on board the satellite Spektr-R, which launched on 18 July 2011. Spektr-R has a highly elliptical orbit with a period of 8 to 9 days and an apogee reaching up to 350 000 km. The SRT operates at frequencies 0.32, 1.6, 5, and 22 GHz (Kardashev et al. 2013). The SRT delivers data in both left (LCP) and right (RCP) circular polarization at 0.32, 1.6, and 22 GHz. At 5 GHz, the SRT provides only LCP data. The RadioAstron observations are correlated at three different facilities: the ASC RadioAstron correlator in Moscow (Likhachev et al. 2017), the JIVE 1 correlator in Dwingeloo (Kettenis 2010), and the DiFX correlator (Deller et al. 2007(Deller et al. , 2011 of the MPIfR 2 in Bonn, which uses a dedicated DiFX correlator the latter (see Mizuno et al. 2012Mizuno et al. , 2014, but it would require observations at frequencies higher than 22 GHz to uncover the relation between the two types of instability in S5 0836+710. The jet ridge line measured in the VSOP and ground VLBI images can indeed be reconciled with the presence of KH instability in the flow Perucho & Lobanov 2007;Perucho et al. 2012a;Vega-García et al. 1998). The oscillations observed in the ridge line with respect to the overall jet direction were associated with the helical, H s , and elliptical, E s , surface modes of Kelvin-Helmholtz instability (Hardee & Stone 1997;Hardee 2000), which introduce ellipticity in the cross-section of the flow (E s mode) and force the entire flow to oscillate around its main propagation direction (H s mode). Investigations of more intricate body modes of the instability that affects the jet interior (see Hardee 2000) have not been feasible so far for S5 0836+710, owing to insufficient transverse resolution of the jet structure in VLBI images made for this object.
Observations of the jet made with MERLIN at arcsecond scales showed an emission gap between 0.2 and 1.0 and a large-scale structure beyond 1 . This was explained by the acceleration and expansion of the jet and its disruption due to a helical instability (Perucho et al. 2012b). The presence of shocks has been suggested in the jet on scales up to ∼ 0.5 kpc, which were revealed by multiple regions with flatter spectrum that are separated by 5 mas .
The multiband RadioAstron observations of S5 0836+710 presented here are described in Sect. 2. The resulting images are introduced in Sect. 3, the interpretation of the asymmetric structure of the jet is discussed in Sect. 3.2, and the large-scale jet structure is analyzed in Sect. 4. The results are summarized and discussed in a broader context in Sect. 5.

RadioAstron observations of S5 0836+710
RadioAstron imaging observations of S5 0836+710 were performed at two epochs: on 24 October 2013 at L band (1.6 GHz, wavelength λ=18 cm) and on 10-11 January 2014 at C band (5 GHz, λ=6 cm) and K band (22 GHz, λ=1.3 cm), combining data recorded at the SRT with data provided by a ground array of radio telescopes. At each epoch the source was observed for about 18 hours, starting from a perigee passage of the SRT, in order to ensure a continuous coverage of Fourier domain (uv coverage) by the ground-space baseline data. The resulting uv coverage of the RadioAstron observations are shown in Fig. 1.
The L-band data were recorded at all participating antennas in dual-circular polarisation, comprising the LCP and RCP channels. The C-and K-band data were recorded at the SRT simultaneously in the mixed C/K mode, with the SRT recording in the LCP mode simultaneously at both frequencies (with two IF bands per frequency). At each frequency, the SRT recording was supported by a subset of the ground telescopes (see Table 1 for details of the frequency allocation) recording in the dual-circular polarization mode at the selected frequency band.
General technical parameters of the telescopes participating in the different observations are given in Table 1. The individual observations and observational setups are summarized in Table 2. The data were recorded at a rate of 128 megabits per second (Mbps) split into two intermediate-frequency (IF) bands, resulting in a total bandwidth of 2x16 MHz for each circular polarization channel (Andreyanov et al. 2014). The tracking stations at Puschino (Andreyanov et al. 2014) and Green Bank (Ford et al. 2014) received and recorded the SRT data and telemetry.
The RadioAstron observations were complemented with simultaneous VLBA observations at 15 GHz (λ = 2 cm) and Column designation: D : antenna diameter ( † -equivalent diameter). SEFD: system equivalent flux density (system noise) for the frequencies at which each individual antenna participated in the observations. 43 GHz (λ = 0.7 cm) made during gaps between the SRT scans, which are required for cooling of the onboard high-gain antenna hardware (consult Table 2 for the duration of the on/off cycles at the SRT). The respective uv coverage and images obtained from these observations are shown in Appendix A (Figs. A.3,A.5,A.7,and A.8).
After each of the observing runs, the respective data from all participating telescopes were transferred to the VLBI correlator of the MPIfR. For the correlation and further calibration of RadioAstron data, information on the orbital position, velocity, and acceleration of the space antenna is introduced in the correlator model. The orbit and the momentary state vector of the SRT were reconstructed by the RadioAstron ballistic team by combining information of the radiometric range and radial velocity, measured in Bear Lakes and Ussurijsk (Russia), Doppler-tracking was performed at both tracking stations, and measurements of the sky position of the satellite were obtained from laser ranging and optical astrometry Zakhvatkin et al. 2014).

Data correlation
The data were processed using the DiFX correlator of the MPIfR at Bonn upgraded for RadioAstron data correlation . Fringe-searching for the SRT was performed separately for each individual observing scan, in order to optimize the centering of the correlation window, and to minimize residual acceleration terms (see Lobanov et al. 2015, for a more detailed discussion). An initial search window with 1024 spectral channels per IF band and 0.1 s integration time was used. Whenever Coverage of the Fourier domain (uv coverage) of the Ra-dioAstron observations of S5 0836+710 at 1.6 GHz (top), 5 GHz (middle), and 22 GHz (bottom) plotted in units of Mλ at each respective frequency. Central concentrations correspond to baselines between the ground antennas (ground-ground baselines). Long arcs represent baselines between the ground antennas and the SRT (ground-space baselines).
feasible, the delay and rate solutions for scans with no fringe detections were interpolated from the adjacent scans with successful fringe detections.
In the L-band observations, about half of the RadioAstron scans yielded fringes at the correlation stage, corresponding to the initial ∼50 minutes of observations, performed with the Green Bank tracking station, and a further ∼3 hours (distributed from 00 UT to 07 UT) with the Puschino tracking station. The final ∼3 hours (from 08 to 14 UT) of observing time were tracked again by Green Bank, and initial best-guess delays and rates were extrapolated from the earlier scans. The detection significance progressively diminished with increasing baseline length, Residual fringe delay (top) and delay rate (bottom) solutions at 5 GHz, obtained for the SRT from the data recorded at the tracking stations in Green Bank (light gray, UT 04-12, 10 January and UT 01-08, 11 January) and Puschino (dark gray, UT 16-24, 10 January). The squares and triangles denote the solutions for the first and second IF band, respectively. The time variations of the residual delays and the delay rates likely reflect the effect of the accelerated motion of the SRT, which was not accounted for during the fringe fitting.
dropping below a signal-to-noise ratio (S/N) of 10 at baselines longer than six Earth diameters, progressively resolving out the source structure at longer space baselines.
For the C-and K-band observations, only 4 scans out of 30 gave fringes on space baselines with an S/N>10 at the correlator stage. These scans corresponded to the perigee part of the orbit (∼40 minutes, from 16 UT to 18 UT). Similar to the approach employed for the L-band data, the initial delay and rate for the remaining RadioAstron scans were extrapolated (∼2.5 hours from 4 UT to 11 UT, and ∼4 hours from 18 UT to 08 UT) and supplied to the correlator. The final correlation of the data was made with 64 spectral channels per IF and 0.5 s of integration time for the L band, and 32 spectral channels per IF and 0.5 s of integration time for the C and K bands. These extrapolated solutions were subsequently refined and improved during fringe fitting performed as part of post-processing of the correlated data.

Post-correlation data reduction
The correlated data were reduced using AIPS 5 for the initial calibration and Difmap (Shepherd 1997(Shepherd , 2011 for imaging. The a priori amplitude calibration was applied using nominal values for the antenna gains and system temperature measurements made at each antenna during the observation. For the SRT, the sensitivity parameters measured in 2011-2013 (Kovalev et al. 2014) were used. The data were edited using the flagging information from each station log. Parallactic angle correction was applied to the ground-array antennas to account for the axis rotation of the dual-polarization antenna feeds with respect to the target source.

Fringe fitting
The data were fringe fit in two steps, first only for the ground array, and then also including the space baselines, with stacked solutions for the ground baselines (baseline stacking) used to improve the detectability of space baseline fringes. For the observations at 1.6 GHz, the solution interval for the ground-and space-array fringe fitting was 4 min and the cutoff S/N = 3 was set. Effelsberg was used as reference antenna for the first fringe fitting and Green Bank for the second. For the 5 GHz data, the solution interval was 1 min and 4 min for the ground and space arrays, respectively. The S/N cutoff was set to S/N = 4.3 for the ground VLBI data and S/N = 3 for the space baselines. Effelsberg was the reference antenna for all baselines. For the observations at 22 GHz, the solution intervals for the ground and space arrays were 1 min and 5 min, respectively. The S/N cutoffs and the reference antenna were the same as used for the data at 5 GHz. For the SRT, a two-way running average smoothing over a time interval of 30 minutes was applied to the solutions obtained at 5 and 22 GHz. Figure 2 shows the resulting fringe delay and delay rate solutions obtained for the SRT at 5 GHz. The solutions are consistent between the two IF bands, and the likely effect of the satellite acceleration is reflected both in the time variations of the fringe delays and in the respective delay rates. The overall behavior of the fringe solutions is in agreement with the expected accuracy of the orbital position and velocity determination (Kardashev et al. 2013;Duev et al. 2015;Lobanov et al. 2015).

Bandpass calibration
For the 1.6 GHz data, corrections for the individual receiver bandpasses were introduced. For this calibration, the autocorrelated data were used and phases and amplitudes were corrected. The correction was applied for the whole time interval, and Effelsberg was used as the reference antenna. For the 5 GHz and 22 GHz data, the attempted bandpass corrections did not show any improvement and were not implemented.

RadioAstron images of S5 0836+710
After the a priori amplitude calibration, fringe fitting, and bandpass calibration, the data were averaged in time over 10 s and in frequency over all channels in a given IF band and exported to Difmap for imaging and further analysis. The data were imaged using the CLEAN hybrid imaging algorithm (Högbom 1974) and self-calibration, with amplitude adjustments limited to a single gain correction applied to each antenna over the entire duration of the respective observation. The natural weighting (increasing the sensitivity to extended emission) was applied to the visibility data for obtaining the ground array images and the uniform weighting (maximizing the image resolution) was used to produce the space VLBI images (see Briggs et al. 1999, for a detailed discussion of the data weighting in VLBI). Radial distributions of the final self-calibrated visibility amplitudes and the respective CLEAN models from the space VLBI images are plotted in Fig. A.9 of Appendix A. The figures show that structure is detected up to the longest ground-space baselines and that the CLEAN models agree well with the data.  Table 3.    Table 3. The curved blue line denotes the ridge line we derived that is discussed in Sect. 4. The 1.6 GHz images from the ground-array (Fig. 3) and space VLBI (Fig. 4) data show a remarkable wealth and complexity of the structure, with the jet extending well beyond 40 mas even in the full-resolution RadioAstron image. The comparison of the space-and ground-VLBI images at 5 and 22 GHz in Fig. 5 further demonstrates the remarkable improvement in the angular resolution provided by RadioAstron. The structure traced by the space baselines is located at the brightest and most central part of the jet. The jet seen at 1.6 GHz shows wiggles that could be part of a relativistic spiral with effects from Kelvin-Helmholtz instabilities (as suggested previously in Lobanov et al. 1998;Perucho et al. 2012a,b).
The RadioAstron image at 5 GHz shows that the jet is substantially bent at small scales. Two clear bends are observed in the two brightest regions. A Kelvin-Helmholtz instability developing in the jet may explain this type of structure. The 22 GHz RadioAstron image of S5 0836+710 yields the highest resolution image ever obtained for this source, reaching down to ∼ 15 µas. This corresponds to a linear scale of ∼ 0.13 pc, or ∼ 650 gravitational radii for a black hole mass of 2 · 10 9 M (Tavecchio et al. 2000). Owing to the larger errors of the visibility phase and amplitude as well as poor uv coverage on the space baselines, the full-resolution 22 GHz image can only trace the brightest regions of the jet flow, with their morphology suggesting a limbbrightened structure, which can be caused by strongly asymmetric emission within the jet, or a double-helical pattern, similar to that detected in 3C273 with VSOP observations (Lobanov & Zensus 2001).

Brightness temperature
We estimate the brightness temperature in the central regions of the jet from two-dimensional Gaussian model fitting of the visibility data, using the Levenberg-Marquardt algorithm implemented in Difmap. We compare the brightness temperatures calculated from the Gaussian modelfits with visibility-based estimates (Lobanov 2015) obtained from the data at the longest ground-space baselines of the observations. The Gaussian modelfits of the central regions and the corresponding estimates of brightness temperature are shown in Fig. 6 and listed in Table 4. At 1.6 GHz, a region of about 2.5 mas in extent was modeled. At the two higher frequencies, only the innermost part of the jet was modeled, limited to the structure before the first gap of emission as seen in Fig. 5. The model fitting was performed in Difmap, starting by removing the clean components in the region of interest, and replacing them with a set of fitted Gaussian components representing the respective structure.
At 1.6 GHz (left panel of Fig. 6), five components can describe the central region. For the innermost jet structure observed at 5 GHz and 22 GHz (central and right panels of Fig. 6), nine and four Gaussian components, respectively, are needed to describe the most central part. At 1.6 GHz, the entire inner jet has a brightness temperature of ≈ 10 12 K, with a maximum value of 3.3 × 10 12 K estimated for the component C4 located downstream from the apparent jet origin (see Table 4). At 5 GHz and 22 GHz, the respective maximum brightness temperatures are also achieved in components downstream in the jet, with 1.7 × 10 13 K estimated for C4 at 5 GHz and 3.8 × 10 12 K for C3 at 22 GHz. This result underlines the complex structure of the inner jet, with opacity (see Lobanov 1998;Gómez et al. 2016), jet acceleration (Lee et al. 2016), or fine-scale substructure (e.g., resulting from recollimation shocks; see Fromm et al. 2013;Lobanov et al. 2015;Gómez et al. 2016 (right) superimposed on images obtained at respective frequencies using only the ground-array data (colors). Ground-array images are made using natural weighting and space VLBI images are made using uniform data weighting. The respective restoring beams are plotted in the lower left corner in orange (ground array) and blue (space VLBI). Basic parameters of all images are listed in Table 3.  Table 4. jet base. In the reference frame of S5 0836+710, the respective brightness temperatures are higher by a factor of 1 + z, hence all of them exceed 10 13 K and require Doppler factors, δ IC = 20-100 in order to reconcile them with the inverse Compton limit of ≈ 5 × 10 12 K. Only the lower limit of the estimated δ IC range agrees with the kinematic Doppler factor δ kin ≈ 17 resulting from the existing measurements of the jet speed and viewing angle (Otterbein et al. 1998;Lister et al. 2013).
Estimating the brightness temperature directly from the visibility data (see Table 5) yields even higher values, with the minimum brightness temperature, T b,min exceeding 10 13 K at all three bands. When averaged over 10 % of the longest baselines in each data sets, the values of T b,min approach the modelfit estimates, but still remain higher. This may be viewed as another indication that there is an ultracompact structure in the jet that is not recovered by the Gaussian modelfit. The visibility-based estimates require Doppler factors of up to ∼ 300 to avoid the inverse Compton catastrophe. We therefore conclude that the RadioAs-tron observations of S5 0836+710 are indicative of intrinsic violation of the Compton limit on the brightness temperature.
Similar indications have also been seen in a number of other sources observed with RadioAstron (see Lobanov et al. 2015;Gómez et al. 2016;Giovannini et al. 2018), and they might imply that the innermost regions of the jet are either not in equipartition or subject to more peculiar physical conditions (Kellermann 2002) such as a monoenergetic electron energy distribution (Kirk & Tsang 2006). A more systematic study of the overall Ra-dioAstron measurements of the brightness temperature should give better insights into this matter.

Asymmetries of the jet structure
The improved resolution of the RadioAstron images at 5 GHz and 22 GHz enables detailed studies of the innermost section of the flow, on scales smaller than 1 mas. The 22 GHz image presented in Fig. 5 suggests that the jet is transversally resolved, [Gλ] L 5.0 × 10 13 0.64 3.9 × 10 12 0.66-0.73 C 3.4 × 10 13 2.25 2.9 × 10 12 2.34-2.60 K 1.2 × 10 13 11.2 5.8 × 10 12 10.4-11.5 Notes: T b,min : minimum brightness temperature obtained from the visibility data at a given baseline length. B. T b,min : average value of the minimum brightness temperature over the baselines in the interval B corresponding to the longest 10% of the baselines in the respective data sets.
with a filamentary structure tracing only one side of the flow. This apparent asymmetric transverse structure may in part be due to the limited dynamic range (estimated to be about 75:1) of the images. It may also result from differential Doppler boosting of the flow within the jet. Following Rybicki & Lightman (1979) (see also Aloy et al. 2000;Lyutikov et al. 2005;Clausen-Brown et al. 2011), the bright side of the jet changes from top to bottom at a critical viewing angle given by cos(θ r ) = β j for a helical magnetic field with its maximum asymmetry when the pitch angle is 45 • . In the case of S5 0836+710, this corresponds to θ r = 4 • . Therefore, the top side of the jet should be brighter, with the estimated jet viewing angle of 3.2 • smaller than θ r . This is indeed the case, as observed in Fig. 5, implying a helical field oriented counterclockwise (e.g., Aloy et al. 2000). Maximum asymmetry is obtained at θ = φ, with θ the viewing angle in the jet reference frame and φ the pitch angle in the fluid frame (Aloy et al. 2000). Fixing the viewing angle to 3 • and β j to that corresponding to Lorentz factor γ = 12 (β j = 0.99652), we obtain φ 64 • for a maximum asymmetry. This pitch angle implies the presence of a helical field with a slightly dominating toroidal component. Nonetheless, the asymmetry in emission could also indicate a physical asymmetry in the jet properties, revealing the presence of distinct regions that vary across the jet channel with time. This has been revealed by a number of VLBI observations that show components that only partially fill the jet cross-section, and the jet is only observed in its full width in stacked images over many epochs (e.g., Lister et al. 2013;Pushkarev et al. 2017;Beuchert et al. 2018). We can take the dynamic range of the 22 GHz image as a measure of the difference in brightness across the jet. Asymmetries in the transverse structure of the jet emission could be caused by pressure differences (Perucho et al. , 2012b. However, this effect alone cannot explain this difference in brightness if the perturbation is still developing in the linear regime (Perucho et al. 2005(Perucho et al. , 2006, and the effect of the magnetic field orientation in the jet frame on the synchrotron emissivity needs to be taken into account. Nonlinear effects such as large instability amplitudes are not expected this close to the jet base unless they respond to the development of short-wavelength fast-growing modes (e.g., Perucho & Lobanov 2007). However, in this case, we would expect a persistent and symmetric emissivity. We thus exclude the option that instability patterns generate this strong asymmetry in brightness.
We focus our discussion on the possibility that the asymmetry (of about the dynamic range of the observation) is produced by differential Doppler boosting. The observed flux density is given by where D is the Doppler factor, B and θ B are the intensity of the magnetic field and the angle between the magnetic field and the line of sight in the fluid frame, respectively, and α is the spectral index (S ν ∝ ν −α ). In this expression we ignore the dependence on the particle energy distribution (in particular, on the Lorentz factor limits of this distribution, γ min and γ max , which introduce dependencies on other parameters, such as the particle number density). Recalling that the inner jet region has a flat spectral index, α 0 (see Lobanov et al. 2006), and assuming that the total magnetic field intensity (as given by both the toroidal and poloidal components of the field) does not change much with the toroidal angle (implicitly assuming axisymmetry), we obtain the following expression for the brightness ratio between the top brighter (at least by a factor of ≈75) part of the jet (indicated with the superscript t) and the bottom part of the jet (superscript b) in terms of the differential Doppler factors, D: and in the following discussion, we assume this ratio to be ∼ 100. Taking into account that sin(θ B ) = D sin(θ B ) (the unprimed value refers to the value in the observer's frame), we obtain The terms D t and D b can only be different if there is a toroidal component of the velocity and henceforth the velocity vector changes across the jet. When we assume that the axial velocity component v z c is constant across the jet, the brightness asymmetry can be ascribed to the top-to-bottom change of the angle θ between the velocity vector and the line of sight, so that from Eq. 3, which allows us to derive θ t in terms of θ b if we approximate sin(θ B,b )/ sin(θ B,t ) 1/3 ∼ 1 . The result is shown in the left panel of Fig. 7 for the Doppler factor ratios of 5 (solid line) and 10 (dashed line). In the plot, θ b spans from 0 • to 10 • in the observer's frame. The required changes in the viewing angle can give estimates of the relative values of the toroidal velocity component, v φ , as v φ /v z = tan(ψ), with ψ defined as the deprojection of the half-angle formed by the velocity vectors on both sides of the jet (here we assumed a mean viewing angle of 3 • to deproject the half-angle between the velocity vectors). This ratio is plotted in the right panel of Fig. 7 for the same Doppler factor ratios. The plots show that low toroidal velocities can result in sufficiently large changes in the Doppler factor of different regions in the flow. When we assume in contrast that D t D b , then the brightness asymmetry has to be given by Taking into account that this ratio sin(θ B,t )/ sin(θ B,b ) sin(θ B,t )/ sin(θ B,b ) (D t D b ) and that 0 ≤ sin(θ B ) ≤ 1, the fraction can only give a value 100 if one or both values of θ B are very close to 0 • , where the sine can take values that differ by several orders of magnitude. We find this coincidence that both angles are aligned to such an accuracy less probable than the presence of low toroidal velocities. Furthermore, a toroidal component of the velocity has been shown to naturally arise in RMHD simulations of expanding jets due to the Lorentz force (Martí et al. 2016). Finally, the required azimuthal velocities are compatible with the axial velocity (Lorentz factor) reported for this jet (Perucho et al. 2012a) in terms of causality. We can conclude that a low toroidal velocity has to be present in the jet flow to explain the brightness asymmetry. Recent numerical GRMHD simularions of jet formation (McKinney & Blandford 2009;McKinney et al. 2012) indicated that the rotational speed reaches up to the Keplerian speed, v c (r j ), at a given jet radius, r j . For the reported black hole mass of 2 × 10 9 M (Tavecchio et al. 2000) and a jet radius r j ≈ 0.05 mas (≈ 0.4 pc) estimated from the RadioAstron images, the respective v c ≈ 0.015 c is substantially higher than the v φ estimates made above. This may result from a slowdown due to thermal mixing or shear, but likely uncertainties in all of the measured quantities involved in these calculations preclude us form making any firm conclusion on the matter. We plan to run RMHD simulations to further investigate the potential physical mechanisms that define the estimated rotation in this jet.

Ridge line calculation
Because of the limited dynamic range of the space VLBI images at 5 GHz and 22 GHz (in which the jet is transversally resolved), we base our quantitative analysis of the jet flow largely on the jet ridge line derived from the 1.6 GHz images of S5 0836+710. At the scales sampled by these images, plasma instability is expected to play an important role, inducing regular patterns into the flow. We study these patterns by estimating the ridge line of the jet. We define the jet ridge as the line that connects the peaks of one-dimensional Gaussian profiles fitted to the profiles (slices) of the jet brightness drawn orthogonally to the jet direction (see Lobanov et al. 1998;Perucho et al. 2012a, for similar previous studies), with the step between individual profiles set to be smaller than the beam size (a measure taken in order to ensure continuity of the brightness profiles recovered from adjacent slices). Similarly to the earlier studies of S5 0836+710, a position angle of −162 • was adopted for the jet direction. To obtain the ridge line, a code was developed using Python and the iMinuit package 6 . The fitting algorithm implements continuity conditions for the parameters fitted to adjacent slices, thus ensuring a robust and self-consistent description of the evolution of the ridge line and flow width along the jet. The first slice is taken across the jet so as to cross the peak of brightness in the image at a given frequency. We calculated the uncertainties of the fitted parameters from the S/N in the image using the method described in Schinzel et al. (2012).
We present the ridge lines obtained from the ground-and space-VLBI images of S5 0836+710 at 1.6 GHz in Figs. 3-4. The ridge lines obtained from the images at other frequencies are presented in Appendix A.
The complex ridge line structure suggests the presence of several periodic patterns in the flow. In the following, we assume that these oscillatory modes represent individual modes of Kelvin-Helmholtz (KH) instability developing in the flow (see Lobanov & Zensus 2001;Perucho & Lobanov 2007;Perucho et al. 2012a).

Modeling the KH instability modes
We initially identifed the potential modes using a wavelet scalogram calculated from the ridge line offsets measured in the ground-and space-VLBI images at 1.6 GHz (Fig. 8), with the Marr function as the kernel for the wavelet transform. The scalograms presented in Fig. 8 indicate oscillatory modes at wavelengths of ∼ 90 mas, ∼ 35 mas, ∼ 15 mas, and ∼ 10 mas. The wavelet scalograms also show that the oscillatory patterns may have a slightly increasing wavelength with distance, a typical behavior expected for a KH mode propagating in an expanding flow (Hardee 2000). These inferences are used as initial guesses for further modeling of the ridge lines.
We used the ridge line measured in the 1.6 GHz ground-and space-array images as the basis for our modeling. We first describe the observed offset, ∆ r(z), of the ridge line with a simplified model by fitting it with several different oscillatory terms with constant amplitudes and wavelengths: 6 https://pypi.python.org/pypi/iminuit where a i is the amplitude, ψ i the phase and λ i the wavelength corresponding to the i-th mode. Four modes were needed to minimize the χ 2 in both ridge lines. The fitting was performed in two ways: 1. Fitting the first mode and subtracting it from the original ridge line. Then, a second mode was fit and again subtracted. The process was iterated until the addition of a new mode did not provide any improvement to the fit. 2. Fitting all modes at the same time. This method depends more strongly on the initial guess parameters. The first three modes were fit with the wavelengths observed in the scalogram as an initial guess.
For both methods, the best fit was chosen to be that which reduced the χ 2 . The two methods yield similar results. The wavelengths of the modes are listed in Table 6 for the ground-image ridge line and in Table 7 for the space-VLBI image ridge line. Figure 9 shows the resulting fit for the two ridge lines. The resulting wavelengths obtained from two independent fits are similar in both cases. This corroborates the robustness of our identification of the oscillatory modes and the little role played by jet expansion on the wavelength, taking into account the small jet opening angle (< 1 • Perucho et al. 2012a).  97.0 ± 4.0 3.70 ± 0.20 0 ± 10 2 44.0 ± 8.0 0.42 ± 0.04 82 ± 14 3 16.9 ± 1.4 0.50 ± 0.12 −106 ± 7 4 8.8 ± 0.1 0.23 ± 0.07 54 ± 6 When we compare these results with earlier works, it is interesting to note that the wavelengths identified in the ridge line at 1.6 GHz agree well with those that were obtained previously from the VSOP image of the source , where oscillations with wavelengths of 100, ∼ 8, and ∼ 5 mas were identified. The longest of these wavelengths is reasonably close to the λ = 97 or λ = 102 mas of the mode 1 in our RadioAstron data, and the two shorter wavelengths are close to λ = 6.6 mas from the ground-array image and λ = 8.8 mas of the mode 4 from the space-VLBI image. It should be noted that identification of short wavelengths in a ridge line measured in different VLBI images can be affected by the differences in the image noise and the uv coverage of the observations. It is therefore more difficult to give a robust representation of short-wavelength oscillations in a ridge line. We note that a 80 mas mode was also detected in VLBA images of S5 0836+710 by Perucho et al. (2012a), who suggested the presence of oscillations with shorter and possibly variable wavelengths, and reported oscillations with a wavelength between 10 mas and 20 mas at r ≤ 40 mas distance from the jet origin and growing to 40 mas  For each of these wavelengths, the respective scalogram contains at least one region of statistically significant wavelet amplitudes, thus reaffirming the physical relevance of the fitted oscillatory modes. Fig. 9. Fits by oscillatory modes to the jet ridge lines in the ground-(left) and space-(right) VLBI images of S5 0836+710 at 1.6 GHz. In each plot, the red line is the total multi-mode fit and the blue lines represent contributions from the individual oscillatory modes as described in the legend.
at r ≥ 40 mas. Finally, Perucho et al. (2012a) also reported the detection of short-wavelength oscillations, with wavelengths 5 mas at r < 15 mas increasing to 7-8 mas at r > 20 mas. These can be related to mode 4 in Tables 6 and 7 (wavelengths of 6.6 and 8.8 mas). It is also relevant to stress that in our fits we did not include growth lengths for the modes, which can also introduce small differences with the data at the largest scales. Based on all these arguments, we conclude that our results are consistent with previous studies of the instability patterns detected in the jet of S5 0836+710.
Regarding the possible identification of the instability modes, because the first, longest-wavelength mode displaces the ridge line, it should correspond to a helical mode. It could be tentatively identified as an helical surface mode, H s .
A tentative identification of the other modes can be obtained using the characteristic wavelength, λ * = λ i (n i + 2m i + 1 2 ), where λ * is the characteristic wavelength, λ i the observed wavelength and n i and m i the two transverse wavenumbers of the mode. The characteristic wavelength should have a similar value for the different modes (Lobanov & Zensus 2001), assuming that all the observed wavelengths correspond to modes that are excited at their maximum growth rates. The second-longest wavelength should then be the first helical body mode, H b1 , whereas the third longest wavelength should be either the second, H b2 , or the third helical body mode, H b3 . We exclude the fourth mode because such small wavelengths are difficult to determine and the results for ground and space arrays do not reconcile, even when the errors are accounted for. For the case in which the third mode is identified as the second-order helical body mode, the characteristic wavelength is 119 ± 12 mas, and when the third mode is identified as the third-order body mode, the characteristic wavelength is 130 ± 12 mas.
Following the approach of Hardee (2000) and references therein, the basic physical parameters of the jet can be obtained from the identification of the different modes, using the linear analysis of the Kelvin-Helmholtz instability. With this approach, it is possible to obtain the Mach number of the jet and the ratio of the jet to the ambient density knowing the characteristic wavelength, the jet radius, jet viewing angle, the jet apparent speed, and the apparent pattern speed. The jet radius can be easily calculated. It corresponds to the point where the peak flux density of the jet is reduced to 1% along the jet direction (Wehrle et al. 1992).
In practice, to choose the slice where we measure the jet width, we considered the slice with a peak flux density 1% of the flux density of the first slice, corresponding to the brightest point of the jet. For this slice, we took the full width at half-maximum (FWHM) of the corresponding Gaussian profile as a measure of the jet width. Then we deconvolved it, using the following relation: R j [mas] = 0.5 θ 2 FWHM − b 2 , where θ FWHM is the FWHM of the corresponding slice and b is the beam size. This yields a jet radius of 2 mas or 16 pc. We note that in earlier works (Perucho & Lobanov 2007, 2011 a substantially higher value of 17 mas was estimated for the jet radius. This value corresponds to the apparent jet radius measured directly from the image, while we deconvolved this measurement with the FWHM of the resolving beam. The large radius considered in earlier works led to the identification of the longest observed wavelength with the first body mode (because the observed structures are relatively shortened when they are measured with respect to the jet radius).
For the jet viewing angle and for the speed we used the values given in Otterbein et al. (1998), that is: θ j = 3 • , and a Lorentz factor γ j = 12, which corresponds to an apparent speed of β app = 10.7. The apparent pattern speed was measured by comparing ridge lines for different epochs at 15 GHz in the images of the MOJAVE monitoring program. This calculation was made using only the highest S/N region at 1-4 mas distances from the jet origin. The measured speed was w app = 0.35 ± 0.25c. With these quantities, the Mach number M j and the ratio of the jet to the ambient density η can be calculated as (e.g., Hardee 2000) where the external jet Mach number, M x , and the intrinsic jet pattern speeds are calculated with and β j = β app sin θ j + β app cos θ j .
This yields a Mach number of M j = 12 ± 3, in contrast to the results given in Lobanov et al. (1998), who reported a Mach number of 6 and a density ratio of η = 0.33 ± 0.08. The ratio of the jet to the ambient density we obtained here is close to unity, which may indicate that the jet is surrounded by a very diluted cocoon. This restriction could be alleviated if the pattern speed we used represents an upper limit because the different modes travel at different speeds; the shorter wavelength modes travel faster. Taking into account that the pattern speed was calculated in a rather small region, it therefore probably corresponds to a mode with a small wavelength, which is difficult to measure. The effect of the speed on the ratio of the jet to the ambient density is also relevant. If we look at different values within the range allowed by the error, for example assuming a pattern speed of w app = 0.10c, the ratio of the jet to the ambient density changes by an order of magnitude, to η ∼ 0.02. Nevertheless, the jet Mach number remains unchanged. Another important remark is that the approximation used before is only valid for cases where M j 1, which is verified by our result, and in the scenario where a contact discontinuity (or a narrow shear-layer) separates the jet and the ambient medium. The possibility of a shear layer playing a role in the jet in S5 0836+710 has been suggested in earlier works (Perucho & Lobanov 2007, 2011 after solving the linear stability problem for sheared jets, using a set of jet parameters given by Lobanov et al. (1998). Arising naturally in spinesheath scenarios, such a shear layer will also affect the growth rates and wavelengths of the KH instability modes (see Mizuno et al. 2007;Hardee 2007).

Summary
Multiband VLBI observations of S5 0836+710 with RadioAstron provided images of the jet with an unprecedentedly high angular resolution, reaching down to 15 microarcseconds at 22 GHz, which corresponds to a linear scale of 0.13 pc. The radio source S5 0836+710 was observed with a ground-and space-VLBI array at the frequencies of 1.6 GHz on 24 October 2013, and 5 and 22 GHz on 10 January 2014. The source was observed with tracks of 16.5 hours at 1.6 GHz and over a 30-hour period (with a 4-hour gap) at 5 and 22 GHz.
Non-standard procedures were needed to detect interferometric fringes with the space antenna; the resulting residual delay solutions for the 5 GHz RadioAstron image are shown in Fig. 2. They show a weak time-dependence due to the acceleration of the space antenna. The source was detected on baselines as long as 10 Earth diameters at 1.6 GHz, and 12 Earth diameters at 5 GHz and 22 GHz.
Hybrid imaging was performed for the ground-and spacearray data. The latter yield resolutions of 3 to 10 times the ground resolution, reaching scales of ∼ 15 µas (∼ 0.12 pc) at 22 GHz. The full-resolution RadioAstron images show a much richer structural detail than the ground-array images at the same frequencies. At 5 GHz and 22 GHz, the jet is transversely resolved in the RadioAstron images, revealing a bent and asymmetric pattern embedded into the flow.
The observed brightness temperature of the jet is estimated for the three frequencies and is listed in Table 4. These estimates indicate that the highest brightness temperatures for each frequency are 3.3 × 10 12 K, 1.7 × 10 13 K, and 3.8 × 10 12 K for 1.6 GHz, 5 GHz, and 22 GHz, respectively. These estimates agree well with the minimum brightness temperature estimated from the visibility data, using the longest 10 % baselines. All of the estimates imply T b ≥ 10 13 K in the reference frame of the source and would require Doppler factors of up to ∼ 300 in order to reconcile them with the inverse Compton limit.
The internal structure of the jet was studied by analyzing the ridge line of the jet. For the lowest frequency of 1.6 GHz, the location of the ridge line observed in the ground-array and space-VLBI images were obtained by fitting transverse brightness profiles with a single-Gaussian component and recording the position of its peak with respect to the overall jet axis oriented at a position angle of −162 • . The ridge lines were represented with a simple model as the sum of multiple oscillatory terms. The parameters of these oscillatory modes are modeled and explained in the framework of a Kelvin-Helmholtz instability that develops in the flow. The wavelengths of the oscillations are found to be similar to those reported in previous studies of the source. These oscillations can be interpreted as resulting from the helical modes of the instability developing in the jet. Based on this interpretation, an estimate of the jet Mach number and the ratio of the jet to the ambient density is obtained of ∼ 12 and ∼ 0.33, respectively. This estimate can be further verified and refined using a more detailed numerical analysis of the jet stability.
Appendix A: Auxiliary plots: uv coverage, visibility amplitudes, hybrid images, and jet ridge lines.
This appendix presents auxiliary information and images obtained from the RadioAstron observations. Ridge lines determined from RadioAstron images at 5 GHz and 22 GHz are shown in Figs 2, 2, etc. ) times 7.0 mJy/beam. Image parameters are given in Table 3. A discussion of the jet ridge lines is presented in Sect. 4.   Image of the jet in S5 0836+710 made from the ground-VLBI data at 5 GHz. The contour levels are drawn at (−1, 1, √ 2, 2, etc. ) times 1.5 mJy/beam. The curved blue line denotes the ridge line of the jet we derived that is discussed in Sect. 4.