EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 653 (September 2021) A&A, 653 (2021) A10 Full HTML
Open Access
Issue
A&A
Volume 653, September 2021
Article Number A10
Number of page(s) 19
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/201937241
Published online 01 September 2021

© N. R. MacDonald and K.-I. Nishikawa 2021

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Open Access funding provided by Max Planck Society.

1. Introduction

Relativistic extragalactic jets are among the most persistent energetic objects in the universe. They are composed of collimated beams of magnetized relativistic plasma that can extend up to thousands (and in some cases millions) of parsecs from their host galaxies. The current theoretical paradigm postulates that the ultimate physical mechanism powering these relativistic outflows is the energy released from matter accreting onto spinning supermassive black holes (Blandford & Znajek 1977; Blandford & Payne 1982). The plasma content of these jets (and their central engines) remains an active area of research (see, e.g., Croston et al. 2018; Fan et al. 2018; Thum et al. 2018; Myserlis et al. 2018a,b; Enßlin et al. 2019; Anantua et al. 2020; Sikora et al. 2020; Emami et al. 2021; Mościbrodzka et al. 2021; Ricarte et al. 2021).

With the advent of global (and recently space-based) millimeter-wave very long baseline interferometry (VLBI), we are able to probe the polarized emission emanating from the innermost regions of a number of jets. In particular, the linearly and circularly polarized synchrotron emission from these jets carry imprints of both the strength and orientation of the collimating magnetic fields as well as the plasma content of each jet. Studying the nature of this synchrotron emission can be used to infer the physical conditions both within the jet and in the surrounding environment into which the jet propagates.

In parallel to this observational advance, modern computational resources have allowed for increasingly sophisticated numerical plasma simulations. In particular, 3D particle-in-cell (PIC) simulations (e.g., Nishikawa et al. 2014, 2016a, 2020; Alves et al. 2018; Guo et al. 2020) have enabled, for the first time, a self-consistent treatment of the kinetic effects occurring within relativistic plasma outflows, such as the following: plasma instabilities, jet shear and entrainment, and magnetic reconnection (see Birdsall & Langdon 1995 and Nishikawa et al. 2021 for a summary of PIC methods). These PIC simulations, however, are numerically intensive, and this kinetic precision comes at the cost of small (relative to relativistic magnetohydrodynamic – RMHD) simulation sizes. There exists, therefore, a synergy between PIC and RMHD jet modeling that should be explored in order to gain a more holistic understanding of both the micro and macro physics of relativistic outflows. While RMHD simulations can effectively model the large-scale fluid motions of the jet (e.g., Mizuno et al. 2020; Tchekhovskoy & Bromberg 2016; Fuentes et al. 2018; Fromm et al. 2019; Mukherjee et al. 2020, 2021), PIC can be used to model the microphysics and radiative processes occurring within the jet plasma (e.g., Sironi et al. 2015, 2021; Zhang et al. 2018, 2020; Petropoulou et al. 2019; Hosking & Sironi 2020; Davelaar et al. 2020). These kinetic-scale processes form a direct link to VLBI observations of the polarized synchrotron emission.

In this paper, we compare the radiative differences between two PIC jet simulations (computed using the TRISTAN-MPI code1; see Niemiec et al. 2008 and Buneman 1993), one in which the jet is composed of an electron-proton (e− − p+) plasma (i.e., a normal plasma jet), and the other in which the jet is composed of an electron-positron (e− − e+) plasma (i.e., a pair plasma jet). We make use of a polarized radiative transfer scheme (see MacDonald & Marscher 2018) that has been embedded into a ray-tracing code for post-process imaging of each numerical jet simulation. We use this ray-tracing code to create synthetic full Stokes (I, Q, U, and V) images of each numerical PIC jet simulation. This paper is organized as follows: In Sect. 2 we summarize the scaling relations used in our PIC jet simulations. In Sects. 3 and 4 we outline the radiative transfer theory adopted in our study. In Sect. 5 we present the results of our ray-tracing calculations through the e− − p+ plasma and e− − e+ plasma jet simulations. Finally, in Sect. 6 we present our summary and conclusions. We adopt the following cosmological parameters: Ho = 71 km s−1 Mpc−1, Ωm = 0.27, and ΩΛ = 0.73.

2. Scaling

PIC calculations are carried out in dimensionless grid units which must be scaled (via scaling factors) into physical units (i.e., cgs) as a post-process step for use in our ray-tracing calculations. The PIC simulations presented in this paper are small (120 × 120 × 240 cells) segments of the e− − p+ and e− − e+ PIC jet simulations published in Nishikawa et al. (2017).

The simulation cell size (Δl), the time step (Δt), and the speed of light (c) are computed in the following manner:

(1)

(2)

(3)

where by definition the grid values: Δlgrid = Δtgrid = cgrid = 1. The quantities: Δlscale, Δtscale, and vscale represent scaling factors that set the physical dimensions of each jet simulation. In both jet simulations presented in this paper, we set vscale = 3 × 1010 cm s−1.

The plasma skin depth (≡c/ωpe) parameterizes the length scale across which low-frequency electromagnetic radiation can propagate within a collisionless plasma. The quantity ωpe denotes the relativistic electron plasma frequency and is given by:

(4)

where ne, e, and me denote the number density [cm−3], charge [cm3/2 g1/2 s−1], and mass [g] of the electrons within the plasma and γe is the electron Lorentz factor. The plasma frequency is the fundamental frequency of oscillation of the electrons and together with the plasma skin depth sets a fundamental length scale within our PIC jet simulations.

In order to model the kinetic scales within our jet simulations, the plasma skin depth of the jet plasma is resolved across ten computational grid cells by the TRISTAN code:

(5)

Equation (5), therefore, defines the characteristic length scale of the numerical cells within our PIC jet simulations. We define a fiducial jet electron number density as follows:

(6)

to which we scale our PIC values. By combining Eqs. (1), (3)–(6) we obtain the following PIC cell size scale factor:

(7)

having set me = 9.1094 × 10−28 g, e = 4.8066 × 10−10 cm3/2 g1/2 s−1, and γe ≃ 15 (the limiting value at the jet base in our models). This scaling implies that the jet axis of both of our PIC simulations (i.e., 240 cells in length) corresponds to a physical size of ∼150 km. This physical size demonstrates the vastly different spatial scales probed by PIC simulations in comparison to RMHD and represents an inherent challenge PIC codes face when attempting to model the plasma of parsec scale astrophysical jets.

As in Eqs. (1)–(3), we scale the magnetic field within our PIC simulations in the following manner:

(8)

where the dimensionless grid value at the base of our jet models is Bgrid ≃ 2.7 (and varies throughout each simulation box), and Bscale is a magnetic field scaling factor [Gauss]. We define a fiducial electron plasma beta (βe) at the jet base in each of our PIC simulations to which we scale the magnetic field strength:

(9)

The value of βe has been chosen to be indicative of a highly-magnetized synchrotron emitting relativistic plasma. By combining Eqs. (3), (8), and (9) we arrive at the following magnetic field scale factor:

(10)

again having set ne ≃ 101 cm−3 and γe ≃ 15. The above scale factor results in magnetic field strengths of B ≃ 100 Gauss at the base of each jet (see Fig. C.1), which in combination with our fiducial electron number density of ne ≃ 101 cm−3 produces numerically tractable levels of synchrotron emission at radio/mm wavelengths within our PIC simulations. This choice of physical scaling also helps to ensure that the electron Larmor radii do not vastly exceed the scaled PIC cell sizes (see Appendix A).

The electron number density [cm−3] within any PIC computational cell is by definition:

(11)

The dimensionless electron grid value at the base of our jet models is ne grid ≃ 64 (and varies throughout each simulation box). Combining Eqs. (1), (7), and (11) results in a scaled electron number density of ne ≃ 2.3 × 10−13 cm−3 at the jet base in each of our models. This value is many orders of magnitude below the levels inferred from observational and theoretical modeling of relativistic jets (i.e., our fiducial value of ne ≃ 101 cm−3). This difference in particle number again highlights the challenge faced by PIC codes when attempting to model the plasma of parsec scale astrophysical jets.

To help mitigate the numerical limitation in our ability to simulate astrophysically plausible numbers of electrons in our PIC jet calculations, we introduce a PIC “super particle” parameter (fp) which we apply to each computational cell within our jet models: , where is a proxy for the number of ‘real’ electrons represented by our ‘simulation’ electrons. As in Eqs. (1)–(3), and (8), we formulate the following electron number density scaling relation which we apply throughout our simulations:

(12)

where:

(13)

Since Δl ≃ 6.5 × 104 cm (Eqs. (1) and (7)) and ne grid ≃ 64, in order to scale the dimensionless grid electron number densities at the jet base to our fiducial jet value of ne ≃ 101 cm−3, we need to tune our PIC ‘super particle’ parameter to fp ≃ 4.6 × 1013. This immense value highlights a numerical limit PIC codes face in simulating plasmas on astrophysical scales and also sets a numerical benchmark for future simulations with larger particle populations.

