# X-ray Emission from the Jets and Lobes of the Spiderweb

Christopher L. Carilli11affiliation: National Radio Astronomy Observatory, P.
O. Box 0,Socorro, NM 87801, USA<EMAIL_ADDRESS>ORCID: 0000-0001-6647-3861
, Craig S. Anderson22affiliation: Jansky Fellow of the National Radio
Astronomy Observatory, P. O. Box 0,Socorro, NM 87801, USA , Paolo
Tozzi33affiliation: INAF - Osservatorio Astrofisico di Arcetri, Largo E. Fermi
5, 50127, Firenze, Italy , Maurilio Pannella44affiliation: Astronomy Unit,
Department of Physics, University of Trieste, via Tiepolo 11, I-34131 Trieste,
Italy , Tracy Clarke55affiliation: Naval Research Laboratory, Code 7213, 4555
Overlook Ave SW, Washington, DC 20375, USA , L. Pentericci1414affiliation:
INAF-OAR, Via Frascati 33, Monte Porzio Catone (RM), Italy , Ang
Liu33affiliation: INAF - Osservatorio Astrofisico di Arcetri, Largo E. Fermi
5, 50127, Firenze, Italy , Tony Mroczkowski66affiliation: Max Planck Institute
for Extraterrestrial Physics, Giessenbachstrasse 1, 85748 Garching, Germany ,
G.K. Miley77affiliation: Leiden Observatory, Leiden University, P.O.Box 9513,
NL-2300 RA, Leiden, The Netherlands , H.J. Rottgering77affiliation: Leiden
Observatory, Leiden University, P.O.Box 9513, NL-2300 RA, Leiden, The
Netherlands , S. Borgani88affiliation: Dipartimento di Fisica dell’ Università
di Trieste, Sezione di Astronomia, via Tiepolo 11, I-34131 Trieste, Italy 9
9affiliationmark:  10 10affiliationmark:  11 11affiliationmark:  , Colin
Norman1212affiliation:  Space Telescope Science Institute, 3700 San Martin
Dr., Baltimore, MD 21210, USA 1313affiliation: Johns Hopkins University, 3400
N. Charles Street, Baltimore, MD 212 , A. Saro88affiliation: Dipartimento di
Fisica dell’ Università di Trieste, Sezione di Astronomia, via Tiepolo 11,
I-34131 Trieste, Italy 99affiliation: INAF-Osservatorio Astronomico di
Trieste, Via G. B. Tiepolo 11, I-34131 Trieste, Italy 1010affiliation: INFN,
Sezione di Trieste, Via Valerio 2, I-34127 Trieste TS, Italy 1111affiliation:
IFPU, Institute for Fundamental Physics of the Universe, via Beirut 2, 34151
Trieste, Italy , M. Nonino99affiliation: INAF-Osservatorio Astronomico di
Trieste, Via G. B. Tiepolo 11, I-34131 Trieste, Italy , L. Di
Mascolo88affiliation: Dipartimento di Fisica dell’ Università di Trieste,
Sezione di Astronomia, via Tiepolo 11, I-34131 Trieste, Italy

###### Abstract

Deep Chandra and VLA imaging reveals a clear correlation between X-ray and
radio emission on scales $\sim 100$ kpc in the Spiderweb radio galaxy at
z=2.16. The X-ray emission associated with the extended radio source is likely
dominated by inverse Compton up-scattering of cosmic microwave background
photons by the radio emitting relativistic electrons. For regions dominated by
high surface brightness emission, such as hot spots and jet knots, the implied
magnetic fields are $\sim 50~{}\mu$G to $70~{}\mu$G. The non-thermal pressure
is these brighter regions is then $\sim 9\times 10^{-10}$ dyne cm-2, or three
times larger than the non-thermal pressure derived assuming minimum energy
conditions, and an order of magnitude larger than the thermal pressure in the
ambient cluster medium. Assuming ram pressure confinement implies an average
advance speed for the radio source of $\sim 2400$ km s-1, and a source age of
$\sim 3\times 10^{7}$ years. Considering the lower surface brightness, diffuse
radio emitting regions, we identify an evacuated cavity in the Ly$\alpha$
emission coincident with the tail of the eastern radio lobe. Making reasonable
assumptions for the radio spectrum, we find that the relativistic electrons
and fields in the lobe are plausibly in pressure equilibrium with the thermal
gas, and close to a minimum energy configuration. The radio morphology
suggests that the Spiderweb is a high-$z$ example of the rare class of hybrid
morphology radio sources (or HyMoRS), which we attribute to interaction with
the asymmetric gaseous environment indicated by the Ly$\alpha$ emission.

Galaxy Formation; Radio Galaxies; X-ray Clusters

66affiliationtext: European Southern Observatory, Karl-Schwarzschild-Str. 2,
D-85748 Garching b. Munchen, Germany

##  1 Introduction

It is well documented that powerful radio AGN act as beacons to massive galaxy
formation in protoclusters in the early Universe (Miley & de Breuck, 2008;
Carilli et al., 2001). The powerful radio emitting jets are driven by mass
accretion onto supermassive black holes, and are thought to play a number of
important roles in galaxy formation. Radio-mode feedback, meaning heating of
the protocluster medium around the host galaxy by the radio jets on tens to
hundreds of kpc scales, is often invoked to suppress late-time gas accretion,
thereby limiting high mass galaxy growth (Fabian, 2012; McNamara et al.,
2005). Conversely, jet-induced star-formation, induced by compression of dense
gas in the host galaxies by shocks driven by the advancing radio source, may
enhance star formation locally in some systems (Gaibler et al., 2012). The
radio jets may also remove angular momentum from the AGN accretion disk,
thereby enhancing accretion rates and early SMBH formation (Blandford & Payne,
1982). These same jets can also transport magnetic fields and possibly heavy
elements to large scales in the protocluster (Carilli & Taylor, 2002).

Among the best studied powerful radio galaxies in a dense environment at large
look back time is the Spiderweb radio galaxy, J1140-2629, at $z=2.16$
(Pentericci et al., 1997, 1998; Kurk et al., 2000; Miley et al., 2006; Miley &
de Breuck, 2008), with high resolution imaging and spectroscopy at optical
through radio wavelengths. The radio source is characterized by two jets
oriented roughly east-west, extending about $7^{\prime\prime}$ ($\sim 60$ kpc)
from the AGN. The protocluster shows a dense distribution of galaxies
extending to radii $>100$ kpc, as well as a Ly$\alpha$ halo on a similar
scale, with a major axis along the jet axis (Miley et al., 2006). Many of
these galaxies have been detected in CO and thermal emission from warm dust,
including the host galaxy of the radio source itself (Gullberg et al, 2016;
Dannerbauer et al., 2014; Jin et al., 2021). There also appears to be a
diffuse CO emitting component on $\sim 100$ kpc scale (Emonts et al., 2016).

A 28 ksec (on-source) observation with Chandra revealed a dominant X-ray AGN,
plus a faint extension of the X-ray emission along the radio jet (Carilli et
al., 2002). However, a clean separation of the faint extended emission from
the AGN was difficult with the limited exposure, and physical conclusions
speculative. To date, a detailed understanding of the physical conditions
(pressures, magnetic fields), in the radio source and the ambient medium in
the Spiderweb complex remains limited.

We have embarked on a deep study of the Spiderweb galaxy using Chandra and the
Jansky Very Large Array. We have obtained a 700 ksec Chandra exposure, as well
as broad band radio observations with arcsecond to sub-arcsecond resolution
from 340 MHz to 36 GHz, and low spatial resolution observations with the ALMA
Small Baseline Array (Saro et al. in prep). In this paper, we present the
extended X-ray emission associated with the radio jets and lobes on $\sim 100$
kpc scales. We derive magnetic fields and pressures in the source under the
standard assumptions of inverse Compton up-scattering of CMB photons by the
relativistic synchrotron emitting electrons in the outer regions of the source
$>20$ kpc from the core. A future analysis will consider contributions of
other photon fields over the full source (Sobolewska et al., 2022). In
parallel papers we explore the polarized radio emission and, in particular,
study the magnetic field orientation in the radio source, and the extreme
rotation measures indicative of a magnetized ambient thermal gas (Anderson et
al., 2021), the proto-cluster AGN population based on the deep X-ray exposure
(Tozzi et al., 2021), and the nuclear X-ray source and diffuse emission away
from the radio jets (Tozzi et al., 2022).

##  2 Observations and data reduction

The Spiderweb Galaxy was the target of a Chandra Large Program observation
with ACIS-S (AO20, PI P. Tozzi). Data reduction is performed starting from the
level-1 event files with CIAO 4.13, with the latest release of the Chandra
Calibration Database at the time of writing (CALDB 4.9.3). The details of data
reduction, including removal of spurious events and residual cosmic rays, and
high-background filtering, are described in Tozzi et al. (2021). The total
effective exposure time is 715 ks.

We obtained deep VLA observations in September, 2019 at X band (8 GHz to 12
GHz) at $0.4^{\prime\prime}\times 0.2^{\prime\prime}$ resolution (major axis
oriented north-south, in all cases), and S band (2 GHz to 4 GHz), at
$1.3^{\prime\prime}\times 0.6^{\prime\prime}$ resolution. We use ’3 GHz’ to
designate the full S band, and ’10 GHz’ to designate the full X band. A total
of 8 hours was observed in each band.

