The Spiderweb Protocluster is Being Magnetized by Its Central Radio Jet

Supermassive black holes (SMBH) in active galactic nuclei (AGNs) interact with the larger cosmos via production of magnetized, radio-emitting synchrotron jets and lobes ¹³ (e.g., Gaspari et al. 2019; Marsden et al. 2020). This can help drive cosmic ecology and evolution (e.g., Hardcastle & Croston 2020) by controlling the flow of SMBH-bound gas on scales ranging from megaparsecs (i.e., in galaxy cluster gaseous halos; Gaspari et al. 2019 and references therein) down to milliparsecs (i.e., in SMBH accretion disks; e.g., Event Horizon Telescope Collaboration et al. 2021; Narayan et al. 2022), by heating galactic and intergalactic gas to limit star formation and the growth of massive galaxies (e.g., Croton et al. 2006; Weinberger et al. 2017a, 2017b), or by conversely generating localized star formation via compression of gas clouds around the host galaxy (e.g., Croft et al. 2006; Gaibler et al. 2012; Fragile et al. 2017; Mukherjee et al. 2018), and by enriching the universe in metals (Reuland et al. 2007), chemical compounds (e.g., Russell et al. 2017; O'Sullivan et al. 2021), magnetic fields (e.g., Furlanetto & Loeb 2001; Arámburo-García et al. 2021), and cosmic rays (e.g., Hardcastle et al. 2009; Abdo et al. 2010; Eichmann et al. 2018; Vazza et al. 2021a).

Such interactions may have had their greatest impact on galactic evolution in the approximate redshift range 2 < z < 4 (e.g., Nesvadba et al. 2011; Hatch et al. 2014; Nesvadba et al. 2017; Falkendal et al. 2019), the nearside of which sees galactic star formation rates beginning to decline after cosmic noon (Förster Schreiber & Wuyts 2020), and a nascent red sequence of quenched galaxies already established (Kriek et al. 2008; Brammer et al. 2009). However, the precise roles that radio jets play in driving high-z galactic evolution remain unclear (Hardcastle & Croston 2020), and directly probing the jet-gas interaction regions at these redshifts remains challenging. Moreover, better-studied low-z radio galaxies are probably poor analogs for their high-z counterparts, since they operate in very different cosmic environments. For example, the mean density of cold gas in high-z galaxies is several times higher than in the present-day universe (e.g., Walter et al. 2020). The energy density of the cosmic microwave background scales as (1 + z)⁴ (meaning that inverse-Compton cooling of radio lobes may scale similarly, e.g., Wu et al. 2017; Hodges-Kluck et al. 2021; Carilli et al. 2022). The star formation rate is also an order of magnitude greater than present-day levels (Madau & Dickinson 2014), galaxy clusters are still in the process of forming (e.g., Muldrew et al. 2015; Tozzi et al. 2022), and magnetic fields may not yet fully permeate the intergalactic gas (Donnert et al. 2018), which can affect its viscosity, pressure support, and thermal conductivity. Therefore, to better understand how SMBHs have shaped the cosmic ecology, it is desirable to observe radio galaxies in nascent galaxy clusters at high redshift, using techniques that illuminate the locus of interaction between the jet and ambient gas, and probe the physics occurring therein.

Broadband radio spectropolarimetry represents a singular set of such techniques. Exploiting analysis of Faraday rotation and depolarization, provides exquisite probes of magnetic field strength and structure around radio-emitting plasma, and by extension, the physical processes operating therein (e.g., Cooper & Price 1962; Burn 1966; Conway et al. 1974; Kronberg & Simard-Normandin 1976; Taylor & Perley 1993; Farnsworth et al. 2011; O'Sullivan et al. 2012, 2013; Anderson et al. 2015, 2018a, 2018b, 2021; Gaensler et al. 2015; Johnston-Hollitt et al. 2015; Anderson 2016; Pasetto et al. 2018). Consider that the linear polarization state of radio emission can be described by a complex vector P , related to the Stokes parameters Q and U, the polarization angle ψ, the fractional polarization p, and the total intensity I as

$$\begin{eqnarray}&&{\boldsymbol{P}}=Q+{iU}={{pIe}}^{2i\psi }.\end{eqnarray} \tag{ 1 }$$

In traveling from source to observer, linearly polarized radiation will be Faraday rotated by magnetized thermal plasma along the line of sight (LOS) to an observer by an amount equal to

$$\begin{eqnarray}&&{\rm{\Delta }}\psi ={RM}{\lambda }^{2},\end{eqnarray} \tag{ 2 }$$

where λ is the observing wavelength, and RM is the Faraday rotation measure, which is related to the thermal electron density n_e [cm⁻³] and magnetic field B [μG] along the LOS as

$$\begin{eqnarray}&&{RM}=0.812{\int }_{{z}_{s}}^{0}\displaystyle \frac{{n}_{e}(z){B}_{| | }(z)}{{\left(1+z\right)}^{2}}\displaystyle \frac{{dl}}{{dz}}\mathrm{dz}\,\,\,\mathrm{rad}\ {{\rm{m}}}^{-2}\,,\end{eqnarray} \tag{ 3 }$$

where z_s is the redshift of the radio source, and the comoving path increment per unit redshift, dl/dz, is in parsecs. The RM thus provides a means to detect magnetized thermal plasma, and to help deduce its properties.

In terms of high-z protocluster/radio galaxy systems to target, there are none more compelling nor better-studied than J1140-2629 ("The Spiderweb Galaxy") at z = 2.16 (Miley & De Breuck 2008; Miley et al. 2006). The radio source is characterized by two powerful jets oriented roughly east–west, extending about 60 kpc from the AGN (Carilli et al. 1997; Pentericci et al. 1997). The jets appear to be interacting with gas in the system over a range of scales in multifaceted ways: The protocluster contains a dense agglomeration of galaxies extending to radii >100 kpc, whose projected long axis appears to align with that of the radio jets, and whose chain and tadpole morphologies may indicate the radio galaxy has impacted their star formation history (e.g., Zirm et al. 2005; Miley et al. 2006). Observations of CO emission with the Jansky Very Large Array (VLA) and Australia Telescope Compact Array has revealed a remarkable molecular gas halo, which is also aligned with the radio jet, and which extends to a projected radius of 70 kpc (Emonts et al. 2016). Atacama Large Millimeter/submillimeter Array observations of H₂O, CI, and CO emission suggest that the passage of the radio jet can induce condensation of cold molecular gas clouds, intrinsically linking the action of the radio jet with star formation activity in the protocluster halo (Gullberg et al. 2016). The entire protocluster is enveloped in a spectacular Lyα halo out to a projected radius of >200 kpc, which is again broadly aligned with the radio jet (Kurk et al. 2002; Miley et al. 2006), and which is strongly displaced by a newly-detected eastern radio lobe (Carilli et al. 2022). An alignment between the radio jet and bright X-ray emission in the system, first described by Carilli et al. (2002), has been confirmed and delineated in great detail with a new ∼700 ks Chandra observation (Carilli et al. 2022). The new images show a remarkably close correlation between the radio and the extended X-ray emission from the jet. The data are consistent with the extended X-ray emission arising mostly from inverse-Compton up-scattering of the local photon field by the radio-emitting relativistic electrons, with a minor thermal contribution from hot gas.

High resolution, broadband, full-polarization radio observations did not previously exist, but represent a singular probe of magnetized gas physics in and around the radio jet, and may thereby shed new light on jet-gas interactions in this archetypal protocluster system. We therefore undertook such observations using the Jansky VLA, which possesses subarcsecond resolution and spans 340 MHz to 36 GHz, as well as a 700 ks Chandra exposure, as part of a definitive new study of the Spiderweb system (Carilli et al. 2022; Tozzi et al. 2022). In this paper, we present an analysis of polarization and Faraday rotation in the Spiderweb system, in order to better understand the properties of the radio jet, the surrounding protocluster medium, and interactions between the two, at an epoch where such interactions may peak in their cosmological importance. The paper is laid out as follows: Section 2 describes our observations and their calibration, and Section 3 describes our polarimetric imaging and analysis, as well as details of ancillary data sets. Section 4 describes the results derived from our radio observations, and Section 5 describes further results obtained via analysis of ancillary data sets. We discuss our findings and conclude in Section 6. We adopt the current Planck cosmology (Planck Collaboration et al. 2020) for our calculations, such that one arcsecond on the sky corresponds to 8.49 kpc at the z = 2.16 redshift of the source.

We observed The Spiderweb Galaxy using the Jansky VLA in multiple bands and array configurations in full polarization for project 19A-024: 2–4 GHz at 13 × 06 resolution (S band in A array), 8–12 GHz at 04 × 02 resolution (X band in A array), 29.2–36.8 GHz at 04 × 02 resolution (Ka band in B array). A total of 8 hr was observed in each band. Additionally, in 2020 December the Jansky VLA was temporarily left in an unscheduled "BnA"-like hybrid configuration due to the COVID-19 pandemic. ¹⁴ We sought and obtained 8 hr of director's discretionary time in this configuration in the Ka band, again covering 29.2–36.8 GHz. These supplemental data primarily help improve our Ka-band sensitivity, but also help increase the resolution we would otherwise achieve across the east–west-oriented jet at these southerly latitudes. The final full-band sensitivities of the S-, X-, and Ka-band observations are 5.5, 3, and 5 μJy beam⁻¹, respectively. The S-band data has been presented in detail by Carilli et al. (2022); their comparatively low spatial resolution and long wavelength make them unsuitable for our polarization analysis, ¹⁵ so they will not be discussed in detail in this work. Further details of our X- and Ka-band observations are listed in Table 1.

Table 1. Observing Parameters

Target	J1140-2629 ("The Spiderweb Galaxy")
Project codes	19A-024, 20B-448
Execution block IDs	36425705, 37201837, 37213317, 39133187
Dates of observation	2019-02-28, 2019-09-06, 2019-09-08, 2020-12-06, 2020-10-06
Field center (J2000); [l, b]	11h40m4840, − 26d29m900; [283872 , +33758]
Frequency bands	X, Ka
Frequency range (X, Ka)	8–12, 29.2–36.8 GHz
Total integration time (X, Ka)	8, 16 hr
Full-band sensitivity (X, Ka) a	3, 5 μJy beam−1
Array configuration(s) (X, [Ka])	A, [B, BnA-like (see footnote 14)]
Angular resolution (X, Ka; robust = 0.25)	032 × 020, 031 × 017
Recorded polarizations	RR, RL, LR, LL
λ2 range (X, Ka)	6.25 × 10−4–1.4 × 10−3, 6.65 × 10−5–1.06 × 10−4 m2
RM spread function width (X; Ka; combined) b , c	4400; 89,000; 2600 rad m−2
Largest recoverable ϕ-scale b , c (X; Ka; combined)	5000; 47,000; 47,000 rad m−2
Largest recoverable ∣ϕ∣ b , c , d (X; Ka; combined)	97,000; 1,400,000; 1,400,000 rad m−2

Notes.

^a Measured per Stokes parameter in MFS images generated with a Briggs' robust weighting value of +0.25. ^b Calculated from equations in Section 6 of Brentjens & de Bruyn (2005), in the observer frame, to two significant figures. Source-frame values can be derived by multiplying tabulated values by a factor of (1 + z)² ≈ 10. ^c At greater than 50% sensitivity. ^d Calculated at the center frequency of the band.

Download table as:  ASCIITypeset image

We used the standard VLA pipeline ¹⁶ to flag the data for radio frequency interference (RFI), apply online calibration measurements and flags, and calibrate the delay, bandpass, flux, and gains. Observations of 3C286 were used to calibrate the bandpass, absolute flux density scale, RL-polarization phase, and absolute polarization angle. J1146-2447 was observed to calibrate the time-dependent gains, and to solve for the on-axis polarization leakage. For the Ka-band observations, J1331+3030 was observed to calibrate the telescope pointing.

To calibrate the absolute polarization angle and on-axis polarization leakage, we first had to de-apply the pipeline-derived instrument-to-sky-frame Stokes Q and U rotation. We then:

• 1.  
  flagged the cross-hand polarization data manually
• 2.  
  used the CASA task setjy to specify a model for the polarization intensity and angle of 3C286, with Stokes I,Q,U,V set to 1.88, 0.076, 0.233, and 0 Jy beam⁻¹ (respectively) at a reference frequency of 32 GHz with α = −0.8, and noting that RM ∼ 0 rad m⁻² for this source (Perley & Butler 2013), as well as the small λ² range spanned by the X and Ka bands.
• 3.  
  used the above model of 3C286 to solve for the cross-hand delays.
• 4.  
  applied the above cross-hand delay calibration to our data, then solved for the instrumental D-terms and polarization (simultaneously) every 4 MHz using observations of J1146-2447 repeated over a range of parallactic angles (Sault et al. 1996).
• 5.  
  applied the above cross-hand delay and leakage calibrations to our data, and then used our model of 3C286 to solve for the instrumental RL-polarization phase every 2 MHz.
• 6.  
  applied all of the above corrections, which complete the calibration of (on-axis) polarization, to the target data.

To verify the accuracy of the D-term leakage calibration, we applied it to 3C286, whose frequency-dependent polarized fraction is tabulated in Perley & Butler (2013). Their measured average P/I values for the source in the X and Ka bands (respectively) are 12.1% and 13.3% (±0.2% systematic), after extrapolation to the present epoch using their Table 5. Our corresponding measurements are 12.4% and 13.2%, indicating that our band-averaged leakage calibration accuracy is generally better than ∼0.3% of Stokes I.

The off-axis polarization response was not calibrated for these observations. However, the angular extent of the source is small when compared to the FWHM of the primary beam (PB) at our upper frequency limit (∼10'' compared to 12, respectively), so the source remains squarely in the regime where on-axis polarization effects dominate, and are corrected for.

The pipeline- and polarization-calibrated data were threshold flagged in all polarizations, and self-calibrated in phase (two rounds) and amplitude (one round) using CASA. We then imaged the data spectropolarimetrically with WSClean (Offringa et al. 2014). For all Stokes parameters, we generated channelized images across the full frequency range having 3500 × 3500 pixels, a pixel scale of 0008, and a Briggs (1995) robust weighting value of 0.25. We imaged and cleaned Stokes I multifrequency synthesis (MFS) maps independently, and then the Stokes I, Q, and U data cubes using WSClean's join polarizations and squared channel joining modes, with 100 and 200 MHz channelization in the X band and Ka band, respectively, and using automatic clean thresholding and masking (at 1 σ and 3σ of the full-band noise floor, respectively) based on local noise estimation.

We then smoothed the spectral images to the spatial resolution of our lowest frequency channel—05 × 05—then re-gridded to a common pixel grid, and concatenated together to form Stokes I, Q, and U data cubes with dimensions R.A., decl., λ².

We note that the image fluxes were not corrected for the effect of the PB, again because the small extent of the source in comparison to the primary beamwidth means that the magnitude of the correction (∼2%) is negligible.

We calculated the Faraday dispersion spectrum (FDS) in the observer frame over the range −1 × 10⁴ to +1 × 10⁴ rad m⁻² in increments of 250 rad m⁻² using RM synthesis, ¹⁷ (Burn 1966; Brentjens & de Bruyn 2005) applied to the Stokes Q and U data cubes with equal weighting per image channel. The result is a complex-valued FDS datacube with dimensions R.A., decl., and ϕ.

We generated a map of the peak polarized intensity (henceforth "peak P") and the Faraday-depth-at-peak-polarized intensity (which for convenience, we henceforth simply refer to as the RM ¹⁸ ) across the field from the FDS cube using Miriad's (Sault et al. 1995) moment function, which uses a three-point quadratic fit to the dominant peak of the FDS on a per-pixel basis to derive the amplitude and position of the peak. ¹⁹ We verified beforehand that the FDS did not, for the most part, appear multi-peaked or Faraday-complex (see, e.g., Anderson et al. 2015, 2016) throughout the source. We then masked both the peak-P and peak-ϕ maps at a full-band polarized signal-to-noise threshold of 7σ, which is required for reliable RM measurements (e.g., Macquart et al. 2012).

Except where otherwise noted, observer-frame RMs have been "K corrected" to the emitting (source) frame through multiplication by a factor of (1 + z)².

3.1. Ancillary Data

4.1. Total Intensity

Figure 1 shows S-, X-, and Ka-band MFS total intensity radio maps of the system. These are described in detail by Carilli et al. (2022). For this work, the features of interest in the total radio intensity maps are (a) the radio core, which is the eastern-most radio component visible in the western jet, and which is so identified due to its spatial coincidence with the bright X-ray-emitting AGN (see Figure 2 of Carilli et al. 2002), (b) a deflection and subsequent bifurcation of the western jet, which is co-located with bright knots of both radio emission and Lyα emission (also described by Kurk et al. 2002), and see Section 4.2.3), and (c) the relatively complex morphology of the eastern hotspot complex, which is resolved in the X and Ka bands. We note that the latter is embedded in a newly-detected eastern lobe described by Carilli et al. (2022), which is faintly visible in the S-band map.

Zoom In Zoom Out Reset image size

4.2. Polarization

Figure 2 presents maps of peak P (row (a)), the source-frame RM (row (b)), intrinsic polarization angle (row (c)), and fractional polarization (row (d)) across the system, derived from a mixture of X- and Ka-band spectral data (specified in the caption). Figure 3 presents raw spectropolarimetric data for selected locations in the system, further described in Section 4.2.2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The appearance of the system is similar in both polarized and total intensity. An obvious difference is that the radio core is not detected in polarization. The upper limit on the core's Ka-band fractional polarization is ∼0.5%, but as we move westwards down the jet, the fractional polarization rapidly increases to ∼30% or more, after which the jet shows a complex mix of high and low fractional polarization values ranging up to 50%. The (X plus Ka band) RMs show similarly complex structure on a variety of scales throughout the eastern hotspot complex and western jet, down to scales we describe as interfaces (see, e.g., Anderson et al. 2018a), where the RM shows large changes in both magnitude and sign over scales smaller than the synthesized beamwidth. Conversely, the sky-projected magnetic field changes only slowly, and is broadly aligned with the jet axis over the entire length of the system. We investigate these basic features of the system in more detail below.

4.2.1. Magnetic Structure of a High-redshift Radio Galaxy

Our observations provide the most detailed ever maps of the magnetic field structure in a high-redshift source. The emission-weighted, sky-projected magnetic field orientation (Figure 2 panel (c)) is broadly aligned with the jet axis over the entire length of the system (∼70 kpc) where polarized emission is detected. At the beginning of the western jet (between locations 4 and 5 in Figure 2 panel (b)), the position angles of the jet and magnetic field are 261° versus 260° (respectively, on average, increasing eastwards from north, and with an uncertainty of ∼5° on the quoted measurements). Post-bifurcation, the jet consists of both a southern and western branch. In the southern branch (roughly between positions 5 and 6 in Figure 2 panel (b)) the measured orientation angles are 227° versus 210°, and likewise 281° versus 258° in the western jet extension. The eastern hotspot lies at a position angle of 95° from the nucleus. While the projected magnetic field orientation is relatively more structured in this location, its mean position angle is 120°. Thus, for both the eastern and western jets, the magnetic field orientation lies predominantly along the jet, but with a minor cross-jet component in some locations (of between 0° and 30°).

The southern jet branch appears brighter and broader in total intensity than the western jet branch, but terminates at a bright component or hotspot after less than an arcsecond. The latter is fainter and narrower, but appears to continue for two or more arcseconds, and terminates in a faint lobe that is visible in our maps at the S band, but only faintly so at the X band (see Figure 1). The western branch is characterized by positive RM and higher fractional polarization compared to the southern branch (Figure 2, panels (b) and (d), respectively), the latter of which suggesting that the magnetic field is more ordered here, and perhaps stronger. The emission-weighted magnetic orientation vectors appear to change smoothly into the brighter post-bifurcation southern branch from the pre-bifurcation main branch, but show an abrupt change moving from the former into the western branch. This is likely because the emission-weighted magnetic orientation vectors are dominated by the brighter southern branch, until the two branches separate completely on the plane of the sky.

Assuming that Faraday depolarization is negligible in the Ka band, the observed fractional polarization of the jet carries information about the intrinsic degree of ordering of the plane-of-sky-projected magnetic field, with p_obs related to the degree of field ordering as ${p}_{\mathrm{obs}}={p}_{\max }{(1+{B}_{{\rm{r}}}^{2}/{B}_{{\rm{o}}}^{2})}^{-1}$ (Burn 1966), where p_max ∼ 0.75 is the theoretical maximum fractional polarization generated by synchrotron radiation, B_r is the strength of the random component of the decomposed magnetic vector field, and B_o is likewise the strength of the ordered or regular component. The observed fractional polarization at the Ka band typically ranges between 25% and 40% in the western jet, and 10%–25% in the eastern hotspot complex, so that the ratio B_r/B_o ranges from 0.9–1.4 in the western jet, and from 1.4–2.5 in the eastern hotspot complex. Generally then, the magnetic field structure of the jet is quite well-ordered, with a substantial poloidal (along-jet) component. There are two classes of models that purport to explain the presence of this poloidal component—i.e., those that invoke (i) helical fields and a roughly cylindrical jet (Lyutikov et al. 2005), and (ii) velocity shear (e.g., Laing 1980; Hughes et al. 1989a, 1989b, 1991; Attridge et al. 1999; Kharb et al. 2005). Our analysis in Section 4.2.4 suggests that in the present case, scenario (i) applies in a sheath-like region around the jet, but is agnostic about scenario (ii).

4.2.2. Updated Search for Extreme Faraday RMs

A map of the intrinsic (source-frame) RMs calculated from the combined X- and Ka-band data is shown in Figure 2 (row (b)). Selected locations are marked with green crosses; corresponding measurements are presented in Table 2, and FDS and raw QU spectra in Figure 3.

Table 2. Source-Frame Faraday RMs for Selected Locations in the Spiderweb Radio Jet

Location Index (Figure 2, Panel (b))a	Coordinates	RM (Source Frame)	Projected Distance from AGN
	(J2000)	(rad m−2)	(kpc)
1	11h40m48945, −26d29m0886	+11,100 ± 400	66.2
2	11h40m48895, −26d29m0926	+5800 ± 180	60.5
3	11h40m48895, −26d29m0959	−4400 ± 360	60.7
4	11h40m48294, −26d29m0864	+270 ± 30	9.0
5	11h40m48135, −26d29m0937	−1200 ± 100	27.3
6	11h40m48077, −26d29m1032	−4400 ± 390	35.8
7	11h40m48002, −26d29m1004	+3500 ± 270	43.5

Note.

^a Note that the ∣RM∣ values shown in Figure 2, panel (b) saturate at 2000 rad m⁻² as per the figure caption.

Download table as:  ASCIITypeset image

We detect intrinsic ∣RM∣ values consistent with those reported by Carilli et al. (1998) for the eastern hotspot complex, but also detect similarly large ∣RM∣ values throughout several fainter regions of the jet. In the north-most part of the eastern hotspot complex, we measure an ∣RM∣ value of 11,100 ± 400 rad m⁻², which is roughly double the previously largest reported value. It is noteworthy that the most extreme ∣RM∣ values in the system tend not to occur in the central regions of the protocluster, but are found (a) toward the eastern and western terminal points of the radio jet (see Table 2), and (b) toward the edges of the jets in general. The latter is relevant to results in Section 4.2.4.

We also detect interfaces in RM, where the value changes dramatically over spatial scales smaller than the synthesized beamwidth—from +5800 ± 180 rad to −4400 ± 360 rad m⁻²in the most extreme case visible in the eastern hotspot complex. The RM values often cross through zero over these interfaces, which implies an associated change in the LOS-projected magnetic field direction. The interfaces may also be associated with some of the most heavily depolarized regions visible in the lobes at the Ka-band fractional polarization map (Figure 2, row (d)). Given the frequency, the depolarization is unlikely to be due to Faraday effects per se, but rather geometric depolarization arising from associated complex magnetic field structures in the emitting region. Similar structures have been described in Fornax A, where they were found to be associated with complex internal lobe structure and entrained thermal gas (Anderson et al. 2018a).

Because extreme RM values often coexist with extreme RM dispersion, the most extreme ∣RM∣ values are often best sought by restricting the analysis to the highest available frequency bands. Otherwise, the signal can be overlooked due to strong depolarization, or overwhelmed by strongly polarized emission emerging from lower-RM regions falling along the same sight lines. We therefore conducted a search for even higher ∣RM∣ values than reported above using the Ka-band data only. This analysis did not reveal RM values significantly different from those shown in Figure 2 (row (b)), though we note that the Ka-band-only source-frame RM uncertainties are typically thousands of rad per square meter for the sensitivity attained in these observations.

4.2.3. A Co-spatial Lyα Knot, Jet Bifurcation, and RM Enhancement

Early studies of the Spiderweb system (e.g., Pentericci et al. 1997; Kurk et al. 2002) highlighted that the western radio jet changes course and bifurcates at the location of a bright knot of Lyα emission. Figure 4 shows that the RM is also enhanced at this location, jumping from ∣RM∣ = 200 ± 40 to 1200 ± 50 rad m⁻² at the position of the Lyα knot (moving WSW down the jet), and then immediately dropping back to ∣RM∣ = 0 ± 80 rad m⁻² —corresponding to an ∼1100 rad m⁻² enhancement (relative to the average of the ∣RM∣ values measured immediately upstream and downstream of the knot) spanning ∼15 (12.75 kpc). While larger ∣RM∣ enhancements occur in other locations in the jet, this particular enhancement distinguishes itself by virtue of (a) being spatially coincident with the bright knots of both Lyα and polarized radio emission, as well as the jet bifurcation and deviation as previously noted, (b) being immediately bracketed by two regions with much lower absolute RM values, and (c) the direction of steepest ∣RM∣ gradient to and from this location being oriented along the jet axis, rather than perpendicular to the jet axis, as it generally is for other notable RM enhancements (see Section 4.2.4). We conclude that the jet impinges on the Lyα emission region, compresses the gas and magnetic fields in this location, and thereafter deviates from its course and bifurcates. In Section 5.1, we use ancillary data (Section 3.1) to constrain the nature of the Faraday-active medium associated with this Lyα knot, and calculate the conditions in the gas.

Zoom In Zoom Out Reset image size

4.2.4. Transverse Faraday RM Gradients

Figure 5 (top panel) shows the same hue-intensity map (tracing source-frame RM and peak P, respectively) as Figure 4, but presented with a diverging binary color map, and the ∣RM∣ enhancement/Lyα knot region masked, to emphasize the structure that we will now discuss. It is evident that regions of the jet northward of the main east–west jet axis show predominantly positive RMs, whereas the converse is true southward of the jet axis. Moreover, larger ∣RM∣ values in the system tend to be found toward the edges of the jet (Section 4.2.2). To illustrate these effects, we extracted the RMs along the cross sections labeled A–E in the top panel of Figure 5, and plot the results in the bottom panel of the same figure. In each case, the RM along the cross sections show large, monotonic decreases in RM, which typically cross through zero. While the cross sections only extend over between 1 and 3 total intensity synthesized beamwidths, we are confident that the RM gradients are genuine because:

• 1.  
  polarization is a vector quantity, meaning that structures smaller than the FWHM of the total intensity synthesized beam can and will be revealed in RM maps (in fact, the same is true of total intensity, given sufficient signal strength—e.g., Dabbech et al. 2018).
• 2.  
  Monte Carlo simulations in the VLBI literature show that RM gradients can be reliably detected over even a small fraction of a beamwidth—it is the statistical significance and monotonicity of the change in RM that are most important in this key respect (Murphy & Gabuzda 2012; Hovatta et al. 2012; Mahmud et al. 2013).
• 3.  
  the RM gradients are not confined to a single poorly-resolved knot of emission (as has been the case for some of the RM gradients claimed to exist in the VLBI literature), but rather exist along almost the entire length of the western jet.
• 4.  
  the extracted RM cross sections plotted in Figure 5 are recovered with similar maximum magnitude (∼2000–4000 rad m⁻² ), RM span (∼4000–8000 rad m⁻²), and orientation (increasing most rapidly from south to north across the jet axis), and in three cases span two or more synthesized beamwidths. None of this is expected if the gradient structure is in fact produced by random unresolved noise (e.g., Hovatta et al. 2012).

Finally, we note that we performed experiments where we intentionally offset the X-band and Ka-band polarization maps by several pixels in both R.A. and decl., to assess whether self-calibration-induced astrometric offsets could spuriously cause the observed RM gradients. The RM gradients remained unchanged.

Zoom In Zoom Out Reset image size

Thus we claim to have robustly detected these RM gradients. Their presence represents additional evidence to constrain the nature of the Faraday rotation observed throughout the system. The RM reversal across the gradients implies a change in the LOS-projected magnetic field direction. The effect is observed along at least ∼20 kpc of the western jet, discounting the local change in RM attributed to the Lyα knot in Section 4.2.3. It follows that the RMs are generated in the immediate vicinity of the jet—it seems implausible that the unassociated ICM in the foreground knows about the axis of the radio jet, and is configured to generate zero RM precisely along this locus of points, with predominantly positive and negative RM above and below this line, respectively. On the other hand, the synchrotron-emitting particles in the jet are probably not well mixed with the Faraday-rotating medium, since relevant depolarization models (like the Burn slab or internal Faraday dispersion; see, e.g., O'Sullivan et al. 2012) generally predict polarized fractions of less than 10% at the X band for measured ∣RM∣ values greater than ∼1000 rad m⁻² . This is lower than what we observe throughout much of the jet. We conclude that the Faraday rotation probably originates in a magnetized cocoon around the jet.

The VLBI literature argues that significant and monotonic RM gradients typically imply the presence of helical or toroidal magnetic fields that are directly coupled with the radio jets (e.g., Mahmud et al. 2013; Gabuzda 2018; Gabuzda et al. 2018). While this type of analysis has most often been performed on parsec scale jet emission, there is emerging evidence that RM gradients are also observed on kiloparsec scales in some sources (Christodoulou et al. 2016; Knuettel et al. 2017). The sense of kiloparsec-scale RM gradients may be an important probe of accretion disk physics, as we discuss in Section 6. We therefore note that in the present case, the observed RM gradient increases from south to north in the western jet, indicating that the toroidal component of the magnetic fields has a clockwise orientation when viewed along the jet from the central AGN (Equation (3); also see Figure 6).

Zoom In Zoom Out Reset image size

5.1. Magnetic Fields in the Lyα Knot

5.2. A Magnetized ICM?

Cluster ICMs are generally observed to be magnetized in the local universe (e.g., Heald et al. 2020; Anderson et al. 2021 and references therein), which is often invoked to explain the RM structure of radio jets embedded therein (e.g., Dreher et al. 1987; Murgia et al. 2004; Vogt & Enßlin 2005; Bonafede et al. 2010, but see Anderson et al. 2018a). But can we assume this equally applies to the Spiderweb system and its (proto-)ICM at z = 2.16 (Section 3.1.1) specifically, or indeed, high-redshift protocluster environments more generally? Certainly, the apparent depolarization of the radio core (Figure 2 panel (d); Section 4.2) is likely due to large RM dispersion in this region (e.g., Burn 1966; Gardner & Whiteoak 1966). But only ∼1'' (8.5 kpc) to the west of the radio core, we observe the fractional polarization increase to ∼30%, with comparatively low net ∣RM∣ values of ∼+100 rad m⁻² . This suggests that the radial extent of the gaseous structure responsible for the strong depolarization is ≲8 kpc, which is similar to that measured for gas bound to embedded cD-type galaxies, rather than the ICM proper (e.g., Paolillo et al. 2002).

Since the X-ray data do not yield a detailed profile of n_e , we must resort to heuristics to understand the possible contribution of the ICM to the observed Faraday rotation measures, and by extension, reasonable values for the ICM magnetic field strength. To do so, we assume a simple model where the LOS depth D to each part of the radio source is identical (reasonable, given the radial extent of the ICM is a factor of three larger than that of the radio source). We approximate the electron density as being constant throughout this volume. We further assume that any given LOS passes through a series of cells of size l in which the orientation of the random magnetic field component (having strength B_r ) is effectively independent of surrounding cells, and so likewise with the RM contribution of any sequence of cells falling along a given sight line. In such a case, the cluster gas will produce a dispersion in the RMs observed across independent sight lines of Gaensler et al. (2001):

$$\begin{eqnarray}&&\displaystyle \frac{{\sigma }_{\mathrm{RM}}}{\mathrm{rad}\ {{\rm{m}}}^{-2}}=\displaystyle \frac{812}{\sqrt{3}}\left(\displaystyle \frac{{n}_{e}}{{\mathrm{cm}}^{-3}}\right)\left(\displaystyle \frac{{B}_{r}}{\mu {\rm{G}}}\right)\sqrt{\left(\displaystyle \frac{D}{\mathrm{kpc}}\right)\left(\displaystyle \frac{l}{\mathrm{kpc}}\right)}.\end{eqnarray} \tag{ 4 }$$

Taking n_e = 1.5 × 10⁻² cm⁻³ and D = 100 kpc, σ_RM becomes a function of only B and l, which we present as a contour plot in Figure 7. We then use bootstrap resampling to estimate the actual RM dispersion in the radio jet, and its 90% confidence interval, obtaining ${\sigma }_{\mathrm{RM},\mathrm{obs}}={3120}_{-920}^{+890}$ rad m⁻². From Figure 7, supposing for the sake of argument that our estimate of B_hot ∼ 9 μG in the jet-Lyα gas interaction region (Section 5.1) extends to the entire ICM, and furthermore that D = 100 kpc is not a gross underestimate of the ICM scale, then σ_RM values of the required magnitude can only be produced if the magnetic reversal scale satisfies l ≳ 25 kpc, or 25% of the cluster radius. This value would seem too large for an actively forming cluster. If the magnetic field strength is weaker than the assumed value, the magnetic reversal scale constraint increases, though we caution that a substantially stronger magnetic field could reduce the constraint derived above to more plausible values. Nevertheless, we conclude that RM structure observed on scales smaller than ∼25 kpc in the jets is probably not associated with mere RM dispersion from the ambient ICM, but rather gas in the vicinity of the radio jet itself.

Zoom In Zoom Out Reset image size

So again we are left with the question: What is the most plausible range of values for B and l in the cluster gas? Consider the following. While we measure ∣RM∣ = 1100 rad m⁻² for the RM enhancement associated with the Lyα knot, the ∣RM∣ values upstream and downstream of the Lyα knot are either consistent with zero, or a factor of 4–12 times lower in absolute value. The same is true of the ∣RM∣ values found along most of the spine of the western jet before it bifurcates. In view of this and the other evidence we have supplied thus far suggest that the RMs are predominantly generated in the vicinity of the jets (transverse RM gradients, higher ∣RM∣ values near the edges and ends of the jet than near the cluster center and radio core), we propose a new prior ²¹ be adopted: Namely, that it is the relatively low ∣RM∣ values of ${ \mathcal O }(200)$ rad m⁻² that are generated by the ambient cluster gas, and over a path length some 100 kpc/10 kpc ≈ 10 times greater than the Lyα RM enhancement. Since the gas pressure derived from our adopted volume-averaged electron density for the ICM (from Tozzi et al. (2022); see Section 3.1.1) is roughly equal to the pressure required to confine the Lyα knot (Section 5.1), we contend that the lower cluster-based RMs cannot be easily attributed to lower n_e elsewhere in the cluster volume. The remaining possibilities are that (a) the magnetic field strength is at least a factor of 10 lower in the broader cluster gas than in the hot gas confining the Lyα knot, or that (b), numerous magnetic field reversals occur along any given sight line through the ICM, which will tend to reduce the observed RMs by a factor of $\sqrt{D/l}$. This latter effect is weak however, and would require a reversal scale that is smaller than the typical inner scale of turbulence in Galaxy clusters—thought to be ${ \mathcal O }(1)$ kpc (e.g., Zhuravleva et al. 2018)—to limit the RM produced by a field of only a few microgauss to only ${ \mathcal O }(200)$ rad m⁻² . For example, a pervasive 10 μG field requires 1000 reversals along the LOS (i.e., l = 0.1 kpc) to produce sufficiently small ∣RM∣ values, and a 3 μG field requires 100 reversals (l = 1 kpc). While the inner scale of turbulence in high-z protoclusters might be expected to differ from that of nearby and mature galaxy clusters, small reversal scales on the order cited above would also give rise to a much patchier RM map than is observed. Thus, we claim that a realistic upper limit on the mean magnetic field strength in the ambient ICM of the Spiderweb system is ${ \mathcal O }(1)$ μG.

5.3. Magnetic Fields in the Interaction Region

We have presented the most detailed map of the magnetic field structure in and around a high-z radio galaxy to date. The RM, fractional polarization, and projected magnetic field orientation all show smooth changes in structure on scales larger than the synthesized beamwidth. The projected magnetic field orientation is oriented almost exactly along the axis of the jet along its entire length, indicating a substantial poloidal component to the fields. Our detection of large-scale cross-jet RM gradients highlights a coherent toroidal field component (discussed in more detail below). The ratio of the ordered to random components in the plane-of-sky-projected magnetic field is of order unity. The magnetic fields in the radio jet are therefore well-ordered on scales as small as ∼a few kiloparsec, up to the size scale of the radio jet itself (∼a few tens of kiloparsec), only ∼3 Gyr after the big bang.

Contrary to previous work (e.g., Pentericci et al. 1997) in which large ∣RM∣ values observed toward the radio jet were assumed to arise in the ambient cluster gas, we have presented new evidence to suggest the Faraday rotation principally occurs in the immediate vicinity of the jet, and in the gas interacting with it, as follows. The most extreme ∣RM∣ values in the system do not occur toward the central regions of the protocluster, but toward the eastern and western terminal points of the radio jet, and the edges of the jet (Table 2; also Figure 2 (b)). The eastern hotspot complex shows RM interfaces (that is, large changes in RM over angular scales smaller than the synthesized beamwidth; Figure 2 (b)), which are correlated with structure in the total intensity emission (Figure 1). There are transverse gradients in RM across the western jet (Figure 5), which cross zero rad per square meter and thereby indicate a change in the LOS-projected magnetic field direction. This orientation flip is spatially coincident with the spine of the radio jet. We also observe a spatial correlation between a bright Lyα knot, the deviation and bifurcation of the radio jet in total intensity, and a local enhancement in the ∣RM∣ of ∼1100 rad m⁻² (Figure 4). Finally, we argued that the magnetic field strength in the cluster ICM is low, and cannot explain the RM structure observed on ∼arcsecond scales within the jet (Section 5.2). Taken together with some key findings of Carilli et al. (2022)—i.e., that the synchrotron plasma is inverse-Compton-up-scattering the cosmic microwave background, and that the radio source has a hybrid FRI/II morphology which appears to be coupled to its one-sided displacement of Lyα-emitting gas in the system—a consistent picture emerges of a radio jet undergoing vigorous interaction with the surrounding gaseous environment, which is likely still in the process of falling into the system and equilibrating. The strong jet-gas interaction signatures may mean that the radio jet drives outflows observed in the system, rather than radiation pressure and winds from the central AGN, though these mechanisms are all energetically capable (Seymour 2012, Nesvadba 2017). It is significant that this is observed at z ∼ 2, where the importance of quasar-mode AGN feedback wanes, the importance of jet-mode or radio-mode feedback begins to dominate, and galaxies transition from more active to passive modes of evolution (Section 1).

The RM gradients identified in Section 4.2.4 increase anti-clockwise when viewed down the jet from the central AGN, presumably tracing toroidal magnetic field components with a clockwise orientation. Combined with our results that the sky-projected magnetic field lies broadly parallel to the jet along its length (Section 4.2.1), and that this field is well-ordered (Section 4.2.1), it is likely that the overall magnetic field structure of the jet is helical. Christodoulou et al. (2016) and Gabuzda (2018), Gabuzda et al. (2018) have cited a statistical preference for clockwise toroidal magnetic field components observed at decaparsec–kiloparsec scales (coupled with counterclockwise toroidal magnetic field components observed at parsec scales; Contopoulos et al. 2009) as evidence of a predicted manifestation of the Poynting–Robertson cosmic battery (Christodoulou et al. 2016) generating nested helical magnetic fields (e.g., see Figure 6; see also Figure 1 of Contopoulos et al. 2009 and Figure 2 of Gabuzda 2018). The sense of the RM gradient we detect is consistent with this picture. The results imply that the net current is flowing outward in the cylindrical region enclosed by the magnetic fields that dominate the RM signal. The cosmic battery is a means of generating strong, ordered, large-scale magnetic fields in the accretion disks of AGN, through which radio jets can be seeded with helical magnetic fields, and which in turn may transport magnetic flux from AGN into the broader universe. We propose that the Spiderweb radio galaxy may show the cosmic battery in action.

The hierarchy of estimated magnetic field strengths in the system is of further relevance to the findings just discussed. The magnetic field strength is ∼90 μG the terminal regions of the jets, where interaction with the surrounding gas is obviously maximal. The field strength drops to 50–70 μG in the radio lobe volume (Carilli et al. 2022), then drops further to typically 10 μG in the gas interacting with the jets, and finally drops again to perhaps 1 μG in the ambient ICM. We note that the magnetic pressure in and around the jet does not appear to dominate the jet's dynamical interactions, but it may not be entirely negligible given the substantial uncertainties. Considering (a) the cascade of decreasing magnetic field strengths as we go from the lobes, to gas interacting with the lobes, to the ambient ICM; (b) the proposed operation of a cosmic battery in the jets; (c) that the radio galaxy is embedded in a protocluster environment at high redshift; and (d) the evidence for strong jet-gas interactions throughout the system, we further propose that the Spiderweb radio jet may be in the process of magnetizing the cluster gas, at a redshift which is relevant to broader questions around how the universe became pervasively magnetized, as it appears to be (e.g., see Neronov & Vovk 2010; Chen et al. 2015; Govoni et al. 2019; Vazza et al. 2021b). Recent results from the IllustrisTNG simulations (Arámburo-García et al. 2021) appear to back this up, not only demonstrating the probable cosmological importance of AGN-driven magnetization of the universe generally, but also the escalating importance of these processes at z ∼ 2 specifically—again placing the z = 2.16 Spiderweb system near the dawn of an important cosmic epoch as far as these questions go.

Finally, our clear detections of (i) ∣RM∣ values up to 10,000 rad m⁻², (ii) kiloparsec-scale transverse RM gradients, and (iii) radio-mode feedback, particularly in the form of an ∣RM∣ enhancement co-located with a bright Lyα knot and jet path deviation, are all novel for high-z radio galaxies, and demonstrate the value of modern spectropolarimetric analysis in this regime. In making these observations, however, even the Jansky VLA is near to its effective limits of sensitivity, resolution, and object-frame λ² coverage, such that it is probably impossible to study statistical samples of high-z radio galaxies in comparable detail at the present time. Such studies must wait for the combination of sensitivity, resolving power, and high-frequency bandwidth coverage that will be provided by the Next Generation VLA (Murphy et al. 2018), and the Square Kilometer Array (Heald et al. 2020, and references therein). In the interim, we plan to test the broader importance of our findings by targeting a select few other high-z radio galaxies using spectropolarimetric Jansky VLA observations, including 4C41.17 at z = 3.6 (see, e.g., Miley & De Breuck 2008).

We thank the anonymous referee for their time and helpful comments, which have significantly improved the content and presentation of this paper. This work employed VLA data from projects 19A-024 and 20B-488. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under a cooperative agreement by Associated Universities, Inc. C.S.A. is a Jansky Fellow of the NRAO. C.L.C. acknowledges support through CXC grant G09-2-1103X. S.B. acknowledges financial support from the INFN INDARK grant. P.T. and S.B. acknowledge financial support from the agreement ASI-INAF n.2017-14-H.0. A.S. is supported by the ERC-StG ClustersXCosmo' grant agreement 716762, by the FARE-MIUR grant "ClustersXEuclid" R165SBKTMA, and by INFN INDARK grant.