With our dimensionless PIC grid values scaled into physical units (i.e., cgs), we are now in a position to apply polarized radiative transfer via ray-tracing, in order to infer both the level and the morphology of the polarized synchrotron emission emanating from the normal plasma jet and the pair plasma jet (both of which are illustrated in Fig. 1).

thumbnail Fig. 1.

Left panel: 3D visualization of the magnetic field within the normal plasma (e− − p+) PIC jet simulation. Each vector highlights the magnetic field strength within an individual computational cell. Right panel: 3D visualization of the magnetic field within the pair plasma (e− − e+) PIC jet simulation. The same convention (normal plasma on the left and pair plasma on the right) is used in the other figures of this paper.

3. Polarized radiative transfer

Jones & O’Dell (1977a,b), Jones (1988) present solutions to the full Stokes equations of polarized radiative transfer for synchrotron emission emanating from a homogeneous and an inhomogeneous magnetized plasma containing isotropic distributions of electrons (see Appendix B for further discussion). We use these solutions to compute the levels of linear polarization (LP) and circular polarization (CP) emanating from our PIC jet simulations at radio frequencies. In particular, we solve the matrix presented in Eq. (14) along individual rays passing through our PIC simulations.

(14)

Here, Iν, Qν, Uν, and Vν denote the frequency-dependent Stokes parameters, while (κI, κQ, κU, κV) and represent, respectively, the synchrotron absorption and emission coefficients for each Stokes parameter. The terms and account for the effects of Faraday conversion and rotation (see Appendix C), respectively, within our PIC plasma. The term l denotes the path length of each ray through the computational cells of our PIC calculations. Jones & O’Dell (1977a) present an analytic solution to this matrix which we apply along each ray. The radiative transfer is carried out in the ‘co-moving’ frame of the plasma with: (i) a relativistic abberation correction being applied, cell-by-cell, to obtain the angle between each cell’s local magnetic field vector and the observer’s inclination to the jet axis, and (ii) a rotation correction being applied, cell-by-cell, to transform the linear polarization ellipse from the local co-moving frame onto the plane of the sky (see Appendix D for further discussion). Once the Stokes parameters have been generated for each sightline (i.e., for each pixel/ray in our synthetic maps) a Doppler factor is applied (which incorporates the velocity/angle dependence of the larger scale jet) to obtain the flux levels in the observer’s frame. This Doppler boosting/de-boosting is apparent upon comparison of the relative flux levels between images generated when viewing each simulation edge-on and at right angles to the jet axis. We have embedded this polarized radiative transfer scheme into the ray-tracing code RADMC-3D2. For the images presented in this paper, RADMC-3D casts 640 000 individual rays through our cartesian PIC grids forming 800 × 800 square-pixel images of the polarized emission produced by our jet models. The resultant emission maps from these ray-tracing calculations for the normal plasma and pair plasma jet simulations are presented in Figs. 3, 4, and 6.

4. Slow-light interpolation

In order to resolve (both spatially and temporally) the kinetic-scale processes occurring within a relativistic jet plasma, the time step (Δt) of the TRISTAN code is set such that:

(15)

Therefore, it takes ten numerical time steps for a light ray to traverse the length of an individual plasma cell within each jet calculation. Equation (15) has profound implications for our ray-tracing calculations which typically assume/invoke the fast-light approximation. Fast-light ray-tracing makes the assumption that the light-crossing time of the jet is far smaller than the dynamical time step of the jet simulation being imaged (i.e., that we are taking a snapshot of static jet plasma). This assumption is not in general valid for PIC simulations. The fact that the plasma is evolving in time as each ray propagates through the computational grid must be taken into account within our ray-tracing calculations. In particular, we have constructed a slow-light interpolation scheme that accounts for this effect and is illustrated in Fig. 2. In essence, our interpolation scheme builds up a hybrid computational grid composed of the stratified components of successive time steps within each jet simulation. This scheme ensures that the plane parallel rays of our ray-tracing calculations encounter the correct upstream plasma values as each ray propagates through the jet. Within our PIC jet simulations:

thumbnail Fig. 2.

Schematic representation of our PIC slow-light interpolation scheme. A ray () enters the jet and encounters successive upstream plasma cells as the ray propagates across the jet (similar to stepping through a fast moving stream). These successive encounters contribute to the emission stored in each ray and result in the final total intensity value (Iν) that exits the jet and produces one pixel in our synthetic radio maps.

– The time scale for light crossing is ∼500 simulation time steps.

– The time scale for jet plasma evolution is ∼10 simulation time steps.

– The time scale for particle acceleration is < 10 simulation time steps.

Since the time scale for light crossing ≫ than the time scale for jet plasma evolution, slow light interpolation has a noticeable effect on the resultant synchrotron emission (see Figs. 3 and 6). Our radiative transfer scheme, however, is unable to incorporate the smallest particle acceleration time scales within our emission calculations, since we implicitly assume that the plasma properties of a given cell remain constant during the 10 simulation time steps it takes a light ray to traverse each cell. We also emphasize that when scaled into physical units, a simulation time step within our jet models corresponds to ∼2.2−7 s. This infinitesimal value is again a reminder of the drastically smaller spatial and temporal scales probed by PIC models in comparison to RMHD jet calculations.

thumbnail Fig. 3.

Upper left panel: fast-light Stokes I image (at an observing frequency of νobs = 1 GHz) of the normal plasma (e− − p+) jet. The internal structure of the normal plasma jet’s magnetic field is illustrated in the left panel of Fig. 1. Upper right panel: corresponding fast-light Stokes I image of the pair plasma (e− − e+) jet. The internal structure of the pair plasma jet’s magnetic field is similarly illustrated in the right panel of Fig. 1. Middle left panel: rendering of the LP intensity of the normal plasma jet. White line segments indicate the electric vector position angles (EVPAs) as projected onto the plane of the sky. The effects of relativistic aberration (see Lyutikov et al. 2005) on the orientation of these EVPAs have been included in these calculations. Middle right panel: corresponding LP intensity image of the pair plasma jet. Lower left panel: stokes V image of the normal plasma jet highlighting the different regions within the jet that produce positive and negative CP. Lower right panel: corresponding Stokes V image of the pair plasma jet. All six images have been convolved with a circular Gaussian beam of FWHM 7.5 mas (shown in the lower left of each panel). The projected PIC cell size is ∼0.9 mas at a source distance of 1 AU.

5. Results

To gain a better understanding of the effects of slow-light interpolation within our ray-traced images, we first construct fast-light images for the purposes of comparison between the two ray-tracing methods. Both the normal plasma jet and the pair plasma jet simulations were run for a sufficient time span to allow both jets to propagate across the full extent of both computational grids (shown in Fig. 1). At each time step within our PIC simulations we apply the scaling relations presented in Sect. 2 to our PIC grid values in order to create ray-tracing output files. These output files consist of three dimensional arrays containing:

  • Magnetic field strength (and orientation):

  • Electron number densities: ne[i, j, k]

  • Minimum electron Lorentz factors: γmin[i, j, k]

where in general E = γemec2, and the indices [i, j, k] denote the [x, y, z] directions within our cartesian computational grids, respectively. In contrast to magnetohydrodynamic calculations where prescriptions for the electron values must be applied in order to infer these quantities from the thermal fluid variables (e.g., Porth et al. 2011; Fromm et al. 2016), in our PIC calculations, we are able to compute these values directly from our jet simulations.

The output files from the TRISTAN particle-in-cell code then become input files for the RADMC-3D ray-tracing code. Given the fact that the entire computational grid in each PIC calculation (when scaled into physical units) only spans roughly ∼150 km (as discussed in Sect. 2) we have arbitrarily placed each jet at a distance of 1 Astronomical Unit (AU) from the ‘observer’ (i.e., the distance to the Sun). This distance scale, in combination with our simulation box sizes, generates radio maps with an angular extent of ∼200 milli-arcseconds (mas). Clearly, 150 cubic kilometers of relativistic jet does not exist at such a close distance to the Earth, but, in the spirit of theoretical study, we proceed with these calculations to infer what: (i) the morphology and (ii) the fractional level of polarized synchrotron emission would be if we could image each of these plasma jets with an idealized interferometric array at such close proximity.

5.1. Fast-light images

The results of our fast-light ray-tracing calculations for both simulations are presented in Figs. 3 and 4. We explore the nature of the polarized synchrotron emission when viewing each jet at right angles to the jet axis (i.e., θobs = 90°, shown in Fig. 3) and when viewing each jet edge-on to the jet axis (i.e., θobs = 0°, shown in Fig. 4). For both sets of images we have performed our ray-tracing calculations through the last time step of each simulation (illustrated in Fig. 1).