Standard pipeline calibration was employed at the
VLA111https://science.nrao.edu/facilities/vla/data-processing/pipeline, using
3C 286 as the absolute flux density scale calibrator. The calibrated data were
then threshold flagged, and self-calibrated using both AIPS and CASA. The
final imaging process employed the CLEAN algorithm, with Briggs weighting and
a robust factor of zero (Briggs, Schwab, Sramak x, 1999). Multiscale clean was
employed, cleaning down to 2$\sigma$ residuals. Bandwidth synthesis was also
employed, across four 0.5 GHz bands spanning the 3 GHz band, and eight 1 GHz
bands across the 8 GHz band. In each case, the images were restored with CLEAN
Gaussian beams corresponding to the resolution at mid-band. The final images
were summed to produce a single image for each band, while the sub-band images
were also analyzed for spectral index across a given band. The resulting rms
noise on the full band images was 3$\mu$Jy beam-1 at 10 GHz, and 5$\mu$Jy
beam-1 at 3 GHz.

As part of the VLITE system (Clarke et al., 2016), commensal observations at
340 MHz (P-band) with a resolution of $8.0^{\prime\prime}\times
4.3^{\prime\prime}$ (major axis north-south), were obtained in parallel with
our higher frequency program. The VLITE 340 MHz data were processed using a
dedicated flagging and calibration pipeline developed within the Obit (Cotton,
2008) and AIPS (van Moorsel et al., 1996) data reduction packages. The
pipeline uses automated flagging routines to remove radio frequency
interference followed by standard delay, gain, and bandpass calibrations. The
calibrated VLITE data were imaged using the wide-band imaging algorithm
MFimage in Obit to map a region of radius 1∘ from the pointing center with
small facets. Two rounds of phase-only self calibration were applied to the
data. While all four days of observations were recorded with VLITE, the period
of time in early September 2019 was known to be significantly impacted by
ionospheric fluctuations, thus we concentrate on the lowest noise day (rms
noise 640 $\mu$Jy beam-1, Sept. 13, 2019) for our analysis herein.

When comparing the Chandra X-ray and VLA radio images, we employ the standard
astrometry for the Chandra observatory, as well as the astrometry set by the
VLA calibrators, using J2000 coordinates. A check of the relative astrometry
is given by the position of the active nucleus of the Spiderweb, for which we
find agreement between the radio and X-ray to within $0.1^{\prime\prime}$. The
X-ray and radio alignments we consider herein are for structures on
significantly large scales, $\sim 1"$ or greater. In particular, the
quantitative spectra considered average over areas with few arcsec scale (see
Figure 8).