thumbnail Fig. 4.

Upper left panel: fast-light Stokes I image (at an observing frequency of νobs = 230 GHz and θobs = 0°) of the normal plasma (e− − p+) jet. Upper right panel: corresponding Stokes I image of the pair plasma (e− − e+) jet. Middle left panel: rendering of the LP intensity of the normal plasma jet. White line segments indicate the electric vector position angles (EVPAs) as projected onto the plane of the sky. Middle right panel: corresponding LP intensity image of the pair plasma jet. Integrated levels of fractional linear polarization () are listed to the lower right in each panel. Lower left panel: stokes V image of the normal plasma jet. Lower right panel: stokes V image of the pair plasma jet. The red dots in the lower panels highlight individual sightlines through each simulation which are illustrated in Fig. 5. Integrated levels of fractional circular polarization () are listed to the lower right in each panel. All images have been convolved with a lower resolution circular Gaussian beam of FWHM 3.5 mas for comparison to Fig. 3.

We point out that in the case of a ‘pure’ pair plasma, Stokes V would be exactly zero, since the equal numbers of electrons and positrons would cancel out their respective contributions to the circularly polarized synchrotron emission. Similar to the methods presented in Wardle & Homan (2003), Homan et al. (2009), and more recently Anantua et al. (2020), we parametrize within each plasma cell the proton-to-position ratio (). This ratio is incorporated into the radiative transfer coefficients, including the CP specific terms: , κV, and (see Appendix C; see also MacDonald & Marscher 2018). We have initially fixed this ratio to r+ = 100 (i.e., 100 protons for every positron) in the normal plasma jet and r+ = 0.01 (i.e., 100 positrons for every proton) in the pair plasma jet. We plan on running a larger set of PIC plasma simulations in the future in which we explore a wider range of jet plasma compositions than the two cases investigated here.

5.1.1. Radio jet orientation

In Fig. 3, we present ray-traced images of each jet simulation at θobs = 90° for an observing frequency of νobs = 1 GHz (see Appendix E for further discussion). In the upper panels we present Stokes I maps of both jets (normal plasma on the left and pair plasma on the right), in the middle panels we present the corresponding maps of linearly polarized intensity () for each jet with electric vector position angles overlaid in white. In the lower panels of Fig. 3 we present the corresponding Stokes V images for both simulations (again normal plasma on the left and pair plasma on the right). To mimic the effects of interferometric resolution, all images have been convolved with a circular Gaussian beam of FWHM 7.5 mas. Upon inspection of Fig. 3, one can see that the pair plasma jet exhibits a slightly more filamentary emission morphology in contrast to the normal plasma jet. This morphological difference in emission is a reflection of the distinct jet dynamics occurring within the two plasmas: the e− − e+ plasma jet is prone to larger plasma instabilities (such as the Kelvin-Helmholtz and Weibel instabilities) and is, as a result, less stable than an e− − p+ plasma jet of similar speed and magnetic field strength (see, e.g., Nishikawa et al. 2016b, 2017, 2019).

5.1.2. Blazar jet orientation

In Fig. 4, we present ray-traced images of each jet simulation at θobs = 0° for an observing frequency of νobs = 230 GHz (see Appendix E for further discussion). In the upper panels we present Stokes I images of both jets (normal plasma on the left and pair plasma on the right), in the middle panels we present LP images, and in the lower panels we present the corresponding CP images for both simulations. We have convolved these higher frequency images with a beam size similar to the lower frequency 1 GHz images for ease of comparison with Fig. 3. The effect of Doppler beaming is apparent upon comparison of the Stokes I flux levels in Figs. 3 and 4. The fractional circular polarization (mc ≡ −V/I) is minimal in each simulation (≪1%) but is many orders of magnitude larger in the e− − p+ plasma jet in comparison to the e− − e+ plasma jet (integrated values are listed to the lower right in each Stokes V image). We present spectropolarimetry of the integrated levels of fractional polarization for each jet in Appendix F.

In addition to mimicking the resolution of an interferometric array (via beam convolution), our ray-tracing algorithm can also mimic the sensitivity of an interferometric array by including a synthetic Gaussian noise floor within our ray-traced images. This is discussed further in relation to the detectability of our PIC jets in Appendix G.

5.1.3. Individual ray-properties

In Fig. 5, we illustrate the variations of the Stokes parameters (upper panels) and the Faraday rotation (τF) and conversion (τC) depths (lower panels) of both jet plasmas along the sightlines indicated by the red dots in the lower panels of Fig. 4. The Faraday depths (which are written out explicitly in Appendix C and are themselves functions of B, ne, and γmin) parameterize the ability of the plasma cells along each sightline to attenuate both the linearly and circularly polarized synchrotron emission within the jet. It is evident upon comparison of the lower panels of Fig. 5 that, for these particular sightlines, the Faraday rotation depth is larger in the e− − p+ plasma jet. Both jets, however, are optically (and Faraday) thin and the Stokes V images presented in Fig. 4 are intrinsic in origin (see the lower panels of Fig. F.1).

thumbnail Fig. 5.

Upper left panel: variations of the Stokes parameters for the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0) and then advancing through the jet plasma toward the near side of the jet (cell 240). The inset shows the variation in the fractional circular polarization along the same sightline. Upper right panel: corresponding ray properties along the pair plasma (e− − e+) jet sightline indicated by the red dot in the lower right panel of Fig. 4. Lower left panel: variations of the Faraday rotation (τF) and conversion (τC) depths along the normal plasma jet sightline. Lower right panel: corresponding Faraday depths along the pair plasma jet sightline.

5.2. Slow-light images

The results of our slow-light interpolated ray-tracing calculations are presented in Fig. 6. The slow-light interpolation results in an averaging of the emission along each sightline. In contrast to the fast-light images, the morphology of the resultant emission becomes blurred and the distinct morphologies present between the normal plasma and pair plasma jet compositions (evident in the fast-light images shown in Fig. 3) are less pronounced in the corresponding slow-light images illustrated in Fig. 6.

thumbnail Fig. 6.

Upper left panel: slow-light interpolated Stokes I image (at an observing frequency of νobs = 1 GHz) of the normal plasma (e− − p+) jet. Upper right panel: corresponding slow-light interpolated Stokes I image of the pair plasma (e− − e+) jet. Hybrid computational grids (constructed using the scheme illustrated in Fig. 2) were used in each ray-tracing calculation. Middle left panel: rendering of the LP intensity of the normal plasma jet. White line segments indicate the electric vector position angles (EVPAs) as projected onto the plane of the sky. The effects of relativistic aberration (see Lyutikov et al. 2005) on the orientation of these EVPAs have been included in these calculations. Middle right panel: corresponding LP intensity image of the pair plasma jet. Lower left panel: stokes V image of the normal plasma jet highlighting the different regions within the jet that produce positive and negative CP. Lower right panel: corresponding Stokes V image of the pair plasma jet. All six images have been convolved with a circular Gaussian beam of FWHM 7.5 mas (shown in the lower left of each panel).

In particular, the hybrid computational grids (constructed via the algorithm illustrated in Fig. 2), through which we ray-trace, are composed of ∼50 distinct/stratified jet epochs (each separated in time by 10 code time steps). This number of epochs corresponds to the number of distinct plasma cells that a plane-parallel ray encounters when each jet is imaged at right angles to the jet axis. As discussed in Sect. 4, as each ray propagates across the jet, it encounters newer upstream values of jet plasma due to the finite light crossing time across each PIC plasma cell (i.e., ten time steps – see Eq. (15)). Further computational time (i.e., more epochs) is required in order to produce slow-light images of these simulations when each jet is viewed edge-on (i.e., with the orientation of a blazar).

We finally point out that features/blobs in our slow-light images, in contrast to the fast-light images, do not necessarily map to individual plasma structures within the jet flow and are instead a mixture of emission from multiple plasma components along various sightlines through the jet.

6. Summary and conclusions

We have carried out full Stokes polarized radiative transfer calculations (via ray-tracing) through 3D relativistic PIC simulations of a normal plasma (e− − p+) and of a pair plasma (e− − e+) jet. We generate two sets of images of each jet simulation, one in which we invoke the fast-light approximation (see Fig. 3), and the other in which we implement slow-light interpolation (see Fig. 6). It is clear, upon comparison of Figs. 3 and 6, that slow-light interpolation has a discernible effect on the emission emanating from within each jet. In particular:

  • The finite light-crossing times through our relativistic jet simulations results in ‘blending’ of various plasma emission features along each sightline.

It is also clear that there are differences both in the morphology and in the fractional level of polarization emanating from the two jet plasma compositions. Specifically:

  • The pair plasma jet exhibits a more filamentary synchrotron emission morphology in comparison with the normal plasma jet.

  • The normal plasma jet has a much larger value of integrated fractional circular polarization in comparison with the pair plasma jet.

These differences are a reflection of the distinct plasma dynamics occurring within each jet: the pair plasma jet is prone to larger kinetic instabilities (e.g., Kelvin Helmholtz and Weibel instabilities) in comparison with the normal plasma jet and is, as a result, less stable, resulting in a more filamentary emission structure. The lower levels of circular polarization in the pair plasma jet are in keeping with the synchrotron theory we have incorporated into our ray-tracing calculations.

We emphasize that the applicability of our calculations to jets on parsec scales is quite limited. In particular:

  • Both jet simulations, when scaled into physical units, are only hundreds of kilometers in extent and are extremely tenuous in nature.

We present these calculations, despite this limitation in scale, for two main purposes: (i) to provide a point of comparison for future polarimetric imaging of larger-scale PIC simulations, and (ii) to emphasize that kinetic scale jet dynamics can produce distinct morphologies in the jet’s synchrotron radiation. These calculations reveal how the micro physics of the jet affect the macro emission we detect in the radio. This relationship between the micro and the macro is not commonly addressed in most relativistic jet simulations. Refinement of: (i) the jet injection scheme, (ii) the grid size, and (iii) the jet particle content are planned for future analysis.

This work clearly emphasizes the dire need for vastly larger computational grid sizes and particle populations when attempting to model kinetic scale effects in relativistic plasmas on astrophysical length scales. Hybrid techniques (e.g., Mignone et al. 2018; Vaidya et al. 2018; Davelaar et al. 2019; Parfrey et al. 2019; Bacchini et al. 2020) will be crucial for future jet simulations. We also highlight the care that must be taken when interpreting features in interferometric radio maps of relativistic jets, especially if the light crossing time of the jet exceeds the relevant dynamical timescales of the emitting plasma.

Acknowledgments

We are grateful to I. Myserlis for fruitful discussions regarding the intricacies of radiative transfer, E. Ros for invaluable feedback regarding VLBI, and to the anonymous referee for a thorough review of this manuscript. We also acknowledge J. L. Gómez for the initial idea to pursue this line of research. K.-I. Nishikawa has been supported by NASA research grants NNG05GK73G and NNX13AP14G. The PIC simulations were performed on Pleiades at the NASA Advanced Supercomputing Division (NAS) as well as on Gordon and Comet at The San Diego Supercomputer Center (SDSC), and on Bridges at the Pittsburgh Supercomputing Center (PSC), which are supported by the NSF. We have made use of elements of radmc3dPy: a python package/user-interface to the RADMC-3D code developed by A. Juhasz.

References

  1. Alves, E. P., Zrake, J., & Fiuza, F. 2018, Phys. Rev. Lett., 121, 245101 [NASA ADS] [CrossRef] [Google Scholar]
  2. Anantua, R., Emami, R., Loeb, A., & Chael, A. 2020, ApJ, 896, 30 [CrossRef] [Google Scholar]
  3. Bacchini, F., Ripperda, B., Philippov, A. A., et al. 2020, ApJS, 251, 10 [NASA ADS] [CrossRef] [Google Scholar]
  4. Birdsall, C. K., & Langdon, A. B. 1995, Plasma Physics via Computer Simulation, 2nd edn. (McGraw-Hill) [Google Scholar]
  5. Björnsson, C.-I. 2019, ApJ, 873, 55 [Google Scholar]
  6. Björnsson, C.-I. 2020, ApJ, 894, 47 [Google Scholar]
  7. Blandford, R. D., & Payne, D. G. 1982, MNRAS, 199, 883 [NASA ADS] [CrossRef] [Google Scholar]
  8. Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433 [Google Scholar]
  9. Buneman, O. 1993, in TRISTAN The 3-D Electromagnetic Particle Code, in Computer Space Plasma Physics, Simulation Techniques and Software, eds. H. Matsumoto, & Y. Omura (Tokyo: Terra Scientifc), 67 [Google Scholar]
  10. Croston, J. H., Ineson, J., & Hardcastle, M. J. 2018, MNRAS, 476, 1614 [Google Scholar]
  11. Davelaar, J., Olivares, H., Porth, O., et al. 2019, A&A, 632, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Davelaar, J., Philippov, A. A., Bromberg, O., & Singh, C. B. 2020, ApJ, 896, L31 [NASA ADS] [CrossRef] [Google Scholar]
  13. Dexter, J. 2016, MNRAS, 462, 115 [NASA ADS] [CrossRef] [Google Scholar]
  14. Emami, R., Anantua, R., Chael, A. A., & Loeb, A. 2021, ApJ, submitted [arXiv:2101.05327] [Google Scholar]
  15. Enßlin, T. A., Hutschenreuter, S., & Gopal-Krishna, 2019, JCAP, 2019, 035 [Google Scholar]
  16. Fan, X.-L., Wu, Q., & Liao, N.-H. 2018, ApJ, 861, 97 [NASA ADS] [CrossRef] [Google Scholar]
  17. Fromm, C. M., Perucho, M., Mimica, P., et al. 2016, A&A, 588, A101 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Fromm, C. M., Younsi, Z., Baczko, A., et al. 2019, A&A, 629, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Fuentes, A., Gómez, J. L., Martí, J. M., et al. 2018, ApJ, 860, 121 [Google Scholar]
  20. Guo, F., Li, X., Daughton, W., et al. 2020, ApJ, submitted [arXiv: 2008.02743] [Google Scholar]
  21. Homan, D. C., Lister, M. L., Aller, H. D., et al. 2009, ApJ, 696, 328 [NASA ADS] [CrossRef] [Google Scholar]
  22. Hosking, D. N., & Sironi, L. 2020, ApJ, 900, L23 [NASA ADS] [CrossRef] [Google Scholar]
  23. Janett, G. 2019, A&A, 622, A162 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Jones, T. W. 1988, ApJ, 332, 678 [NASA ADS] [CrossRef] [Google Scholar]
  25. Jones, T. W., & O’Dell, S. L. 1977a, ApJ, 214, 522 [NASA ADS] [CrossRef] [Google Scholar]
  26. Jones, T. W., & O’Dell, S. L. 1977b, ApJ, 215, 236 [NASA ADS] [CrossRef] [Google Scholar]
  27. Kim, J.-Y., Krichbaum, T. B., Marscher, A. P., et al. 2019, A&A, 622, A196 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Lyutikov, M., Pariev, V. I., & Gabuzda, D. C. 2005, MNRAS, 360, 869 [Google Scholar]
  29. MacDonald, N. R., & Marscher, A. P. 2018, ApJ, 862, 58 [NASA ADS] [CrossRef] [Google Scholar]
  30. Marscher, A. P. 2014, ApJ, 780, 87 [NASA ADS] [CrossRef] [Google Scholar]
  31. Mignone, A., Bodo, G., Vaidya, B., et al. 2018, ApJ, 859, 13 [NASA ADS] [CrossRef] [Google Scholar]
  32. Mizuno, Y., Gómez, J. L., & Nishikawa, K.-I. 2020, ApJ, 809, 38 [Google Scholar]
  33. Moscibrodzka, M., & Gammie, C. F. 2018, MNRAS, 475, 43 [Google Scholar]
  34. Mościbrodzka, M., Janiuk, A., & De Laurentis, M. 2021, MNRAS, submitted [arXiv: 2103.00267] [Google Scholar]
  35. Mukherjee, D., Bodo, G., Mignone, A., Rossi, P., & Vaidya, B. 2020, MNRAS, 499, 681 [Google Scholar]
  36. Mukherjee, D., Bodo, G., Rossi, P., Mignone, A., & Vaidya, B. 2021, MNRAS, 505, 2267 [NASA ADS] [CrossRef] [Google Scholar]
  37. Myserlis, I., Komossa, S., Angelakis, E., et al. 2018a, A&A, 619, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Myserlis, I., Angelakis, E., Kraus, A., et al. 2018b, A&A, 609, A68 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. Niemiec, J., Pohl, M., Stroman, T., & Nishikawa, K.-I. 2008, ApJ, 684, 1184 [NASA ADS] [Google Scholar]
  40. Nishikawa, K.-I., Hardee, P. E., Duţan, I., et al. 2014, ApJ, 793, 60 [NASA ADS] [CrossRef] [Google Scholar]
  41. Nishikawa, K.-I., Frederiksen, J. T., Nordlund, A., et al. 2016a, ApJ, 820, 94 [NASA ADS] [CrossRef] [Google Scholar]
  42. Nishikawa, K.-I., Mizuno, Y., Niemiec, J., et al. 2016b, Galaxies, 4, 38 [NASA ADS] [CrossRef] [Google Scholar]
  43. Nishikawa, K.-I., Mizuno, Y., Gómez, J. L., et al. 2017, Galaxies, 5, 58 [NASA ADS] [CrossRef] [Google Scholar]
  44. Nishikawa, K.-I., Mizuno, Y., Gómez, J. L., et al. 2019, Galaxies, 7, 29 [NASA ADS] [CrossRef] [Google Scholar]
  45. Nishikawa, K.-I., Mizuno, Y., Gómez, J. L., et al. 2020, MNRAS, 493, 2652 [CrossRef] [Google Scholar]
  46. Nishikawa, K. I., Duţan, I., Köhn, C., et al. 2021, Liv. Rev. Comput. Astrophys., 7, 1 [Google Scholar]
  47. O’Sullivan, S. P., McClure-Griffiths, N. M., Feain, I. J., et al. 2013, MNRAS, 435, 311 [NASA ADS] [CrossRef] [Google Scholar]
  48. Pacholczyk, A. G. 1970, Radio Astrophysics. Nonthermal Processes in Galactic and Extragalactic Sources (W. H. Freeman) [Google Scholar]
  49. Parfrey, K., Philippov, A., & Cerutti, B. 2019, Phys. Rev. Lett., 122, 035101 [Google Scholar]
  50. Petropoulou, M., Sironi, L., Spitkovsky, A., & Giannios, D. 2019, ApJ, 880, 37 [Google Scholar]
  51. Porth, O., Fendt, C., Meliani, Z., et al. 2011, ApJ, 737, 42 [Google Scholar]
  52. Ricarte, A., Qiu, R., & Narayan, R. 2021, MNRAS, 505, 523 [NASA ADS] [CrossRef] [Google Scholar]
  53. Ruszkowski, M., & Begelman, M. C. 2002, ApJ, 573, 485 [NASA ADS] [CrossRef] [Google Scholar]
  54. Sikora, M., Nalewajko, K., & Madejski, G. M. 2020, MNRAS, 499, 3749 [NASA ADS] [CrossRef] [Google Scholar]
  55. Sironi, L., & Narayan, R. 2015, ApJ, 800, 88 [NASA ADS] [CrossRef] [Google Scholar]
  56. Sironi, L., Petropoulou, M., & Giannios, D. 2015, MNRAS, 450, 183 [Google Scholar]
  57. Sironi, L., Rowan, M. E., & Narayan, R. 2021, ApJ, 907, L44 [CrossRef] [Google Scholar]
  58. Tchekhovskoy, A., & Bromberg, O. 2016, MNRAS, 461, 46 [Google Scholar]
  59. Thum, C., Agudo, I., Molina, S. N., et al. 2018, MNRAS, 473, 2506 [NASA ADS] [CrossRef] [Google Scholar]
  60. Vaidya, B., Mignone, A., Bodo, G., et al. 2018, ApJ, 865, 144 [NASA ADS] [CrossRef] [Google Scholar]
  61. Wardle, J. F. C., & Homan, D. C. 2003, Ap&SS, 288, 143 [NASA ADS] [CrossRef] [Google Scholar]
  62. Zhang, H., Li, X., Guo, F., & Giannios, D. 2018, ApJ, 862, L25 [NASA ADS] [CrossRef] [Google Scholar]
  63. Zhang, H., Li, X., Dimitrios, G., et al. 2020, ApJ, 901, 149 [CrossRef] [Google Scholar]