We also employ an archival image of the Ly$\alpha$ emission from the
Spiderweb. Narrow band imaging was carried out on April 1999 with ESO VLT Antu
(UT1) using FORS1 in imaging mode (ESO programme 63.O-0477(A), P.I. Miley
(Kurk et al., 2000). The narrow band filter has a central wavelength of 381.4
nm and a FWHM of 6.5 nm. The diffuse emission in this filter is dominated by
the rest frame Ly$\alpha$ emission from the Spiderweb, although there remains
some continuum emission in the filter as well. Details of these VLT
observations and results will be presented in Nonino et al. (2022). Herein, we
present a morphological comparison of the Ly$\alpha$ and radio emission.

##  3 Results

###  3.1 Total Intensity Images

Figure 1 shows the total X-ray emission from 0.5 to 7 keV from the Chandra
Observatory as both contour and color scale. The X-ray counts image has been
convolved with a $1^{\prime\prime}$ Gaussian. The range 7 keV to 10 keV is not
used for imaging since the signal is negligible and the background very high.
However, the X-ray fluxes we derive refer to the range 0.5 keV to 10 keV (or
separately 0.5 keV to 2 keV and 2 keV to 10 keV), since they are obtained from
the full-band spectral analysis. Figure 2 shows the same X-ray image as color
scale, with the VLA 3 GHz contour image at $1.3^{\prime\prime}\times
0.6^{\prime\prime}$ resolution (major axis north-south).

The net counts (background subtracted) for the full source in the 0.5-7 keV
X-ray band is 10400, which is dominated by emission from the AGN, which has
$\sim 9600$ counts. The net 0.5 keV to 7 keV X-ray counts for the western jet,
outside of $2^{\prime\prime}$ radius from the core, is 200, while that for the
eastern jet is 272. There is a point X-ray source $5^{\prime\prime}$ northwest
of the AGN which appears to be associated with an optical galaxy in the
protocluster at a large velocity offset from the main Spiderweb host galaxy.
Likewise, there appears to be faint, diffuse emission north and south of the
AGN, possibly indicating thermal emission from a hot ambient medium on tens of
kpc-scales. We consider the AGN and other X-ray point sources, and thermal
emission, in parallel papers (Tozzi et al., 2021, 2022). The focus of this
paper is the extended X-ray emission associated with the radio jets.

The rms noise on the 3 GHz image is $5.5~{}\mu$Jy beam-1. The radio emission
to the west shows a prominent jet, defined by a series of hot spots and more
diffuse emission, extending at least $7^{\prime\prime}$ ($\sim 60$ kpc), from
the AGN, and becoming fainter with distance from the AGN. The emission to the
east shows a bright radio hot spot at the source extremity, about
$7^{\prime\prime}$ from the AGN, and fainter emission from a radio lobe
between this hot spot and the AGN. The width of the eastern radio lobe, and
the more diffuse regions of the western jet, is about $3^{\prime\prime}$
($\sim 25$ kpc), down to the $5\sigma$ level.

The extended X-ray emission follows closely the radio jets to both the east
and west. There is an X-ray hot spot toward the end of the eastern lobe, close
to, but not exactly aligned with, the radio hot spot, with an offset of $\sim
0.7"$. Again, we estimate the relative X-ray and radio astrometry is accurate
to $\sim 0.1"$, based on the strong AGN emission at both wavebands. We
consider this X-ray – radio hot spot offset in Sobolewska et al. (2022).

Figure 3 shows the 3 GHz radio contours, and the ESO Ly$\alpha$ image (Nonino
et al., 2022). Extended Ly$\alpha$ emission is seen from the Spiderweb, which,
at lower surface brightness, extends over more than 100 kpc, enveloping the
protocluster (Miley & de Breuck, 2008). Herein, we focus on the brighter
Ly$\alpha$ emission regions toward the center of the protocluster. Note that
the Ly$\alpha$ emission does not peak at the location of the AGN, likely the
result of dust obscuration toward the galaxy center. In general, obscuration
will be an issue when interpreting UV measurements.

There are two features of note with the comparison of the Ly$\alpha$ with the
radio jets. First, on the western side, the radio jet deflects to the
southwest just at a position of the brightest Ly$\alpha$ knot. This spatial
coincidence suggests a physical association of the ambient thermal gas and the
material within the radio jet, such as a strong shock where the fluids
collide. This region will be discussed at length in a future paper, including
the polarized radio emission (Anderson et al., 2021).

The second Ly$\alpha$ – radio (anti-) association is in the tail of the
eastern lobe. There is a clear deficit of Ly$\alpha$ emission at the position
of the diffuse radio lobe, and the Ly$\alpha$ emission follows the inner edge
of the radio lobe. Again, this anti-coincidence suggests a physical
association, in which the radio lobe somehow excludes the thermal gas, such as
has been seen for eg. radio cavities in lower redshift X-ray clusters
(McNamara et al., 2005).

Figure 4 shows the 10 GHz radio emission as black contours, plus a saturated
color scale image of the 3 GHz radio emission, with the lower 3 GHz contours
in white. The 10 GHz image has been convolved with the same Gaussian restoring
beam as the 3 GHz image, and has an rms noise of $5~{}\mu$Jy beam-1. This
image shows the close association of the brighter radio knots all along the
western jet, and for the eastern radio hot spot. However, the eastern radio
lobe is undetected at 10 GHz down to very low surface brightness levels,
implying a very steep radio spectrum, $\alpha<-2.5$.

###  3.2 Spectra and Spectral Index

Figure 5 shows the integrated radio spectrum of the Spiderweb from 150 MHz to
36 GHz. The data are roughly consistent with single power-law index of
$\alpha=-1.3\pm 0.1$, across the full frequency range, where spectral index is
defined as $S_{\nu}\propto\nu^{\alpha}$. There is evidence for modest
steepening at high frequency, with the spectral index changing from $-1.15\pm
0.15$ between 150 MHz and 1450 MHz, to $-1.4\pm 0.1$ from 3 GHz to 32 GHz. The
integrated spectrum mixes emission components across the source, including the
core, flatter spectrum knots and steep spectrum lobes, and hence is of limited
use for extrapolation of spectral properties in particular regions.

Figure 6 shows the spectral index between images at 2.2 GHz and 3.8 GHz, along
with the 3 GHz total intensity image. The radio core, and first western jet
knot, have spectral indices $\leq-1$, while the brighter knots further along
the western jet, and the radio hot spot at the end of the eastern lobe, have
steeper spectral indices $\sim-1.3\pm 0.2$. The outer jet to the west shows
spectral steepening to $\alpha<-2$ with distance from the core, while the
diffuse tail of the eastern lobe is very steep spectrum, with $\alpha\leq-2.5$
in regions where it could be measured. For comparison, the spectral index
study in Pentericci et al. (1997) shows steeper spectra in most regions
between 4 GHz and 8 GHz, indicative of significant curvature (steepening),
above 3 GHz.

Figure 7 shows the spectral index between 3 GHz and 340 GHz at lower spatial
resolution. In this case, the resolution is adequate to only get two
independent measurements: one for the western Jet, and one for the eastern hot
spot plus radio lobe. Both are consistent with the integrated spectral index
for the full source of $-1.3\pm 0.1$.

##  4 Analysis

###  4.1 Regions and Flux Densities

Figure 8 shows the regions chosen for analysis. These regions were chosen to
avoid emission from the bright AGN (radii $>2.5^{\prime\prime}$ from the AGN),
while capturing the outer eastern and western radio jet emission. We discuss
the choice of regions in more depth below. Table 1 lists the 340 MHz and 3 GHz
radio flux densities for the regions in columns 2 and 3. Column 4 lists the
spectral indices. For the larger regions, the spectral index is between 340
MHz and 3 GHz. For the smaller regions (K3 and K4), the spectral index is
between 2 GHz and 4 GHz. The errors include both off-source rms noise, plus
gain uncertainties estimated during the flux bootstrap calibration process.
For the west jet at 340 MHz, we subtracted the contribution of the AGN itself,
and the innermost west jet knot, since these contribute to the low resolution
340 MHz flux density, but are outside the region used for the X-ray summation.
The contributions (about 10%), were estimated by extrapolating from 3 GHz
using the measure spectral indices $\sim-1.2$ of these features (see Figure
6).

For the X-ray analysis of the jet emission, we have subtracted the background,
and we have subtracted any residual contribution from the bright, compact AGN
emission due to the far-out wings of the Chandra point spread function. We
have also subtracted possible diffuse thermal emission (Tozzi et al., 2022).
These components contribute $\sim 25\%$ to the measured totals in the regions.
The jet X-ray fluxes (in erg cm-2 s-1) in the soft band (0.5 keV to 2.0 keV),
and hard band (2.0 keV to 10 keV), for the regions are listed in columns 5 and
6. Also listed are the power-law indices for the X-ray emission, derived from
fitting power-law spectra jointly across the soft and hard X-ray bands (column
7). The X-ray spectral index in frequency, $\alpha$, relates to the X-ray
photon index as: $\alpha$ = 1.0 - photon index. The X-ray flux errors are
based on photon statistics, while the errors for the spectral indices are
derived during the fitting process. For the magnetic field and pressure
analysis below, the dominant uncertainties arise from the assumptions in the
model calculations, and we quote the model uncertainties explicitly where
relevant.

The X-ray frequency relates to photon energy, $E_{keV}$, as:
$\nu_{x}=2.44\times 10^{17}\times E_{keV}$ Hz. We adopt an effective
frequency, $\nu_{x}$, equal to the geometric mean of the low and high energy
limits of a given Chandra band (appropriate for a power-law spectrum of index
$\sim-1$, consistent with the X-ray measurements). The corresponding values
for the soft and hard bands are: $2.44\times 10^{17}$ Hz, and $1.11\times
10^{18}$ Hz, respectively.

To derive X-ray flux densities required for the inverse Compton calculation
below (in Jy = $10^{-23}$ erg s-1 cm-2 Hz-1), we adopt a spectral index of
$\alpha\sim-1.0$, for which the flux density becomes:
$S_{\nu}=S_{band}/(\nu_{x}\times ln(\nu_{high}/\nu_{low}))$, where $S_{band}$
is the measured band-integrated flux (in erg cm-2 s-1), and $\nu_{high}$,
$\nu_{low}$ are the high and low frequencies of the band. Again, $\nu_{x}$ is
the geometric mean between $\nu_{high}$ and $\nu_{low}$. This equation is
accurate to $\sim 15\%$ for $\alpha=-1.0\pm 0.2$. The X-ray flux densities at
the effective frequencies for the soft and hard Chandra bands are given in
columns 8 and 9 in Table 1. We include the $15\%$ uncertainty in flux density
due to spectral index uncertainty in the table.

###  4.2 X-ray Mechanisms

The very close association between X-ray and radio emission along both the
eastern and western jets in the Spiderweb suggests physically related emission
mechanisms. One possible mechanism to explain this spatial coincidence is
X-ray synchrotron radiation from an extrapolation of the radio emitting
relativistic electron population to very high energy. However, even if we
assume no spectral steepening between radio and X-ray frequencies,
extrapolating the radio flux densities into the X-ray under-predicts the
observed X-ray flux densities by a factor $\sim 40$. Hence, the observed
X-rays are much too bright to be an extrapolation of the radio synchrotron
spectra to very high energy.

A second mechanism to explain the relationship could be thermal emission from
shocked gas surrounding the expanding radio source. This mechanism was favored
by Carilli et al. (2002) based on shorter Chandra observations, where the
extended X-ray emission was detected at low significance and difficult to
separate from the AGN itself. The deep X-ray image presented herein implies a
substantially lower X-ray flux for the emission along the jet (the difference
arising from a poor subtraction of the AGN in the Carilli et al. (2002)
analysis). The new observations show that the X-ray emission is immediately
coincident with the radio jets, and not along the edges of the radio source –
the latter being expected for thermal X-ray emission from a shocked cocoon
(Snios et al., 2018).

Based on the morphology, flux, and spectral index (see Section 4.3), we posit
the most physically plausible mechanism to explain the Radio-X-ray association
along the jet is inverse Compton up-scattering of the ambient photon field by
the synchrotron emitting relativistic electrons. This mechanism has long been
considered in the context of X-ray emission from cosmic radio sources (Harris
& Grindlay, 1979), and it has been pointed out that up-scattering of the
cosmic microwave background (CMB) photons could become a dominant X-ray
emission mechanism in high redshift sources, due to the strong increase in the
CMB energy density with redshift (Hodges-Kluck et al., 2021; Schwartz, 2002;
Tavecchio et al., 2000). Numerous studies have now shown resolved X-ray
emission from high redshift radio jets is best interpreted as inverse Compton
emission (Scharf et al., 2003; Cheung et al., 2012; Worrall et al., 2020;
Marhsall et al., 2018; Snios et al., 2021; Hodges-Kluck et al., 2021; Ighina
et al., 2019), and even in some sources without extended radio emission
(Simionescu et al., 2016; Schwartz et al., 2019).

###  4.3 Magnetic Fields derived from CMB inverse Compton Emission

X-ray inverse Compton emission and the radio synchrotron emission arise from
the same power-law distribution of relativistic electrons: the X-rays
correspond to up-scattering of the ambient photon field, and the radio
emission is, in a quantum mechanical sense, up-scattering of the virtual
photons that comprise the static magnetic field. Hence, if one knows the
energy density of the ambient photon field, then the X-ray emission
effectively ’counts’ the relativistic electrons, which, using the radio
luminosity, then dictates the magnetic field strength (Longair, 2011).

An important consistency check is agreement between the X-ray and radio
spectral indices, given the emission arises from the same power-law
distribution of electrons. For all the regions analyzed herein, the spectral
indices in the radio and X-ray agree, within the errors, with a value of
$-1.2\pm 0.2$

Considering the ambient photon field, for high $z$ radio galaxies it has been
shown that the ambient photons relevant for comparison of the observed radio
synchrotron and X-ray inverse Compton emission are in the submm through far IR
(Simionescu et al., 2016; Schwartz et al., 2019; Scharf et al., 2003). While
the Spiderweb is known to be luminous in the far IR (Seymour et al., 2012),
the dominant FIR emission is isolated to the inner $\leq 1"$ radius (Gullberg
et al, 2016). Gullberg et al (2016) estimate an FIR luminosity of the host
galaxy of $4.4\times 10^{12}$ L⊙. For a centrally concentrated emitting
region, the photon density behaves as $3L/(4\pi r^{2}c)$, where $L$ is the
luminosity and $r$ the radius (eg. Scharf et al. (2003)). The CMB energy
density behaves as: $U_{CMB}=4.2\times 10^{-13}\times(1+z)^{4}$ erg cm-3.
Comparing the CMB to the centrally-sourced FIR, one finds the CMB dominates
the energy density beyond $\sim 20$ kpc, or $2.4^{\prime\prime}$, from the
nucleus. Note that our high resolution images imply that the synchrotron
photon energy density is below the CMB, even in the higher surface brightness
regions.

For the preliminary analysis presented herein, we will focus on the outer
regions of the source (radii $>20$ kpc from the core), and adopt the CMB as
the ambient photon field. A more complete treatment of the full extent of the
source, including spectral contributions from various photon fields will be
given in Sobolewska et al. (2022).

The magnetic field calculation assuming up-scattering of the CMB was first
presented by Harris & Grindlay (1979), and is given by:

$B_{G}^{(1-\alpha)}=5.75\times 10^{-17}\times(5.05\times
10^{4})^{-\alpha}\times(1+z)^{(3-\alpha)}\times[S_{r}\times\nu_{r}^{-\alpha}/S_{x}\times\nu_{x}^{-\alpha}]$

where $S_{x}$ and $S_{r}$ are the X-ray and radio flux densities at observed
frequencies, $\nu_{x}$ and $\nu_{r}$ (all in matching units), and spectral
index is defined as: $S=\nu^{\alpha}$. Again, $\alpha$ is dictated by the
power-law relativistic electron energy distribution, and should be common to
both the X-ray and radio emission.

We examine four regions for which adoption of a spectral index of $-1.2\pm
0.2$ in the X-ray and radio is reasonable. The regions are delineated in
Figure 8. Two large regions capture the integrated X-ray emission and radio
emission at 340 MHz, for the east and west jets beyond $\sim
2.5^{\prime\prime}$ from the AGN, respectively. These larger regions average
over the source structure east and west of the nucleus in large boxes that
contain all the emission, hence, the results are weighted toward the brighter
radio hot spots and jet knots. We also consider spatially resolved
measurements at 3 GHz of the two brightest knots beyond $2.5^{\prime\prime}$
from the AGN, for which the spectral index between 2 GHz and 4 GHz is
consistent with $\alpha=-1.2\pm 0.2$.

Table 1 lists the X-ray and radio flux densities for the four regions. Using
the Harris & Grindlay (1979) equation then leads to the magnetic fields listed
in Table 2. For the larger regions, we derive inverse Compton magnetic field
strengths of $\sim 50~{}\mu$G. The implied magnetic energy density
($B^{2}/8\pi$ erg cm-3, for $B$ in Gauss), is $\sim 1.0\times 10^{-10}$ erg
cm-3. For the two bright western jet knots, we obtain field values $\sim
70~{}\mu$G, and correspondingly a factor two higher magnetic energy densities.
The field values derived from the Soft and Hard bands are consistent to within
20%. An uncertainty in the spectral index of $0.2$ implies an uncertainty in
the magnetic field of $\sim 20$%.