Appendix A: Comparison of electron Larmor radii to the PIC grid cell sizes

Due to the very small length scales within our PIC simulations, care must be taken in the implementation of the physical scaling relations (Sect. 2) to ensure that the individual electron Larmor radii do not vastly exceed the plasma grid cell sizes. This is a necessary criterion for modeling isotropic synchrotron emission in the radiative transfer (Eq. (14)) that we apply along each ray in our transfer calculations. The radiative transfer is extremely sensitive to the lower cutoff in the electron power-law distribution (γmin). For a power-law distribution of synchrotron emitting electrons: n(γ)∝γ−s, the synchrotron emissivity is a function of the integral . Provided that s > 1 (s ≃ 2.3 in our calculations) and γmax ≫ γmin, the contribution from the high-energy end of the power-law is minimal. Values of γmin along sightlines through the e− − p+ and e− − e+ jet plasmas are illustrated in the lower panels of Fig. C.1. In particular, at the base of the normal plasma jet: γmin ≃ 15 and B ≃ 100 Gauss. This implies that the low-energy electrons in this region will emit synchrotron radiation at a critical frequency (νc) of:

(A.1)

which corresponds to Doppler beamed GHz emission in the observer’s frame when a bulk Lorentz factor Γ = 15 jet is oriented edge-on to our line of sight (i.e., with the orientation of a blazar). The corresponding Larmor radii of these low-energy electrons is:

(A.2)

The above Larmor radius is (by design) smaller than the scaled PIC plasma cell sizes (Δlscale ≃ 6.5 × 104 cm, Eq. (7)). Our choice of physical scaling (i.e., the fiducial values of ne ≃ 101 cm−3; Eq. (6) and βe ≃ 10−3; Eq. (9)) is motivated by: (i) the ability of our PIC simulations to produce numerically tractable levels of synchrotron radiation, and (ii) to ensure that rgyro(γmin) ≤ Δlscale. Equation (A.2), however, parameterizes the gyration of electrons about a single ordered magnetic field line. The magnetic fields within our PIC simulations, in contrast, are highly turbulent in nature. A more applicable expression for the electron Larmor radii within our PIC simulations can be recast in terms of the plasma beta (βe) and the electron plasma frequency (ωpe):

(A.3)

(see, e.g., Sironi & Narayan 2015). Equation (A.3) can be derived using the ratios of the electron (e) to proton (p) cyclotron frequencies: , the plasma frequencies: , and the ratio of the Alfvén speed to the speed of light: vA/c = ωcp/ωpp, assuming relativistic electrons and nonrelativistic protons. Combining these ratios allows one to re-formulate the electron Larmor radius (rgyro ≡ ve/ωce) into the form of Eq. (A.3) in the relativistic limit where: . In Fig. A.1 we plot individual electron Larmor radii computed using Eq. (A.3) (for the normal plasma jet – left panel, and the pair plasma jet – right panel) along rays highlighted by the red dots in the lower panels of Fig. 4. The red dashed horizontal line in each panel highlights the scaled PIC grid cell size (Δlscale). These ray profiles verify that the Larmor radii of our simulated PIC electrons are largely contained within our numerical cells. The regions where the electron Larmor radii exceed the plasma cell size occur where the magnetic field drops precipitously (see the top panels of Fig. C.1). These zones contribute minimally to the integrated synchrotron emission along each ray since the synchrotron emissivity is ∝B2.

thumbnail Fig. A.1.

Left panel: variation of the electron Larmor radii (rgyro) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0), and then advancing through the jet plasma toward the near side of the jet (cell 240). The red dashed line in each panel demarcates the scaled PIC grid cell size (Δlscale). Right panel: corresponding variation of the pair plasma (e− − e+) jet Larmor radii along the sightline indicated by the red dot in the lower right panel of Fig. 4.

Appendix B: Electron phase-space and energy distributions

The synchrotron emission/absorption coefficients (contained in Eq. (14)) assume that the underlying electron distribution within the jet is isotropic. While we have attempted to scale our PIC simulations in a manner that ensures the individual electron Larmor radii are largely contained within our grid cells (see Appendix A), this scaling does not guarantee that the particle distribution is indeed isotropic in nature. To explore this point further, we have generated 2D phase-space distribution plots of the electrons contained within the jet regions we generate synchrotron emission from (see the left and middle panels of Fig. B.1). Clearly, anisotropy is present within the underlying particle distribution used in our emission calculations (especially along the jet axis – Fig. B.1 left panel). Despite this anisotropy, we compute the synchrotron emission using the fully isotropic emission coefficients of Eq. (14) in order to gain a ‘first approximation’ of the synchrotron emission emanating from our PIC simulations.

thumbnail Fig. B.1.

Left panel: an electron phase-space distribution plot along the jet axis – individual electrons are plotted as points. The underlying color scheme highlights the corresponding probability density function (PDF). Middle panel: an electron phase-space distribution plot perpendicular to the jet axis. Right panel: electron energy spectrum for the particles contained within our PIC simulations. The high-energy tail exhibits an exponential cut-off in energy (highlighted in red).

The synchrotron emission is calculated in a post-process fashion. We first obtain (within each cell) the local 3D fluid velocity and magnetic field components which are averaged over 27 neighboring cells with weighting to help ensure good statistical properties within each cell. These values are then used to compute the cell’s ‘local’ magnetic field strength and ‘local’ electron Lorentz factor. In addition, the number of particles are counted and scaled (via Eq. (12)) to obtain a ‘local’ electron number density.

As a further simplification in our emission calculations, we have arbitrarily fixed a constant power-law index (s = 2.3) within each plasma cell. Future refinement of our algorithm will be needed in order to connect the spectral indices used in our synchrotron emission calculations more closely to the ‘local’ particle energy spectrum within each cell. To explore this point further, we have generated a particle energy spectrum of the electrons contained within the jet regions we compute synchrotron emission from (see the right panel of Fig. B.1). Clearly, the particle spectrum (especially the high-energy tail) is much steeper than the power-law index of s = 2.3 used in our emission calculations. This remains a limitation/discrepancy between our PIC models and our synchrotron emission calculations.

Appendix C: Faraday rotation and conversion depths

The Faraday rotation (τF) and conversion (τC) depths within an individual plasma cell may be written as: and (see Wardle & Homan 2003; Homan et al. 2009; O’Sullivan et al. 2013). The term τ denotes the synchrotron opacity of the plasma. The terms and denote normalized plasma absorption coefficients: (i.e., the κ’s in the matrix shown in Sect. 3) and may be written as:

(C.1)

(C.2)

where and are physical constants (of order unity) that are tabulated in Jones & O’Dell (1977a) and depend on the value of the spectral index α. The quantity , where νB⊥ is given by:

(C.3)

The angle ϑ denotes the angle that each sightline within our ray-tracing calculations makes with respect to the local magnetic field vector in the co-moving frame of each plasma cell. Under the assumption of charge neutrality (), the term in Equation C.2 may be re-written as (1 + 2/r+)−1, where and parameterizes the plasma composition of each cell (see Sect. 5.1). In all of the computations presented in this paper, we have assumed a constant optically thin spectral index of α = 0.65. In future, we plan on refining our algorithm to vary α cell-to-cell and compute it directly from the local electron power-law index (s): α = (s − 1)/2. The synchrotron opacity (τ) of an individual plasma cell may be written as:

(C.4)

where l is the path length of the ray through the plasma cell (computed using RADMC-3D) and where ne is the electron number density within the plasma cell. The length scale of the PIC plasma cells in the simulations presented in this paper is ≃6.5 × 104 cm (see Eq. (7)). The parameter κα is a physical constant (of order unity) that is tabulated in Jones & O’Dell (1977a) and depends on the value of the spectral index α. The term g(ϑ), represents the electron pitch angle distribution. It is evident, upon inspection of equations: (C.1)–(C.4), that τF, τC, and τ are together themselves functions of each plasma cell’s: B, ne, and γmin. The ray profiles (similar to those presented in Fig. 5) of these three variables (for the sightlines indicated by the red dots in the lower panels of Fig. 4) are illustrated in Fig. C.1. Upon inspection of these panels one can see that our choice of physical scaling (i.e., Sect. 2) results in magnetic field strengths that range from B ≃ 10−2 − 100 Gauss and electron number densities that range from ne ≃ 100 − 102 cm−3. These values are in rough agreement with the plasma conditions that have been inferred in blazar jets from theoretical modeling of shock acceleration and turbulence (see, e.g., Marscher 2014).

thumbnail Fig. C.1.

Upper left panel: variation of the scaled magnetic field strength (B) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0), and then advancing through the jet plasma toward the near side of the jet (cell 240). Upper right panel: corresponding variation of the scaled magnetic field strength along the pair plasma (e− − e+) jet sightline indicated by the red dot in the lower right panel of Fig. 4. Middle left panel: variation of the scaled electron number density (ne) along the normal plasma jet sightline. Middle right panel: corresponding variation of the scaled electron number density along the pair plasma jet sightline. Lower left panel: variation of the minimum electron Lorentz factor (γmin) along the normal plasma jet sightline. Lower right panel: corresponding variation of the minimum electron Lorentz factor along the pair plasma jet sightline. These panels highlight the physical scaling (i.e., Sect. 2) that we have applied to our dimensionless PIC grid values.

We present here an alternate physical scaling method (commonly used in magnetohydrodynamic simulations) which generates a similar magnetic field scale factor for comparison to our PIC approach (Sect. 2). In an RMHD simulation one typically defines the following three scale factors:

(C.5)

(C.6)

(C.7)

(see the PLUTO3 code manual for further discussion). All other physical scaling quantities can be derived from these three unit values. In particular, the magnetic field strength scale factor (in cgs) may be written as:

(C.8)

As discussed in Sect. 2, we set vo = 3 × 1010 cm s−1. Again, under the assumption of charge neutrality (), and invoking the normal plasma case, in which , one may make the assumption that: . While the electrons dominate the nonthermal synchrotron emission, the protons (in contrast) dominate the thermal fluid dynamics. It follows that: ρo = nomp, where no ≃ 101 cm−3 (i.e., setting the proton number density equal to the electron value from Sect. 2). This results in ρo ≃ 1.7 × 10−23 g cm−3, having used mp = 1.6726231 × 10−23g cm−3. Inserting vo and ρo into Eq. (C.8) yields a magnetic field strength scaling factor of Bo ≃ 0.4 Gauss which is in agreement with the scaling value we obtain in Sect. 2 (Eq. (10)) after invoking a fiducial jet plasma beta of βe ≃ 10−3.

Appendix D: Radiative transfer angles

The synchrotron emission/absorption coefficients (contained in Eq. (14)) are functions of the angle (ϑ) each ray makes with respect to the local magnetic field vector in the co-moving frame of each plasma cell. When computing the Stokes parameters, we execute the following numerical steps (cell-by-cell) to ensure a proper treatment of the angles involved in the radiative transfer:

  • (i)

    We first compute the angle our line-of-sight makes with respect to each cell’s local magnetic field vector accounting for relativistic aberration due to the bulk flow of the jet. In particular, the following Lorentz transformation maps our sightline unit vector in the observer’s frame to the corresponding unit vector in the co-moving frame of the jet plasma:

    (D.1)

    (see Lyutikov et al. 2005), where Γ is the bulk Lorentz factor of the jet flow and . We make the simplifying assumption that there exists a common ‘jet frame’; namely, that all cells within our simulations propagate in a laminar fashion along the jet-axis and share a common bulk Lorentz factor (Γ = 15 in this case). This simplification makes the radiative transfer more tractable and permits us to compute the synchrotron emission along plane parallel rays (depicted in Fig. 2). Future refinement of our algorithm will be required to account properly for cell-to-cell variations in the relativistic motions of the plasma within the jet. As illustrated in Fig. D.1, we orient the jet along the -axis and therefore β = { 0,  0,  β }. From Fig. D.1 it also follows that: , where θobs is the angle the observer’s line-of-sight makes with respect to the jet axis. Recalling the definition of the relativistic Doppler boosting factor: , Eq. (D.1) can be rewritten in component form as:

    (D.2)

    thumbnail Fig. D.1.

    A schematic representation of the angles involved in our radiative transfer calculations. The observer’s sightline (wave vector ) makes an angle θobs with respect to the jet/ axis. The magnetic field unit vector, , of an individual plasma cell is depicted in red. The radiative transfer is carried out in a cell-specific ‘Stokes U = 0’ linear polarization basis (depicted in blue). The cell-specific linear polarization basis is then rotated by an angle ϕrot onto a generalized observer’s plane (shown to the left).

    after combining terms and simplifying. The local magnetic field unit vector is given by . It then follows that:

    (D.3)

  • (ii)

    With ϑ determined (cell-by-cell) the radiative transfer is then computed across each cell. These calculations are done in a cell specific ‘Stokes U = 0’ linear polarization basis. The analytic solution used in our computations is presented in Jones & O’Dell (1977a) (Appendix C) and summarized in MacDonald & Marscher (2018) (Appendix A).

  • (iii)

    Once the radiative transfer across a cell is complete, we then rotate the cell specific linear polarization basis onto a generalized observer’s plane (illustrated in Fig. D.1). In particular, the angle of this rotation (ϕrot) is:

    (D.4)

    and is equal to the angle the projected magnetic field makes with respect to the x-axis (i.e., west) on the observer’s plane. This rotation is then applied to the Stokes parameters in the following manner:

    (D.5)

  • (iv)

    These ‘generalized’ Stokes parameters are then recorded and passed on to the next cell via our slow-light interpolation scheme at which point we rotate back into the jet plasma frame and repeat steps (i–iii) for the next plasma cell along each sightline.

Appendix E: Radiative transfer limits