###  4.4 Minimum Energy Radio Synchrotron Calculations

Assuming a minimum energy configuration of relativistic particles and magnetic
fields to estimate energy densities and magnetic fields in celestial
synchrotron sources was first presented by Burbidge (1956). We adopt equations
1 and 2 in Miley (1980). Note that assuming energy equipartition between
fields and relativistic particles, as opposed to a minimum energy
configuration, leads to field values 8% higher, and a total energy density
within 1% of the minimum (Readhead, 2011; Carilli, 1989).

We assume the standard values of: electron/proton energy density = 1, filling
factor = 1, a low frequency cutoff of 10 MHz, and a high frequency cutoff of
32 GHz (Miley, 1980). The calculation is relatively insensitive to the
frequency cut-offs: raising the low frequency cutoff by a factor five
increases the fields derived from the minimum energy assumption by 30%, and
the high frequency cut-off has an even smaller effect.

We have adopted a low frequency spectral index of $\alpha=-1.2$. Changing
$\alpha$ by $\pm 0.2$ alters the derived fields by about 25%, with steeper
$\alpha$ leading to higher fields. In all cases, we assume a line-of-sight
path length of 25 kpc ($\sim 3^{\prime\prime}$), corresponding to the width of
the eastern lobe and western jet down to the $5\sigma$ surface brightness
level.

The field values and minimum energy densities are listed in Table 2. The
values derived are comparable to the early analysis in Pentericci et al.
(1997). The typical field values derived assuming minimum energy are about a
factor two higher than those derived from the inverse Compton calculation. The
implied minimum energy densities (summed fields and relativistic particles),
are then about eight times higher than the magnetic energy density derived
from the inverse Compton calculation. This factor eight arises from the fact
that the minimum energy fields are $\sim 2\times$ higher, meaning a factor
four higher magnetic energy density, and then there is another factor two when
including the relativistic particles.

###  4.5 Pressures

The minimum energy calculation for the synchrotron emission assumes that a
minimum energy configuration is established by the relativistic particles and
magnetic fields, implying equipartition of energy between fields and
particles. There is no a priori reason to assume such a minimum energy
configuration is established. Hence, in the following analysis, we adopt the
inverse Compton fields, and consider the implications for pressures within the
radio source, in comparison to the ambient medium. The relevant quantity is
the total non-thermal pressure, meaning, the summed pressure in fields and
relativistic particles. Hence, the synchrotron radio emission provides a
critical constraint, since it must be explained in the context of the IC
derived magnetic field strength, and hence dictates the total non-thermal
pressure.

The synchrotron minimum energy density implies a minimum pressure, which, for
an isotropic relativistic gas equals $1/3~{}\times$ energy density. Hence, the
implied minimum pressure based on the synchrotron emissivity is $\sim 3\times
10^{-10}$ dyne cm-2 (Table 2). However, the inverse Compton magnetic field
strengths are a factor of two lower than dictated by the minimum energy
assumption. For fields lower than the minimum configuration, the total non-
thermal energy density (and pressure), increases as $(B/B_{minE})^{3/2}$, in
order to generate the observed radio emissivity. In this case, the pressure is
dominated by the relativistic particles. For reference, for fields larger than
$B_{minE}$, the total non-thermal pressure increases as $(B/B_{minE})^{2}$, as
the total pressure becomes dominated by the field (Readhead, 2011; Carilli,
1989). A field strength a factor two lower than minimum energy then implies a
factor $\sim 3$ higher total non-thermal pressure $\sim 9\times 10^{-10}$ dyne
cm-2.

The first conclusions based on the X-ray and radio emission from the jets of
the Spiderweb are that the non-thermal pressures within the radio source are
larger than dictated by minimum energy by a factor three, and are dominated by
relativistic particles, not magnetic fields.

The Ly$\alpha$ emission provides a rough estimate of the thermal gas pressures
in the Spiderweb protocluster. Pentericci et al. (1997) derive a typical
density for the Ly$\alpha$ emitting knots in the Spiderweb of $\sim 40$ cm-3.
Assuming a temperature of $\sim 10^{4}$ K implies a thermal pressure of
$5.5\times 10^{-11}$ dyne cm-2.

An important point to keep in mind is that the volume filling factor for the
dense Ly$\alpha$ emitting gas is small ($\sim 10^{-5}$; (McCarthy, 1993;
Pentericci et al., 1998)). Hence, there must be a hot, intercloud medium
confining the line emitting clouds. Tozzi et al. (2022) have nominally
identified X-ray emission from this hot, diffuse medium in the Spiderweb, on
the same spatial scale as the Ly$\alpha$ and radio jet emission, by
subtracting the strong X-ray AGN, and by avoiding regions that may be
contaminated by the X-ray jets. They estimate an average gas density and
temperature in a sphere of radius 50 kpc of $n_{e}\sim 0.014\pm 0.0014$ cm-3
and $T\sim 2.4_{-0.8}^{+0.5}\times 10^{7}$ K. The implied thermal pressure is
then $4.6\pm 1.3\times 10^{-11}$ dyne cm-2. Encouragingly, this hot gas
pressure is comparable to the pressure derived from the Ly$\alpha$ emission,
and both imply an ambient thermal gas pressure an order of magnitude lower
than the non-thermal pressures derived above for the bright radio emitting
regions. A more detailed analysis of the cluster density radial profile is
currently in progress.

The second conclusion therefore is that the high surface brightness radio
emitting regions (hot spots and jet knots), are over-pressured relative to the
ambient thermal gas in the cluster, by more than an order of magnitude.

Assuming ram pressure confinement of the radio structures then implies a
typical advance speed for the radio source $v_{ram}=(P_{nt}/\rho)^{1/2}\sim
2400$ km s-1, where $\rho$ is the ambient mass density. Using the source
maximum radius $\sim 60$ kpc implies a source age of $\sim 3\times 10^{7}$
years.

###  4.6 The Tail of the Eastern Lobe

We re-emphasize that the calculations for the regions above are weighted to
the higher brightness regions of the source, such as hot spots and jet knots,
and hence to regions of higher pressure. As a final calculation, we consider
the synchrotron minimum energy density in the very low surface brightness,
steep spectrum tail of the eastern radio lobe (but still beyond 20 kpc from
the core). This is the region that shows the ’cavity’ in the Ly$\alpha$
emission coincident with the radio lobe.

The difficulty in calculating the pressure in the diffuse radio lobe is,
again, that spatially resolving observations are only available at 3 GHz, and
that the spectrum at this high frequency is seen to be very steep, $\sim-2.5$,
but must flatten to lower frequency to explain the integrated radio spectrum
of the source. In this case, we use the typical surface brightness in the tail
of the radio lobe of $\sim 50~{}\mu$Jy beam-1 at 3 GHz, where the beam is
$1.3^{\prime\prime}\times 0.6^{\prime\prime}$. We extrapolate the 3 GHz
surface brightness to 1 GHz using the steep spectral index, implying a surface
brightness at 1 GHz of $\sim 780~{}\mu$Jy beam-1. We then calculate the
synchrotron minimum energy density using this 1 GHz surface brightness, and
the lower frequency spectral index, $\alpha=-1.2$, as above. This calculation
leads to a magnetic field strength of $37~{}\mu$G, a minimum non-thermal
energy density of $1.3\times 10^{-10}$ erg cm-3, and a minimum non-thermal
pressure of $4\times 10^{-11}$ dyne cm-2. Decreasing the turn-over frequency
from 1 GHz to 0.5 GHz increases the minimum non-thermal pressure by 60%.

Hence, if the diffuse tails of the radio lobes have a spectral flattening
between 0.5 GHz and 1 GHz, and if the fields and particles have been able to
reach a minimum energy configuration in these potentially more relaxed and
older regions, then the non-thermal radio lobe pressure is comparable to the
thermal pressures of the ambient medium.

##  5 Discussion

We comment on one curiosity and a general question. Considering the radio
source morphology, the eastern jet of the Spiderweb shows a radio hot spot at
the extremity of the lobe, with clear spectral steepening from the hot spot
back toward the AGN. Conversely, the western jet shows a gradual decline in
surface brightness, and spectral steepening with distance from the AGN. The
eastern morphology, in both intensity and spectrum, is characteristic of high
power radio Fanaroff-Riley Class II radio galaxies (’edge-brightened’,
typically with $P_{178MHz}>10^{35}$ erg s-1 Hz-1), while the western jet is
typical of lower power FRI (’edge-darkened’), radio jets (Fanaroff & Riley,
1974). However, the Spiderweb is one of the most luminous radio galaxies
known, with $P_{178MHz}\sim 4.6\times 10^{36}$ erg s-1 Hz-1, placing it
several orders of magnitude above the canonical FRII - FRI break. Why then the
dichotomy in both morphology and spectral index between lobes?

The Spiderweb is a high-redshift example of a hybrid morphology radio source
(or HyMoRS), which are a known but rare and poorly understood class of radio
galaxy, first described by Gopal-Krishna & Wiita (2000), and for which
correlated differences in morphology and spectral index between lobes have
previously been observed (Gawronski et al., 2006; de Gasperin, 2017). Their
nature is not well constrained: They may arise because of the physical
attributes of the source (e.g. jet power) or its environment (Kapinsk et al.,
2017), or because of observational or orientation effects (Harwood, Vernstrom
& Stroe, 2020), or by some combination of these.