Care must be taken when implementing the analytic solution to the polarized radiative transfer equation presented in Eq. (14). As discussed in Janett (2019), there are many numerical challenges one faces when attempting to implement polarized radiative transfer through plasma simulations. The analytic solution to the matrix in Eq. (14) (which is presented in Appendix B of Jones & O’Dell 1977a) is laced with hyperbolic sine and cosine functions whose arguments consist of the expression χτ. Here, χ is a term (of order unity) related to the Faraday rotation and conversion depths of the plasma (also presented in Appendix B of Jones & O’Dell 1977a) and τ is the synchrotron opacity of the plasma (see Appendix C of this manuscript). Given the extremely small cell sizes within our PIC simulations, the path length l (see Eq. (C.4)) of each ray through our plasma cells is many orders of magnitude smaller than typical path lengths through RMHD grids. This results in extremely small synchrotron opacities along the sightlines through our PIC jets (see the left panel of Fig. E.1). Fortran’s hyperbolic cosh and sinh functions do not go exactly to 1.0 and 0.0 respectively, when there arguments (i.e., χτ) fall below 1.0 × 10−14 (i.e., the red dashed line in the left panel of Fig. E.1). Performing radiative transfer beyond this numerical limit can introduce erroneous values in the Stokes parameters. As discussed in MacDonald & Marscher (2018; see their Appendix A), as τ → 0 within a given plasma cell the Stokes parameters should analytically equal the values from the previous plasma cell along each ray’s path (i.e., ). We enforce this criterion directly in our PIC ray-tracing calculations when a ray encounters a plasma cell in which χτ ≤ 1.0 × 10−14. Also, when viewing our PIC simulations at right angles (i.e., θobs = 90°) sightlines through the jet becomes increasingly optically thin because each ray intersects less jet plasma. Since τ ∝ ν−(α + 5/2), where α = 0.65, observing at lower frequencies helps to increase the synchrotron opacity and mitigate this numerical limitation in the radiative transfer. This is why, when imaging our plasma jets at right angles, we tune νobs = 1 GHz.

thumbnail Fig. E.1.

Left panel: variation of the term χτ (synchrotron opacity) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The red dashed line indicates a numerical limit below which we enforce the criterion directly in our PIC ray-tracing calculations. Right panel: variation of the minimum frequency (νmin) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The red dashed line indicates the emitted frequency νem in the plasma frame corresponding to νobs = 230 GHz. In both panels, the radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0), and then advancing through the plasma toward the near side of the jet (cell 240).

As stated in Jones & O’Dell (1977a), the frequency dependent synchrotron emission/absorption coefficients contained in Eq. (14) are only valid provided that the frequency of emission νem (in the co-moving frame of the plasma) exceeds the minimum frequency νmin of the plasma (defined in Appendix C). We check, cell-by-cell, that this criterion is met in all of our ray-tracing calculations (see the right panel of Fig. E.1). In particular, νobs = δνem, where δ is the Doppler factor and is given by δ ≡ ( Γ(1 − β cos θobs) )−1. Here Γ is the bulk Lorentz factor of the jet, β ≡ v/c = (1 − Γ−2)1/2, and θobs is the angle between the jet axis and our line-of-sight. The minimum frequency within each plasma cell is a function of: (i) the cell’s minimum electron Lorentz factor, (ii) the local magnetic field strength, and (iii) the angle our line-of-sight makes with respect to each plasma cell’s local magnetic field vector, νmin(γmin, B, ϑ). Within our PIC simulations, νmin typically ranges in value from ∼107 − 109 Hz. In order to ensure the validity of the radiative transfer (i.e., νmin < νem) when imaging our plasma jets at θobs = 0°, we tune νobs = 230 GHz (which corresponds to an emitted frequency of νem ≃ 7.7 × 109 Hz in the rest frame of the plasma and is highlighted by the red dashed lined in the right panel of Fig. E.1).

We finally highlight an argument presented in Björnsson (2019) and Björnsson (2020), that questions the applicability of the analytic expressions we use to model the polarization of relativistic jets (namely, by treating our PIC cells as piecewise homogeneous regions of magnetized plasma). It is argued that instead of performing the radiative transfer on the Stokes parameters (i.e., Eq. (14)) one should instead implement the method of characteristic waves and carry out the polarized radiative transfer on the electric and magnetic fields within the plasma (computing the Stokes parameters only as a final step). While both methods should, in principle, be interchangeable, the contributing terms to CP are more easily identifiable in the characteristic wave method. In future, we plan on implementing the characteristic wave method of polarized radiative transfer in our ray-tracing algorithm for comparison to the calculations presented here. We also plan on implementing an integration scheme (similar to the methods presented in Ruszkowski & Begelman 2002; Dexter 2016; Moscibrodzka & Gammie 2018) to solve Eq. (14) numerically, thus providing a further comparison to the analytic expressions of Jones & O’Dell (1977a) used in this manuscript.

Appendix F: Particle-in-cell spectropolarimetry

As a further comparison between the radiative properties of the normal plasma jet and the pair plasma jet we generate spectra of each jet’s: (i) integrated total intensity – I, (ii) integrated fractional linear polarization – ml, (iii) integrated EVPA – χ, and (iv) integrated fractional circular polarization – mc, over the frequency range νobs = 86 − 230 GHz (see Fig. F.1), when each jet is oriented edge-on to the observer as shown in Fig. 4 (i.e., like blazars). In particular, we follow the formalism of Kim et al. (2019) and compute integrated Stokes values for each simulation in the following manner:

(F.1)

(F.2)

(F.3)

(F.4)

thumbnail Fig. F.1.

Left panel: spectropolarimetry of the normal plasma jet: (top) integrated total intensity – I, (upper middle) integrated fractional linear polarization – ml, (lower middle) integrated EVPA – χ, and (bottom) integrated fractional circular polarization – mc, over the frequency range νobs = 86 − 230 GHz. Right panel: corresponding spectropolarimetry of the pair plasma jet. The solid red and dashed red lines in the lower panels highlight theoretical predictions for the frequency dependence of mc(ν) in the limits of intrinsic emission and emission dominated by Faraday conversion in the high-rotation limit, respectively.

where ∑Ii, j, ∑Qi, j, ∑Ui, j, and ∑Vi, j are summations of the map pixel values contained within the outermost Stokes I contour of each ray-traced image (conservatively set at 20% of each map’s peak value). The term Apixel (= 0.4 × 0.4 square mas) denotes the angular extent of each pixel in our ray-traced images and is computed from the ratio of the RADMC-3D pixel size (∼3 × 104 cm) to the source distance (∼1.5 × 1013 cm; 1 AU). The term Abeam (=πψmaj × ψmin/4 ln 2) denotes the angular extent of the convolving beam, where ψmaj and ψmin are the FWHM of the major and minor beam axes, respectively. For this analysis we have used a fixed circular beam of ψmaj = ψmin = 3.5 mas (illustrated in the lower left of each panel in Fig. 4). From these integrated Stokes values we compute the following polarimetric quantities:

(F.5)

(F.6)

(F.7)

(F.8)

The above quantities are then recorded for each ray-tracing calculation at individual frequencies ranging from νobs = 86 − 230 GHz (see Fig. F.1).

The normal plasma jet and the pair plasma jet are both optically thin in this frequency range. Both spectra exhibit a general trend of decreasing intensity (top panels) with increasing fractional linear polarization (upper middle panels) accompanied by decreasing fractional circular polarization (bottom panels) in agreement with synchrotron theory. An anticorrelation between linear and circular polarization across frequency has been detected in the quasar PKS B2126-158 with high-precision CP measurements made using the Australian Compact Telescope Array (ATCA; O’Sullivan et al. 2013). The solid red and dashed red lines in the lower panels highlight theoretical predictions for the frequency dependence of mc(ν) in the limits of intrinsic emission and emission dominated by Faraday conversion in the high-rotation limit, respectively (see Pacholczyk 1970; Jones & O’Dell 1977a). Clearly, the emission produced in our PIC models is intrinsic in origin due to the very small Faraday depths through the jet plasma along each sightline (see Fig. 5).

Appendix G: Interferometric array sensitivity

In addition to being able to mimic the resolution of an interferometric array (via beam convolution), our ray-tracing algorithm can also mimic the sensitivity limit of an interferometric array through the introduction of a Gaussian noise floor into the synthetic images we produce. We illustrate this effect in Fig. G.1, in which we contrast images of the same simulation epoch; one in which we have no resolution/sensitivity limit (Left panel), one in which we apply a resolution limit by convolving our resolved image with a circular Gaussian beam of FWHM 3.5 mas (Middle panel), and the other in which we additionally introduce a infinitesimally small Gaussian noise floor of ∼10−6 Jy beam−1 (Right panel). Clearly, the Doppler de-beamed jets presented in Figs. 3 and 6 would not be detectable by existing ground-based interferometric arrays. We also point out that all of the images presented in this paper assume complete uv-coverage and, as such, represent highly idealized radio images. Future refinement of our ray-tracing algorithm is required in order to include the effects of finite uv-coverage in our synthetic maps.

thumbnail Fig. G.1.

Left panel: resolved fast-light Stokes I image (at an observing frequency of νobs = 230 GHz) of the normal plasma (e− − p+) jet simulation. The orientation of the jet axis to our line-of-sight is θobs = 0°. Middle panel: fast-light Stokes I image of the normal plasma jet, but in which the image has been convolved with a circular Gaussian beam of FWHM 3.5 mas (shown in the lower left of the panel). Right panel: fast-light Stokes I image of the normal plasma jet, but in which an artificial Gaussian noise floor (set at ∼10−6 Jy beam−1) has been included to illustrate the extreme interferometric array sensitivity needed to detect our PIC jets even at such close proximity.