In the case of the Spiderweb system, a potential explanation is seen in the
Ly$\alpha$ emission and galaxy density in the center of the cluster. The
thermal gas and galaxies are asymmetrically distributed in the cluster, with a
higher density toward the west than the east. It is possible that this density
inhomogeneity, in a dynamically active and young protocluster, leads to more
disturbed jet advance toward the west, while the eastern jet remains
relatively unperturbed until it reaches the terminal jet shock (hot spot). The
perturbed nature of the western jet is supported by the bright Ly$\alpha$ knot
at the jet deflection point, discussed above.

We note the very steep spectra for the diffuse emitting regions of the
Spiderweb, such as in the eastern lobe, and the end of the western jet, where
$\alpha<-2$ (see Figure 4 and Figure 6). This may, in part, be due to
increasing inverse Compton energy losses off the CMB for the relativistic
electrons. Indeed, this phenomenon has been invoked to explain the correlation
between redshift and ultra-steep spectrum radio galaxies (for a review, see
Miley & de Breuck (2008)). The VLA snapshot survey of Pentericci et al. (2000)
shows evidence that the majority of $z>2$ powerful radio galaxies have lobe
regions with very steep spectra ($\alpha<-2$ between 5 GHz and 8 GHz),
although the snap-shot sensitivity of the pre-upgraded VLA was inadequate to
properly image these regions in detail.

A general question is: given the rapid rise with redshift of the energy
density in the CMB, at what redshift will cooling of the relativistic
electrons be dominated by inverse Compton emission over synchrotron radiation?
The calculation is straight-forward, since the ratio of radio synchrotron to
X-ray luminosity is given simply by: $L_{r}/L_{x}=U_{B}/U_{CMB}$ (Schwartz,
2002; Carilli & Taylor, 2002). The magnetic energy density in a radio source
is: $U_{B}\sim 1\times 10^{-10}(B/50\mu G)^{2}$ erg cm-3, while the CMB energy
density is: $U_{CMB}=4.2\times 10^{-13}\times(1+z)^{4}$ erg cm-3. Hence,
$L_{x}>L_{r}$ at: $(1+z)>3.9\times(B/50\mu G)^{1/2}$. For typical radio jet
magnetic fields of 50 $\mu$G, inverse Compton dominates over radio synchrotron
losses at $z>2.9$. A recent analysis of 11 radio galaxies from $z=1.3$ to 4.3
shows that, for the diffuse emission, the X-ray and radio data are consistent
with both Inverse Compton scattering of the CMB, and with equipartition
magnetic fields (Hodges-Kluck et al., 2021).

##  6 Summary

Using the observed inverse Compton and radio synchrotron radiation, we derive
typical magnetic field strengths in the jets of the Spiderweb $\sim 50~{}\mu$G
to $70\mu$G, in regions of the source beyond 20 kpc from the core. These are
regionally averaged values, derived over areas of significant radio structure,
hence represent a brightness weighted mean, ie. weighted toward knots and hot
spots, corresponding to regions of particle acceleration due to shocks or
regions of high compression.

The inverse Compton fields are a factor two lower than those derived assuming
a minimum energy ($\sim$ equipartition), configuration for the synchrotron
emitting relativistic particles and fields. The implied non-thermal pressure
(sum of relativistic particles and fields), implied by the combination of
inverse Compton X-rays and radio synchrotron radiation, is $\sim 9\times
10^{-10}$ dyne cm-2. This non-thermal pressure is a factor three higher than
the pressure derived by assuming minimum energy conditions, and the pressure
is dominated by relativistic particles. Comparing the typical non-thermal
pressure in the radio source with the thermal pressure in the ambient cluster
gas derived from Ly$\alpha$ emission and from possible diffuse thermal X-ray
emission, implies a factor 15 higher pressures in the brighter radio emitting
regions than in the ambient thermal gas. Assuming ram pressure confinement of
the radio source by an hypothesized hot, diffuse cluster gas, implies an
advance speed $\sim 2400$ km s-1, and a radio source lifetime of $\sim 3\times
10^{7}$ years.

We see two regions where there is a clear spatial correlation between radio
structures and the Ly$\alpha$ emitting gas. One region is the brightest
Ly$\alpha$ knot in the western jet coincident with the position where the
radio jet deflects toward the southwest. This deflection suggests a strong
shock due to the impact of the radio jet on an ambient density enhancement. We
will discuss this region in detail in our polarization analysis of the radio
emission (Anderson et al., 2021).

The second spatial association of the Ly$\alpha$ and radio emission is seen in
the tail of the eastern lobe, where the radio lobe appears evacuated of
Ly$\alpha$ emission, while the line emission outer edge (away from the AGN),
aligns closely with the inner edge of the radio lobe. This anti-correlated
morphology is reminiscent of cavities in X-ray clusters due to expanding radio
sources (McNamara et al., 2005). A calculation of the minimum pressure in the
non-thermal gas in the faintest radio emitting regions, based on 3 GHz surface
brightness, and assuming a spectral flattening between 0.5 GHz and 1.0 GHz,
leads to approximate pressure equilibrium between the minimum non-thermal
pressures and the ambient thermal cluster gas. In the context of the standard
model for powerful radio galaxies, in which the lobes represent the older,
more relaxed, populations of relativistic particles and fields left behind by
the expanding radio source, it seems physically plausible that a minimum
energy configuration for fields and particles would be more likely in the
relaxed lobes. Spatially resolving observations at low frequency would help to
test this hypothesis.

We show that the Spiderweb radio source is a high redshift example of a
’Hybrid Morphology’ radio source, with the western jet have the morphology in
total intensity and spectral index of an FRI-class radio galaxy (edge-
darkened), and the eastern region having the morphology of an FRII-class radio
galaxy (edge-brightened) (Gopal-Krishna & Wiita, 2000). The origin of such
HyMoRS remains an area of active research. In the case of the Spiderweb, a
clue to the origin may be seen in the asymmetric gas environment, for which
the Ly$\alpha$ emission is predominantly west of the nucleus, potentially
leading to higher disruption of the advancing radio jet.

This program employed VLA data from project 19A-024. The National Radio
Astronomy Observatory is a facility of the National Science Foundation
operated under cooperative agreement by Associated Universities, Inc.. CLC
acknowledges support through CXC grant G09-2-1103X. SB acknowledges financial
support from the INFN INDARK grant and the agreement ASI-INAF n.2017-14-H.0.
MP, AS and LDM are supported by the ERCStG ‘ClustersXCosmo’ grant agreement
716762. AS is also supported by the FARE-MIUR grant ’ClustersXEuclid’
R165SBKTMA, and by INFN InDark Grant.

## References

  * Anderson et al. (2021) Anderson, C. et al. 2021, in prep.
  * Blandford & Rees (1974) Blandford, R. D. & Rees, M. J. 1974, MNRAS, 169, 395
  * Blandford & Payne (1982) Blandford, R.D. & Payne, D.G. 1982, MNRAS, 199, 883
  * Briggs, Schwab, Sramak x (1999) Briggs, D., Schwab, F., & Sramek, R. 1999, in Synthesis Imaging in Radio Astronomy II, eds. Taylor, Carilli, Perley, ASP, San Francisco 180, 127
  * Bertin (2006) Bertin, E. 2006, ASPC, 351, 112
  * Bertin (2002) Bertin, Emmanuel, Mellier, Yannick, Radovich, Mario, Missonnier, Gilles, Didelon, Pierre, Morin, Bertrand 2002, ASPC, 281, 228
  * Brown et al. (2021) Brown, A.G., Vallenari, A., Prusti, T., de Bruijne, J. et al. 2021, A& A, 649, A1
  * Burbidge (1956) Burbidge G. 1956, ApJ, 124, 416
  * Carilli et al. (2002) Carilli, C.L., Harris, D.E., Pentericci, L. et al. 2002, ApJ, 567, 781
  * Carilli (1989) Carilli, C.L. 1989, PhD Thesis, MIT
  * Carilli et al. (2001) Carilli, C. L., Miley, G., Röttgering, H. J. A. et al. 2001, in Gas & Galaxy Evolution, eds. Hibbard, Rupen, van Gorkum, ASP, San Francisco 240, 101
  * Carilli & Taylor (2002) Carilli, C.L. & Taylor, G.B. 2002, ARA&A, 40, 319
  * Carilli & Barthel (1996) Carilli, C.L. & Barthel, P. 1996, A& A Rev, 7, 1
  * Cheung et al. (2012) Cheung, C.C., Stawarz, L., Siemiginowska, A. et a. 2012, ApJ, 756, L20
  * Clarke et al. (2016) Clarke, T. E., Kassim, N. E., Brisken, W., et al. 2016, Proc. SPIE, 9906, 99065B. doi:10.1117/12.2233036
  * Cotton (2008) Cotton, W. D. 2008, PASP, 120, 439. doi:10.1086/586754
  * Dannerbauer et al. (2014) Dannerbauer, H., Kurk, J. D., De Breuck, C. et al. 2014, A& A, 570, 55
  * Emonts et al. (2016) Emonts, B., Lehnert, M., Villar-Martin, M. et al. 2016, Science, 354, 1128
  * Fabian (2012) Fabian, A. C. 2012, ARA&A, 50, 455
  * Fanaroff & Riley (1974) Fanaroff, B. & Riley, J. 1974, MNRAS, 167, 31
  * Gaibler et al. (2012) Gaibler, V., Khochfar, S., Krause, M., Silk, J. 2012, MNRAS, 425, 438
  * de Gasperin (2017) de Gasperin, F. 2017, MNRAS, 467, 2234
  * Gawronski et al. (2006) Gawronski, M., Marecki, A., Kunert-Bajraszewska, M., Kus, A. 2006, A& A, 447, 63
  * Gopal-Krishna & Wiita (2000) Gopal-Krishna, Wiita, P. 2000, A& A, 363, 507
  * Greve, Ivison, & Stevens (2006) Greve, T., Ivison, R., Stevens, J. 2006, AN, 327, 208
  * Gullberg et al (2016) Gullberg, B., Lehnert, M., de Breuck, C. et al. 2016, A& A, 859, 12
  * Harris & Grindlay (1979) Harris, D.E. & Grindlay, J.E. 1979, MNRAS, 188, 25
  * Harwood, Vernstrom & Stroe (2020) Harwood, J.J., Vernstrom, T., Stroe, A. 2020, MNRAS, 491, 803
  * Hodges-Kluck et al. (2021) Hodges-Kluck, E., Gallo, E., Ghisellini, G. et al. 2021, MNRAS, 505, 1543
  * Ighina et al. (2019) Ighina, L., Caccianiga, A., Moretti, et al. 2019, MNRAS, 489, 2732
  * Jin et al. (2021) Jin, S., Dannerbauer, H., Emonts, B. et al. 2021, A& A, 652, 11
  * Kapinsk et al. (2017) Kapinska, A., Terentex, I., Wong, O.I. et al. 2017, AJ, 154, 253
  * Kurk et al. (2000) Kurk, J., Röttgering, H. J. A., Pentericci, L. et al. 2000, ApJ, 351, L1
  * Longair (2011) Longair, Malcolm 2011, High Energy Astrophysics, Cambridge University Press, Cambridge, UK
  * Marhsall et al. (2018) Marshall, H. L., Gelbord, J. M., Worrall, D. M., et al. 2018, ApJ, 856, 66
  * McCarthy (1993) McCarthy, P. 1993, ARA&A, 31, 639
  * McNamara et al. (2005) McNamara, B. R., Nulsen, P. E. J., Wise, M. W., et al. 2005, Sci, 433, 45
  * Miley & de Breuck (2008) Miley, G. & de Breuck, C. 2008, A& ARv, 15, 67
  * Miley (1980) Miley, G. 1980, ARA&A, 18, 165
  * Miley et al. (2006) Miley, G., Overzier, R., Zim, A. et al. 2006, ApJ, 650, L29
  * Nesvadba et al. (2006) Nesvadba, N., Lehnert, M., Eisenhauer, F et al. 2006, ApJ, 650, 693
  * Nonino et al. (2022) Nonino et al. (2022), in prep
  * Pentericci et al. (2000) Pentericci, L., van Reeven, W., Carilli, C.L. et al. 2000, A& A Supp, 361, 25
  * Pentericci et al. (1998) Pentericci, L., Rottgering, H.J.A., Miley, G.K., et al. 1999, ApJ, 504, 139
  * Pentericci et al. (1997) Pentericci, L.; Roettgering, H. J. A.; Miley, G. K.; Carilli, C. L.; McCarthy, P. 1997, A& A, 326, 580
  * Readhead (2011) Readhead, A. 2021, Lecture Notes in Radiative Processes, Caltech, AY121
  * Scharf et al. (2003) Scharf, C., Smail, I., Ivison, R. et al. 2003, ApJ, 596, 105
  * Scheuer (1974) Scheuer, P. A. G. 1974, MNRAS, 166, 513
  * Schwartz (2002) Schwartz, D. A. 2002, ApJ, 571, L71
  * Schwartz et al. (2019) Schwartz, D. A., Siemiginowska, A., Snios, B. et al. 2019, ApJ, 904, 57
  * Seymour et al. (2012) Seymour, N., Altieri, B., De Breuck, C. et al. 2012, ApJ, 755, 146
  * Simionescu et al. (2016) Simionescu, A., Stawarz, L., Ichinohe, Y., et al. 2016, ApJ, 816, L15
  * Snios et al. (2021) Snios, B., Schwartz, Daniel A., Siemiginowska, Aneta, et al. 2021, ApJ, 914, 130
  * Snios et al. (2018) Snios, B., Nulsen, P.E., Wise, M.W., et al. 2018, ApJ, 855, 71
  * Sobolewska et al. (2022) Sobolewska, M. et al. 2022, in prep
  * Tavecchio et al. (2000) Tavecchio, F., Maraschi, L., Sambruna, R. M., & Urry, C. M. 2000, ApJ, 544, L23
  * Tozzi et al. (2021) Tozzi, P. et al. 2021, å, submitted
  * Tozzi et al. (2022) Tozzi, P. et al. 2021, in prep.
  * van Moorsel et al. (1996) van Moorsel, G., Kemball, A., & Greisen, E. 1996, Astronomical Data Analysis Software and Systems V, 101, 37
  * Worrall et al. (2020) Worrall, D. M., Birkinshaw, M., Marshall, H. L., et al. 2020, MNRAS, 497, 988

Table 1: Observed Properties Region |  S340MHz |  S3GHz |  $\alpha_{0.34}^{3}$ or $\alpha_{2}^{4}$ |  Ssoft |  Shard | $\alpha_{soft}^{hard}$ |  S2.4e17Hz |  S1.1e18Hz | Area
---|---|---|---|---|---|---|---|---|---
| Jy | mJy |  |  $10^{-15}$ erg cm-2 s-1 |  $10^{-15}$ erg cm-2 s-1 |  |  $10^{-11}$ Jy |  $10^{-11}$ Jy |  arcsec2
West Jet | $2.6\pm 0.4$ | $170\pm 12$ | $-1.2\pm 0.1$ | $1.2\pm 0.13$ | $1.4\pm 0.14$ | $-1.0\pm 0.2$ | $35\pm 6$ | $7.9\pm 1.4$ | 27
East Jet | $1.8\pm 0.3$ | $111\pm 8$ | $-1.3\pm 0.1$ | $1.5\pm 0.14$ | $1.1\pm 0.12$ | $-1.3\pm 0.2$ | $46\pm 8$ | $6.4\pm 1.1$ | 30
K3 | – | $51\pm 4$ | $-1.2\pm 0.2$ | $0.22\pm 0.05$ | $0.16\pm 0.04$ | $-1.3\pm 0.2$ | $6.5\pm 1.8$ | $0.92\pm 0.24$ | 2.4
K4 | – | $86\pm 6$ | $-1.4\pm 0.2$ | $0.33\pm 0.6$ | $0.41\pm 0.8$ | $-1.0\pm 0.2$ | $9.8\pm 2.4$ | $2.3\pm 0.6$ | 5.1
Table 2: Magnetic Fields and Energy Densities Region |  BIC,soft |  UBIC,soft |  BIC,hard |  UBIC,hard |  BminE |  UminE
---|---|---|---|---|---|---
|  $\mu$G |  erg cm-3 |  $\mu$G |  erg cm-3 |  $\mu$G |  erg cm-3
West Jet | 60 | $1.4\times 10^{-10}$ | 54 | $1.2\times 10^{-10}$ | 103 | $9.9\times 10^{-10}$
East Jet | 45 | $8.1\times 10^{-11}$ | 50 | $1.0\times 10^{-10}$ | 95 | $8.4\times 10^{-10}$
West Jet K3 | 68 | $1.8\times 10^{-10}$ | 76 | $2.3\times 10^{-10}$ | 125 | $1.4\times 10^{-9}$
West Jet K4 | 72 | $2.0\times 10^{-10}$ | 63 | $1.6\times 10^{-10}$ | 125 | $1.4\times 10^{-9}$

Figure 1: Contour and color scale image of the $z=2.16$ powerful radio galaxy,
J1140-2629 (the ’Spiderweb galaxy’), from Chandra Observatory total X-ray
emission (0.5 keV to 7 keV), convolved with a $1^{\prime\prime}$ Gaussian.
This is the total counts image. Further processing, including AGN subtraction,
will be presented in Tozzi et al. (2021). The contour levels are a geometric
progression in the square root two, starting at 0.8 counts beam-1.  Figure 2:
Contour image of the 2 GHz to 4 GHz VLA image of the Spiderweb galaxy, at a
resolution of $1.3^{\prime\prime}\times 0.6^{\prime\prime}$, major axis north-
south. The contour levels are a geometric progression in a factor two,
starting at 20$\mu$Jy beam-1. Negative contours are dashed. The color scale is
the same Chandra Observatory total X-ray emission (0.5 keV to 7 keV), as shown
in Figure 1.  Figure 3: Same contour image as in Figure 1. The color image is
the ESO narrow band filter image, roughly centered on the Ly$\alpha$ emission
line from the Spiderweb galaxy at $z=2.16$ (Nonino et al., 2022).  Figure 4:
The color scale is the 2 GHz to 4 GHz VLA image of the Spiderweb galaxy, as
shown in Figure 1, saturated in intensity to show the low surface brightness
emission. The white contours are the same as Figure 1 for the 2 GHz to 4 GHz
image. The black contours are the 8 GHz to 12 GHz integrated emission,
convolved to the same resolution as the S band image
($1.3^{\prime\prime}\times 0.6^{\prime\prime}$). The contour levels are a
geometric progression in factor two, starting at 20 $\mu$Jy beam-1.  Figure 5:
Integrated radio spectrum of the Spiderweb radio galaxy.  Figure 6: The color
scale is the spectral index across the 2 GHz to 4 GHz band, at the same
resolution as shown in Figure 1, blanked at $8\sigma$. The contours are the
same as Figure 1.  Figure 7: The color scale is the spectral index between 340
MHz and 3 GHz, at the resolution of the 340 MHz image
($8.0^{\prime\prime}\times 4.3^{\prime\prime}$, major axis north-south). The
black contours are the surface brightness at 340 MHz with geometric contour
levels in factor two, starting at 20 mJy beam-1. The blue contours are the
same as in Figure 1.  Figure 8: Same X-ray total counts color image and 3 GHz
contours, with boxes and ellipses indicating the four regions for the inverse
Compton and Synchrotron analysis in Section 4.2.