All Figures

thumbnail Fig. 1.

Left panel: 3D visualization of the magnetic field within the normal plasma (e− − p+) PIC jet simulation. Each vector highlights the magnetic field strength within an individual computational cell. Right panel: 3D visualization of the magnetic field within the pair plasma (e− − e+) PIC jet simulation. The same convention (normal plasma on the left and pair plasma on the right) is used in the other figures of this paper.
In the text
thumbnail Fig. 2.

Schematic representation of our PIC slow-light interpolation scheme. A ray () enters the jet and encounters successive upstream plasma cells as the ray propagates across the jet (similar to stepping through a fast moving stream). These successive encounters contribute to the emission stored in each ray and result in the final total intensity value (Iν) that exits the jet and produces one pixel in our synthetic radio maps.
In the text
thumbnail Fig. 3.

Upper left panel: fast-light Stokes I image (at an observing frequency of νobs = 1 GHz) of the normal plasma (e− − p+) jet. The internal structure of the normal plasma jet’s magnetic field is illustrated in the left panel of Fig. 1. Upper right panel: corresponding fast-light Stokes I image of the pair plasma (e− − e+) jet. The internal structure of the pair plasma jet’s magnetic field is similarly illustrated in the right panel of Fig. 1. Middle left panel: rendering of the LP intensity of the normal plasma jet. White line segments indicate the electric vector position angles (EVPAs) as projected onto the plane of the sky. The effects of relativistic aberration (see Lyutikov et al. 2005) on the orientation of these EVPAs have been included in these calculations. Middle right panel: corresponding LP intensity image of the pair plasma jet. Lower left panel: stokes V image of the normal plasma jet highlighting the different regions within the jet that produce positive and negative CP. Lower right panel: corresponding Stokes V image of the pair plasma jet. All six images have been convolved with a circular Gaussian beam of FWHM 7.5 mas (shown in the lower left of each panel). The projected PIC cell size is ∼0.9 mas at a source distance of 1 AU.
In the text
thumbnail Fig. 4.

Upper left panel: fast-light Stokes I image (at an observing frequency of νobs = 230 GHz and θobs = 0°) of the normal plasma (e− − p+) jet. Upper right panel: corresponding Stokes I image of the pair plasma (e− − e+) jet. Middle left panel: rendering of the LP intensity of the normal plasma jet. White line segments indicate the electric vector position angles (EVPAs) as projected onto the plane of the sky. Middle right panel: corresponding LP intensity image of the pair plasma jet. Integrated levels of fractional linear polarization () are listed to the lower right in each panel. Lower left panel: stokes V image of the normal plasma jet. Lower right panel: stokes V image of the pair plasma jet. The red dots in the lower panels highlight individual sightlines through each simulation which are illustrated in Fig. 5. Integrated levels of fractional circular polarization () are listed to the lower right in each panel. All images have been convolved with a lower resolution circular Gaussian beam of FWHM 3.5 mas for comparison to Fig. 3.
In the text
thumbnail Fig. 5.

Upper left panel: variations of the Stokes parameters for the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0) and then advancing through the jet plasma toward the near side of the jet (cell 240). The inset shows the variation in the fractional circular polarization along the same sightline. Upper right panel: corresponding ray properties along the pair plasma (e− − e+) jet sightline indicated by the red dot in the lower right panel of Fig. 4. Lower left panel: variations of the Faraday rotation (τF) and conversion (τC) depths along the normal plasma jet sightline. Lower right panel: corresponding Faraday depths along the pair plasma jet sightline.
In the text
thumbnail Fig. 6.

Upper left panel: slow-light interpolated Stokes I image (at an observing frequency of νobs = 1 GHz) of the normal plasma (e− − p+) jet. Upper right panel: corresponding slow-light interpolated Stokes I image of the pair plasma (e− − e+) jet. Hybrid computational grids (constructed using the scheme illustrated in Fig. 2) were used in each ray-tracing calculation. Middle left panel: rendering of the LP intensity of the normal plasma jet. White line segments indicate the electric vector position angles (EVPAs) as projected onto the plane of the sky. The effects of relativistic aberration (see Lyutikov et al. 2005) on the orientation of these EVPAs have been included in these calculations. Middle right panel: corresponding LP intensity image of the pair plasma jet. Lower left panel: stokes V image of the normal plasma jet highlighting the different regions within the jet that produce positive and negative CP. Lower right panel: corresponding Stokes V image of the pair plasma jet. All six images have been convolved with a circular Gaussian beam of FWHM 7.5 mas (shown in the lower left of each panel).
In the text
thumbnail Fig. A.1.

Left panel: variation of the electron Larmor radii (rgyro) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0), and then advancing through the jet plasma toward the near side of the jet (cell 240). The red dashed line in each panel demarcates the scaled PIC grid cell size (Δlscale). Right panel: corresponding variation of the pair plasma (e− − e+) jet Larmor radii along the sightline indicated by the red dot in the lower right panel of Fig. 4.
In the text
thumbnail Fig. B.1.

Left panel: an electron phase-space distribution plot along the jet axis – individual electrons are plotted as points. The underlying color scheme highlights the corresponding probability density function (PDF). Middle panel: an electron phase-space distribution plot perpendicular to the jet axis. Right panel: electron energy spectrum for the particles contained within our PIC simulations. The high-energy tail exhibits an exponential cut-off in energy (highlighted in red).
In the text
thumbnail Fig. C.1.

Upper left panel: variation of the scaled magnetic field strength (B) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0), and then advancing through the jet plasma toward the near side of the jet (cell 240). Upper right panel: corresponding variation of the scaled magnetic field strength along the pair plasma (e− − e+) jet sightline indicated by the red dot in the lower right panel of Fig. 4. Middle left panel: variation of the scaled electron number density (ne) along the normal plasma jet sightline. Middle right panel: corresponding variation of the scaled electron number density along the pair plasma jet sightline. Lower left panel: variation of the minimum electron Lorentz factor (γmin) along the normal plasma jet sightline. Lower right panel: corresponding variation of the minimum electron Lorentz factor along the pair plasma jet sightline. These panels highlight the physical scaling (i.e., Sect. 2) that we have applied to our dimensionless PIC grid values.
In the text
thumbnail Fig. D.1.

A schematic representation of the angles involved in our radiative transfer calculations. The observer’s sightline (wave vector ) makes an angle θobs with respect to the jet/ axis. The magnetic field unit vector, , of an individual plasma cell is depicted in red. The radiative transfer is carried out in a cell-specific ‘Stokes U = 0’ linear polarization basis (depicted in blue). The cell-specific linear polarization basis is then rotated by an angle ϕrot onto a generalized observer’s plane (shown to the left).
In the text
thumbnail Fig. E.1.

Left panel: variation of the term χτ (synchrotron opacity) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The red dashed line indicates a numerical limit below which we enforce the criterion directly in our PIC ray-tracing calculations. Right panel: variation of the minimum frequency (νmin) within the normal plasma (e− − p+) jet along the sightline indicated by the red dot in the lower left panel of Fig. 4. The red dashed line indicates the emitted frequency νem in the plasma frame corresponding to νobs = 230 GHz. In both panels, the radiative transfer progresses from left to right: starting on the far side of the jet relative to the observer (cell 0), and then advancing through the plasma toward the near side of the jet (cell 240).
In the text
thumbnail Fig. F.1.

Left panel: spectropolarimetry of the normal plasma jet: (top) integrated total intensity – I, (upper middle) integrated fractional linear polarization – ml, (lower middle) integrated EVPA – χ, and (bottom) integrated fractional circular polarization – mc, over the frequency range νobs = 86 − 230 GHz. Right panel: corresponding spectropolarimetry of the pair plasma jet. The solid red and dashed red lines in the lower panels highlight theoretical predictions for the frequency dependence of mc(ν) in the limits of intrinsic emission and emission dominated by Faraday conversion in the high-rotation limit, respectively.
In the text
thumbnail Fig. G.1.

Left panel: resolved fast-light Stokes I image (at an observing frequency of νobs = 230 GHz) of the normal plasma (e− − p+) jet simulation. The orientation of the jet axis to our line-of-sight is θobs = 0°. Middle panel: fast-light Stokes I image of the normal plasma jet, but in which the image has been convolved with a circular Gaussian beam of FWHM 3.5 mas (shown in the lower left of the panel). Right panel: fast-light Stokes I image of the normal plasma jet, but in which an artificial Gaussian noise floor (set at ∼10−6 Jy beam−1) has been included to illustrate the extreme interferometric array sensitivity needed to detect our PIC jets even at such close proximity.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.