# Non-perturbative renormalization of bilinear operators with Möbius domain-
wall fermions
in the coordinate space

JLQCD Collaboration: a,b, G. Cossub, S. Hashimotoa,b, J. Noakib
E-mail  a School of High Energy Accelerator Science, The Graduate University
for Advanced Studies (Sokendai), Ibaraki 305-0801, Japan, b KEK Theory Center,
Institute of Particle and Nuclear Studies, High Energy Accelerator Research
Organization (KEK), Ibaraki 305-0801, Japan

###### Abstract:

We study the non-perturbative determination of the renormalization constants
of flavor non-singlet quark bilinear operators on the lattice. The
renormalization condition is imposed on correlation functions of bilinear
operators in the coordinate space. The results are converted to the value at 2
GeV in the $\rm\overline{MS}$ scheme by a perturbative matching. The
calculation is carried out on gauge configurations generated with the Mobius
domain-wall fermions at two lattice spacings $a^{-1}=2.4$ GeV and $a^{-1}=3.6$
GeV.

##  1 Introduction

There are several methods to renormalize the lattice operators onto those
defined in the continuum renormalization scheme such as the $\rm\overline{MS}$
scheme. Schrödinger functional [1] and RI/MOM [2] approaches are among the
popular choices. The former provides precise and fully non-perturbative
renormalization at the cost of generating dedicated ensembles. The
renormalization condition in the latter approach is imposed on the vertex
function at a fixed gauge and it involves (continuum) perturbative expansion
which is calculated only at one- or two-loop order.

In this work, we investigate the determination of the renormalization
constants by the X-space method, which was originally suggested in [3] and has
been developed in [4, 5]. In this method, the renormalization condition is
imposed on two-point correlation functions of the operators to be
renormalized. Unlike RI/MOM, the X-space method enables us to renormalize by a
gauge invariant quantity. By keeping the distance finite, one can avoid extra
divergences due to contact terms. Another important advantage is that the
perturbative expansion on the continuum side is available to the four-loop
level for the quark bilinear operators [6]. On the other hand, a potential
problem of the X-space method is the requirement for a ”window” where the
continuum perturbation theory can be applied and the lattice calculation with
minimal discretization effect is possible.

In this report, we present the preliminary results of the determination of
renormalization constants of flavor non-singlet quark bilinear operators using
the X-space method. We obtain good precision on the lattices generated with
the $2+1$ flavors of Möbius domain-wall fermions and Symanzik improved gauge
action [7]. We work on $32^{3}\times 64$ lattices at $a^{-1}=2.4$ GeV and
$48^{3}\times 96$ lattices at $a^{-1}=3.6$ GeV, both having matched physical
volume and similar input physical masses in the window of $M_{\pi}:300\sim
500$ MeV.

##  2 Sketch of the X-space method

We impose the renormalization condition on massless correlation functions of
(flavor non-singlet) quark bilinear operators, which are defined as

$\displaystyle\Pi_{\rm SS}(x)=\big{\langle}S(x)S(0)\big{\rangle},\ \ \ \
\Pi_{\rm PP}(x)=\big{\langle}P(x)P(0)\big{\rangle},$ (1)
$\displaystyle\Pi_{\rm
VV}(x)=\sum_{\mu=1}^{4}\big{\langle}V_{\mu}(x)V_{\mu}(0)\big{\rangle},\ \ \ \
\Pi_{\rm
AA}(x)=\sum_{\mu=1}^{4}\big{\langle}A_{\mu}(x)A_{\mu}(0)\big{\rangle},$ (2)

where $S$ and $P$ are scalar and pseudoscalar densities, $V_{\mu}$ and
$A_{\mu}$ are vector and axial-vector currents. Flavor indices are omitted for
simplicity, but they are understood as isospin triplet operators of light
quarks. These bilinear operators on the lattice are to be renormalized onto
the $\rm\overline{MS}$ scheme at 2 GeV, i.e.

${O}_{\Gamma}^{\rm\overline{MS}}|_{2\rm\ GeV}=Z_{\Gamma}^{\overline{\rm
MS}/lat}(2{\rm\ GeV};a){O}_{\Gamma}^{lat}|_{a},$ (3)

for each quark bilinear operators ${O}_{\Gamma}\in\\{S,P,V_{\mu},A_{\mu}\\}$.

Since these correlation functions involve two bilinear operators, the
renormalization condition is given by

${Z_{\Gamma}^{\overline{\rm MS}/lat}(2{\rm\
GeV};a)}^{2}\Pi_{\Gamma\Gamma}^{lat}(x)|_{a}=\Pi_{\Gamma\Gamma}^{\rm\overline{MS}}(x)|_{\rm
2\ GeV},$ (4)

or

$Z_{\Gamma}^{\overline{\rm MS}/lat}(2{\rm\
GeV};a)=\sqrt{\Pi_{\Gamma\Gamma}^{\rm\overline{MS}}(x)|_{\rm 2\
GeV}\over\Pi_{\Gamma\Gamma}^{lat}(x)|_{a}},$ (5)

where $\Pi_{\Gamma\Gamma}^{\rm\overline{MS}}$ and $\Pi_{\Gamma\Gamma}^{lat}$
are correlation functions calculated in the continuum theory or measured on
the lattice.

We need two calculations to obtain the renormalization constants:

  * •

$\Pi_{\Gamma\Gamma}^{lat}(x)|_{a}:$ Correlation functions calculated on the
lattice and taken the chiral limit.

  * •

$\Pi_{\Gamma\Gamma}^{\rm\overline{MS}}(x)|_{\rm 2\ GeV}:$ Massless continuum
correlation functions renormalized at 2 GeV in the $\rm\overline{MS}$ scheme.

In this process, we have to choose $x$ in the window $a\ll|x|\ll\Lambda_{\rm
QCD}^{-1}$ in order to avoid discretization effect and perturbative ambiguity.

##  3 Correlation functions measured on the lattice

Figure 2 shows $x^{2}$-dependence of short-distance lattice correlator of
pseudoscalar channel plotted as a function of $x^{2}$ (red). In this plot, we
also show the two-point correlator calculated at the tree-level, i.e. no
strong interaction (blue). The data show substantial violation of the
rotational symmetry, but apparently their short-distance behavior is rather
precisely reproduced by the free propagator. The simplest example is the two-
points at $(x/a)^{2}=4$ where $(2,0,0,0)$ and $(1,1,1,1)$ are ten-times
different, but they are well reproduced at tree-level.

Figure 1: Correlation function for P channel which are purely measured on the
lattice in the ensemble $48^{3}\times 96,\ \beta=4.35,\ am_{ud}=0.012,\
am_{s}=0.018$ and that in the free system plotted with that in free continuum
theory.  Figure 2:  Correlation function after applying tree-level correction
to that in Fig. 1 with the discrimination of the range $\theta>30^{\circ}$
(open diamond) from $\theta\leq 30^{\circ}$ (filled diamond).
---

We can thus eliminate the bulk of discretization effect by subtracting the
tree-level contribution of the discretization effect as

$\Pi_{\Gamma\Gamma}^{lat}(x)\longrightarrow\Pi_{\Gamma\Gamma}^{lat}(x)-\left(\Pi_{\Gamma\Gamma}^{lat,free}(x)-\Pi_{\Gamma\Gamma}^{cont,free}(x)\right).$
(6)

The results are shown in Fig. 2. We note that the vertical axis is in the
linear scale rather than logarithmic. Reductions of the Lorentz violation is
quite clear. A similar subtraction was introduced, but multiplicatively, i.e.
multiplying the original correlator $\Pi_{\Gamma\Gamma}^{lat}(x)$ by the ratio
$\Pi_{\Gamma\Gamma}^{cont,free}(x)/\Pi_{\Gamma\Gamma}^{lat,free}(x)$. We find
that (6) works slightly better to eliminate the Lorentz violating effect.

Another way to reduce discretization effect is suggested [5] on top of (6).
Let us define $\theta$ as an angle between $x$ and the direction (1,1,1,1).
One may find that the discretization effect becomes more pronounced for large
$\theta$. By inspecting the results, we decided to drop the points of
$\theta>30^{\circ}$. In Fig. 2, the discarded points are shown by open
symbols. The remaining points (blue diamonds) are quite smooth as a function
of $(x/a)^{2}$.

##  4 Correlation functions computed from continuum perturbation theory

Figure 3:  Massless correlators from perturbation theory with $N_{f}=3,\
\Lambda_{\rm QCD}=340$ MeV, scalar channel at 2 GeV in $\rm\overline{MS}$
scheme (left panel) and vector channel (right panel).

Perturbative expansion of the massless correlators in the $\rm\overline{MS}$
scheme

$\displaystyle\Pi_{\rm SS,PP}^{\rm\overline{MS}}(x,\mu)=\Pi_{\rm
SS,PP}^{\rm\widetilde{MS}}(x,\tilde{\mu})={3\over\pi^{4}x^{6}}\left(1+\sum_{n}\widetilde{C}_{n}^{\rm
S}\tilde{a}_{s}^{n}\right),$ (7) $\displaystyle\Pi_{\rm
VV,AA}^{\rm\overline{MS}}(x)=\Pi_{\rm
VV,AA}^{\rm\widetilde{MS}}(x)={6\over\pi^{4}x^{6}}\left(1+\sum_{n}\widetilde{C}_{n}^{\rm
V}\tilde{a}_{s}^{n}\right),$ (8)
$\displaystyle\tilde{a}_{s}={\alpha_{s}^{\rm\widetilde{MS}}(\tilde{\mu}=1/x)\over\pi}={\alpha_{s}^{\rm\overline{MS}}(\mu=2{\rm{e}}^{-\gamma_{E}}/x)\over\pi},$
(9)

are known up to $n=4$ [6]. $\widetilde{C}_{n}^{\rm S,V}$ are perturbative
coefficients and $a_{s}=\alpha_{s}/\pi$. The running of $\alpha_{s}$ is also
known up to four-loop level [8] in the $\rm\overline{MS}$ scheme.

Since the scale $\mu$ and $\tilde{\mu}$ in these formulae depend on $x$, we
need to perform the scale evolution for correlation functions of scalar and
pseudoscalar correlators to those at 2 GeV in the $\rm\overline{MS}$ scheme.
The scale evolution is calculated as

$\Pi_{\rm
SS,PP}^{\rm\widetilde{MS}}(x,\tilde{\mu}_{1})=\left[{c(a_{s}(\tilde{\mu}_{1}))\over
c(a_{s}(\tilde{\mu}_{0}))}\right]^{-2}\Pi_{\rm
SS,PP}^{\rm\widetilde{MS}}(x,\tilde{\mu}_{0}),\ \ \ \
c(x)\equiv\exp\left[\int^{x}{\rm{d}}x^{\prime}{\gamma_{m}(x^{\prime})\over\beta(x^{\prime})}\right],$
(10)

where $\gamma_{m}$ and $\beta$ are the quark mass anomalous dimension and the
QCD beta function, and $c(x)$ is given in [9, 10]. The Ward-Takahashi identity
guarantees that the correlation functions of vector or axial-vector currents
are scale independent.

Figure 4:  Massless correlators from perturbation theory with improved
convergence and chiral limits of lattice correlators times some constants,
which are squares of renormalization constants. These lattice correlators are
measured on $48^{3}\times 96,\ a^{-1}=3.6$ GeV lattice, while a sea quark is
still massive $am_{s}=0.0180$.

As shown in Fig. 3, the convergence of the perturbative series in (7) and (8)
is not sufficiently good to achieve precise determination of the
renormalization constant already at $x\simeq\rm 1\ GeV^{-1}$. At longer
distances, the expansion is even divergent. The perturbative series in
$N_{f}=3$ theory are given as

$\displaystyle\Pi_{\rm
SS}^{\rm\widetilde{MS}}(x,\tilde{\mu})={3\over\pi^{4}x^{6}}(1+0.67\tilde{a}_{s}-16.3\tilde{a}_{s}^{2}-31\tilde{a}_{s}^{3}+497\tilde{a}_{s}^{4}),$
(11) $\displaystyle\Pi_{\rm
VV}^{\rm\widetilde{MS}}(x)={6\over\pi^{4}x^{6}}(1+\tilde{a}_{s}-4\tilde{a}_{s}^{2}-1.9\tilde{a}_{s}^{3}+94\tilde{a}_{s}^{4}).$
(12)

When $\tilde{a}_{s}=0.1$, the expansion for SS is like
$1+0.067-0.163-0.031+0.050$ and does not seem to converge.

In order to avoid this problem, we expand correlators in terms of the coupling
$a_{s}^{*}$ at another scale $\mu^{*}$ according to the BLM prescription [11].
In this way, one can effectively absorbs higher-order contribution from vacuum
polarization effects into lower-orders. The scale for VV (and AA) is given by
$\mu^{*}=\mu{\rm{e}}^{-11/6+2\zeta(3)}\simeq 1.8\mu$. We use the same scale
also for SS (and PP).

The resulting expansion becomes

$\displaystyle\Pi_{\rm
SS}^{\rm\widetilde{MS}}(x,\mu^{*})={3\over\pi^{4}x^{6}}(1+2.9a_{s}^{*}+1.1a_{s}^{*2}-42a_{s}^{*3}+24a_{s}^{*4}),$
(13) $\displaystyle\Pi_{\rm
VV}^{\rm\widetilde{MS}}(x)={6\over\pi^{4}x^{6}}(1+a_{s}^{*}+0.083a_{s}^{*2}-6a_{s}^{*3}+18a_{s}^{*4}),$
(14)

which obviously has much better convergence property. The scale evolution for
$\Pi_{\rm SS}$ to 2 GeV in $\rm\overline{MS}$ is then calculated using (10).

Figure 4 shows the continuum correlators obtained with the above procedure at
each order of the perturbative expansion. It clearly shows a good convergence
even at lower scales $x^{2}\sim 3\rm\ GeV^{-2}$.

Figure 5:  Renormalization constants with the dependence on the renormalization point $x$, where the parameters are same as Fig. 4.  Table 1:  Preliminary results of computing renormalization constants. First error means statistical error and second one means systematic error. $a^{-1}$ [GeV] | $am_{s}$ |  $Z_{S}^{\overline{\rm MS}/lat}$(2 GeV) | $Z_{V}^{\overline{\rm MS}/lat}$
---|---|---|---
$2.4$ | $0.030$ | 1.092(6)(12) | 1.013(2)(7)
$2.4$ | $0.040$ | 1.093(8)(28) | 1.017(2)(9)
$3.6$ | $0.018$ | 0.973(5)(6) | 0.999(3)(5)
$3.6$ | $0.025$ | 0.965(7)(11) | 0.994(6)(5)

##  5 Preliminary results

Figure 4 also shows the lattice data rescaled by a factor such that they agree
with the perturbation theory in the region of $x^{2}\sim\rm 2\ GeV^{-2}$. This
factor corresponds to the renormalization factor squared as given in (4).

Figure 5 shows the renormalization constant calculated from eq. (5) with
dependence on the renormalization point $x$. Since LHS of eq. (5) is
independent of $x$, we should extract renormalization constants from a range
where $x$-dependence of them is approximately absent. Since we employ the
Möbius domain wall fermions, which have excellent chiral symmetry, we can
assume $Z_{\rm S}=Z_{\rm P},\ Z_{\rm V}=Z_{\rm A}$ in the chiral limit, and
extract $Z_{\rm S}$ and $Z_{\rm V}$ as an average of $\\{Z_{\rm S},Z_{\rm
P}\\}$ or $\\{Z_{\rm V},Z_{\rm A}\\}$, respectively. Here the difference
between $Z_{\rm S}$ and $Z_{\rm P}$, or $Z_{\rm V}$ and $Z_{\rm A}$, coming
from higher orders of the operator product expansion is considered as a part
of systematic errors.

Table 1 shows the preliminary results for the renormalization constants. Since
we calculate only at two strange quark masses for each lattice spacings, we do
not take the chiral limit of $m_{s}$. Systematic errors, which are the second
errors in Table 1, are typically $1\%$ or less except for $Z_{\rm S}$ at
$a^{-1}=2.4$ GeV with $am_{s}=0.040$. This $1\%$ precision for $Z_{\rm S}$ is
better than that of other methods, such as RI/MOM at similar lattice spacings.

##  6 Summary

We investigate the non-perturbative renormalization of flavor non-singlet
quark bilinear operators using the gauge invariant X-space method for the
action of $2+1$ Möbios domain-wall fermions and Symanzik improved gauge
action.

Discretization effect of correlation functions measured on the lattice is
substantially reduced by applying the tree-level correction and the angle cut.
Convergence of perturbation theory is dramatically improved by expanding
correlators as a polynomial of coupling $a_{s}^{*}$. The systematic error we
can reach is roughly $1\%$ or better.

###### Acknowledgments.

Numerical simulations are performed on the IBM System Blue Gene Solution at
High Energy Accelerator Research Organization (KEK) under a support of its
Large Scale Simulation Program (No. 13/14-04). This work is supported in part
by the Grant-in-Aid of the Japanese Ministry of Education (No. 26247043) and
the SPIRE (Strategic Program for Innovative Research) Field5 project.

## References

  * [1] M. Lüscher, R. Narayanan, P. Weisz, and U. Wolff, Nucl. Phys. B 384 (1992) 168, hep-lat/9207009.
  * [2] G. Martinelli, C. Pittori, C.T. Sachrajda, M. Testa, and A. Vladikas, Nucl. Phys. B 445 (1995) 81, hep-lat/9411010.
  * [3] G. Martinelli, G.C. Rossi, C.T. Sachrajda, S.R. Sharpe, M. Talevi, M. Testa, Phys. Lett. B 411 (1997) 141, hep-lat/9705018.
  * [4] V. Giménez, L. Giusti, S. Guerriero, V. Lubicz, G. Martinelli, S. Petrarca, J. Reyes, B. Taglienti, and E. Trevigne, Phys. Lett. B 598 (2004) 227, hep-lat/0406019.
  * [5] K. Cichy, Karl Jansen, and Piotr Korcyl, Nucl. Phys. B 865 (2012) 268, arXiv:1207.0628 [hep-lat].
  * [6] K.G. Chetyrkin and A. Maier, Nucl. Phys. B 844 (2011) 266, arXiv:1010.1145 [hep-ph].
  * [7] J. Noaki, S. Aoki, G. Cossu, H. Fukaya, S. Hashimoto and T. Kaneko, (Now applying to PoS LATTICE 2014 (2014) 069).
  * [8] K.G. Chetyrkin, B.A. Kniehl, and M. Steinhauser, Phys. Rev. Lett. 79 (1997) 2184, hep-ph/9706430.
  * [9] K. G. Chetyrkin, Phys. Lett. B 404 (1997) 161, hep-ph/9703278.
  * [10] J.A.M. Vermaseren, S.A. Larin, and T. van Ritbergen, Phys. Lett. B 405 (1997) 327, hep-ph/9703284.
  * [11] S. J. Brodsky, G. P. Lepage and P. B. Mackenzie, Phys. Rev. D 28 (1983) 228.