Simulating the transport of relativistic electrons and magnetic fields injected by radio galaxies in the intracluster medium
AAstronomy & Astrophysics manuscript no. main ©ESO 2021February 9, 2021
Simulating the transport of relativistic electrons and magneticfields injected by radio galaxies in the intracluster medium
F. Vazza , , , D. Wittor , , G. Brunetti , M. Brüggen Dipartimento di Fisica e Astronomia, Universitá di Bologna, Via Gobetti 93/2, 40122, Bologna, Italy Hamburger Sternwarte, University of Hamburg, Gojenbergsweg 112, 21029 Hamburg, Germany Istituto di Radioastronomia, INAF, Via Gobetti 101, 40122, Bologna, ItalyAccepted ???. Received ???; in original form ???
ABSTRACT
Radio galaxies play an important role in the seeding of cosmic rays and magnetic fields in galaxy clusters. Here, we simulate theevolution of relativistic electrons injected into the intracluster medium by radio galaxies. Using passive tracer particles added tomagnetohydrodynamical adaptive-mesh simulations, we calculate the evolution of the spectrum of relativistic electrons taking intoaccount energy losses and re-acceleration mechanisms associated with the dymamics of the intracluster medium. Re-acceleration canoccur at shocks via diffusive shock acceleration, and in turbulent flows via second-order Fermi re-acceleration. Relativistic electronsfrom radio galaxies are found to fill the intracluster medium over scales of several 100 Myr, and they create a stable reservoir of fossilelectrons which remains available for further re-acceleration by shock waves and turbulent gas motions. In the near future, deep radioobservations (especially at low frequencies) are likely to probe such mechanisms in galaxy clusters. alaxy clusters, ICM, radio galaxies, radio emission
1. Introduction
Galaxy clusters are the largest reservoir of relativistic plasma inthe Universe, in which relativistic particles are subject to ac-celeration and energy losses and are confined by turbulent mag-netic fields (e.g. Berezinsky et al. 1997; Sarazin 1999; Brunetti& Jones 2014; Bykov et al. 2019). Gamma-ray observations (e.g.Ackermann et al. 2014) and radio observations (e.g. van Weerenet al. 2019) have provided important constraints on the complexlife cycle of relativistic plasma in the intracluster medium (ICM).Radio-bright, bipolar outflows from active galactic nuclei(AGN) are commonly found in clusters of galaxies. Radio galax-ies play an important role as they have the potential to releaselarge quantities of relativistic plasma and magnetic fields into theICM (e.g. Völk & Atoyan 1999). Many recent numerical simu-lations have investigated the dynamics of relativistic jets as theyexpand into the ambient medium, and the role of MHD instabili-ties in determining the jet morphology and stability (e.g. Normanet al. 1982; Bodo et al. 1998; Perucho & Martí 2007; Mignoneet al. 2005, 2010; Hardcastle & Krause 2014; Massaglia et al.2016; Borse et al. 2020).Cosmological simulations predict that radio-mode feedbackis crucial for shaping the thermodynamic properties of the gas ingalaxy groups and clusters at low redshifts (e.g. Puchwein et al.2008; Booth & Schaye 2009; McCourt et al. 2011; Fabjan et al.2010; Tremmel et al. 2017). At higher redshifts ( z ≥ −
3) mostof the mass growth of SMBH occurs in the radiatively efficientquasar mode, with an accretion rate larger than a few percent ofthe Eddington rate (e.g. Sijacki et al. 2007).Observationally, the duty cycle of radio galaxies in theICM is still poorly constrained. Observations of nearby massive
Send offprint requests to : E-mail: [email protected] In this paper we will use the terms relativistic particles and cosmic-ray particles interchangeably. galaxy clusters by Bîrzan et al. (2020) concluded that only abouthalf of such systems show clear evidence for radio-mode AGNfeedback. Several deep LOFAR observations have discoveredsteep-spectrum, filamentary and distorted radio structures, oftenconnecting old tails of radio galaxies with diffuse radio emissionacross a wide range of spatial scales (Wilber et al. 2018; Botteonet al. 2020b; Mandal et al. 2020). Very recent deep JVLA ob-servations Gendron-Marsolais et al. (2021) have also unveiledcomplex substructures in the wake of radio galaxy NGC 1272 inthe Perseus clusters. Its morphology and spectrum support thepossibility that fossil plasma ejected by AGN jets is being re-accelerated by shear and compressive motions in the ICM andthat the fossil plasma fuels the the mini-halo in the Perseus clus-ter. The emission from tailed radio galaxies is sometimes linkedto diffuse radio emission from the ICM, both in the form of radiohalos (Wilber et al. 2018) or radio relics (e.g. Markevitch et al.2005; Owen et al. 2014; Bonafede et al. 2014; van Weeren et al.2017; Stuardi et al. 2019). In most cases, a volume-filling dis-tribution of fossil relativistic electrons is required to explain theobserved radio power (e.g. Brunetti et al. 2001; Markevitch et al.2005; Kang et al. 2012; Pinzke et al. 2013; Vazza & Brüggen2014; Vazza et al. 2015; Brunetti & Jones 2014; Kang 2018; Bot-teon et al. 2020a). The recent ASKAP/EMU observation of themerging A3391-3395 pairs of galaxy clusters has also detected athree times higher sky density of giant radio galaxies (Brüggenet al. 2020) than previously assumed. This suggests that radiojets can have a larger role in seeding cosmic rays far from theirsource.Radio galaxies can also play a role in the genesis of mag-netic fields in galaxy clusters (Kronberg et al. 1999; Furlanetto& Loeb 2001; Xu et al. 2009), as they can seed primordial mag-netic fields (e.g. Vachaspati 2020, for a recent review). Outsideof the virial regions of clusters, the impact of galaxies on ex-tragalactic magnetic fields is expected to be detectable with thenext generation of deep radio surveys, both in total intensity and
Article number, page 1 of 21 a r X i v : . [ a s t r o - ph . H E ] F e b &A proofs: manuscript no. main polarisation (e.g. Vazza et al. 2017; Locatelli et al. 2018). Thedynamical interaction between radio jets and the ICM occurs ona wide range of spatial scales (1 kpc ≤ L ≤
100 kpc). A num-ber of simulations have given valuable insights into the interac-tion between AGN and the surrounding ICM, mostly based onhydrodynamics and starting from the point where the relativis-tic plasma has reached approximate pressure equilibrium withthe surrounding ICM and starts to be dominated by buoyancy(e.g. Churazov et al. 2001; Brüggen & Kaiser 2002; Mathews &Brighenti 2007; Heinz et al. 2006; McCarthy et al. 2010; Gaspariet al. 2011; Yang et al. 2012, 2019).On the other hand, a few papers have focussed on the di-rect injection of magnetic fields by AGN, and on their impact onthe evolution of the ICM. In particular, Xu et al. (2009) and Xuet al. (2011) included magnetised outflows from radio galaxiesin cosmological
ENZO simulations, and studied the build-up ofcluster magnetic fields from the injection of individual AGN jets.Mendygral et al. (2012), simulated the effect of "cluster weather"on radio lobes by injecting magnetised jets into a galaxy clusterextracted from a (Smoothed Particle Hydrodynamics) cosmolog-ical simulation, and resimulating a cluster cutout region with agrid-based MHD method, in order to evolve radio jets and theirmagnetic properties at high resolution. Using the AREPO code,Bourne & Sijacki (2020) performed high-resolution simulationsof jets from AGN, and produced realistic X-ray and radio prop-erties (albeit assuming a distribution of magnetic fields and rela-tivistic electrons in post-processing).In this new work we devote more computational resourcesto the simulation of the properties of non-thermal componentsinjected by jets, while neglecting the effect of gas cooling. Ourmain focus is on the evolution of radio jets, as well as their mag-netic fields and relativistic electrons. In particular, we model theevolution of relativistic electrons subject to ageing and accelera-tion processes. We rely on passive tracers, which in turn, do notallow us to study the interplay of a relativistic plasma with theICM, as is done in other work (O’Neill & Jones 2010; Mendygralet al. 2012; Nolting et al. 2019a,b). However, our simulations of-fer an unprecedented view of the large-scale circulation and ad-vection of fossil electrons on Mpc scales, in a regime scarcelyaffected by radiative cooling and other processes neglected inthis work. They also allow us to test several re-acceleration pro-cesses together with a realistic model of the evolving ICM.Our paper is structured as follows: in Sec. 2 we present thecosmological simulations and numerical methods employed inthis paper, while in Sec. 3 we give our results on the evolutionof radio sources (3.1), on the impact of radio galaxies on theICM (3.2), and on the the energy evolution relativistic electrons(3.3). Limitations of our approach are discussed in Sec. 4 andour conclusions are given in Sec. 5.
2. Simulations
ENZO -simulations
We used the cosmological
ENZO -MHD (Bryan et al. 2014) codeto produce realistic simulations of the formation of a galaxygroup and of the thermal and non-thermal feedback from radiogalaxies within it.We use nested grids to provide a uniform resolution at thehighest refinement level. Each simulation covers a root-grid of(50 Mpc / h ) and is sampled with 128 cells. Using MUSIC (Hahn & Abel 2011), four additional nested regions with in-creasing spatial resolution were nested, until the innermost (4Mpc/ h ) region where the cluster form is uniformly covered at a Table 1.
Input parameters for the SMBH models in our simulations, andreference parameters for the AGN jets, measured as an average at theend of each jet’s life ( t jet ).Parameter Run2 Run1 z jet M BH [M (cid:12) ] 10 α Bondi
10 10 ˙ M BH [M (cid:12) / yr] 10 − − t jet [Myr] 10 200 L BH [erg/s] 1 . · . · E kin , jet [erg] 2 . · . · E th , jet [erg] 1 . · . · E mag , jet [erg] 1 . · . · B av , jet [ µ G] 1 . . B max , jet [ µ G] 3 . . v av , jet [km / s] 630 1267 ∆ x = . / h uniform resolution. During the simulation, twoadditional level of mesh refinement were added using a localgas/DM overdensity criterion ( ∆ ρ/ρ ≥ ≈ .
86 kpc, typically inthe cluster core region. As a result of our nested grid approach,the mass resolution for dark matter in our cluster formation re-gion is of m DM = . · M (cid:12) per dark matter particle, forthe highest resolution particles that are used to fill the innermostAMR level since the start of the simulation.The MHD solver employs the Local Lax-Friedrichs (LLF)Riemann solver to compute the fluxes in the Piece-wise LinearMethod (PLM). We initialised a simple uniform magnetic fieldat z =
50 with a value of B = . M ≈ . · M (cid:12) at z = . M ≈ . · M (cid:12) at z = . M ≈ · M (cid:12) at z =
1. The green colors in Fig. 1 helpin visualising its main evolutionary events. At z = / z ∼ . z = . z = .
1, following the accretion of a second massive companion(entering to the right of the panel for z = . h = . Ω Λ = . Ω M = .
308 and Ω b = . ENZO
The simulation of SMBHs follows the formulation by Kim et al.(2011), released in the public version of
ENZO and supplementedwith a few ad-hoc prescriptions to include the magnetic feedbackfrom radio jets. Given the absence of radiative cooling, our runslack a self-regulating mechanism to switch the feedback cycleon and off (e.g. Gaspari et al. 2011; Rasia et al. 2015; Ricarteet al. 2019, for a few examples). The thermal structure of the
Article number, page 2 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 1.
Maps of all tracers evolved in our run Run0 (left) and Run1(right) at z = z = . z = . z = .
1. The colours representthe projected gas density (with a log stretching). The line of sight isparallel to the jet axis, to give a sense of the lateral expansion of tracersrelative to the jet (which points towards the observer in this case). ICM around the SMBH is thus not entirely realistic, which alsoimplies that typical parameters for the launching of jets fromAGN must be adjusted to reproduce realistic radio galaxies.After extensive testing of the possibility of generating ki-netically dominated jets as implemented in other simulations(Dubois et al. 2010; Gaspari et al. 2011), we resorted to purely thermal energy feedback around the SMBH. This is motivatedby the fact that even the interaction of a fast jet with a T ICM ≥ K cluster core will not produce realistic shock heating and dynamics. This is in conflict with usual models for AGN jets incool-core atmospheres (where T ≤ K during the jet launch-ing stage e.g. Brighenti & Mathews 2000; Brüggen & Kaiser2002). While a pure thermal feedback would also be in conflictwith recent Sunyaev-Zeldovich measurements of AGN-blowncavities that show a substantial amount of non-thermal supportinside these cavities (Abdulla et al. 2019), a large fraction of thefeedback energy in our model is also provided by the magneticand kinetic energy produced by our bipolar magnetic feedbackmodel (see below). This approach allowed us to relax the other-wise tight constraint on the time-step in our simulation, set bythe typically large jet velocities ( v j ∼ km / s) in the vicinity ofthe SMBH.In summary, our approach is suitable to investigate the im-pact of radio galaxies on the gas in large-scale structures, pro-vided that we focus on (cid:29)
10 kpc scales and that dynamical fea-tures on smaller scales can be neglected.Following the implementation by Kim et al. (2011), eachSMBH particle is allowed to accrete gas, based on the spher-ical Bondi-Hoyle model. Due to the lack of spatial resolution,it is necessary to simulate the mass growth of SMBH with ad-hoc parametric formulations. For example, the commonly usedEddington-limited Bondi–Hoyle accretion generally has to beboosted by an efficiency parameter α ≥
1, so that˙ M BH = min πα G M ρ B c ; 4 π GM BH m p (cid:15) r σ T c , (1)where M BH is the SMBH’s mass, c s is the sound speed of thegas at the SMBH’s location, ρ B is the local gas density, m p is theproton mass, σ T is the Thomson scattering cross-section, and (cid:15) r is the radiative efficiency of the SMBH, so that L BH = (cid:15) r ˙ M BH c is the bolometric luminosity of the SMBH.Given our resolution, we cannot resolve the Bondi radiusand simulate how the SMBH gravitationally influence their sur-roundings. Neither can we resolve the multi-phase interstel-lar medium around the host galaxy, given that we do not in-clude cooling and star formation. Hence we adopted a factor, α Bondi >
1, (e.g. Booth & Schaye 2009; Gaspari et al. 2012;Tremblay et al. 2016) that parametrises the true mass accretionrate onto the SMBH. After testing we used two different val-ues of α Bondi depending on the epoch of jet injection, which alsodepends on the spatial resolution at the location of SMBH seedparticles. For models in which radio jets start at z = α , i.e. α Bondi = , in order to producejets leading to a realistic radio morphology. For the simulationin which our radio jets are launched at z = . α Bondi =
10) is foundto produce a jet morphology comparable to real radio galaxiesin a dense atmosphere, as in, e.g. Booth & Schaye (2009). Wealways assume that the gas temperature in the accretion disc is T disc = K.The
ENZO model for SMBH feedback by Kim et al. (2011)deposits energy at the maximum numerical resolution in theform of thermal energy. A SMBH particle releases thermalfeedback on the surrounding gas, as an extra thermal energyoutput from each black hole particle, assuming an efficiency (cid:15) r (cid:15) BH = E jet / ( ∆ M ∆ tc ) between the accreted mass, ∆ M , dur-ing the timestep ( ∆ t ) and the feedback energy E jet . The ther-mal feedback energy is distributed to the neighbouring 27 gas Article number, page 3 of 21 &A proofs: manuscript no. main cells around the SMBH. We assumed
ENZO ’s default efficiencyparameters, (cid:15) r = . (cid:15) BH = .
05 (e.g. Kim et al. 2011;Bryan et al. 2014), which also yielded a good match with ob-served galaxy clusters scaling relations as found in previouswork (Vazza et al. 2013; Vazza et al. 2017).The SMBHs also inject magnetic fields in the form of mag-netic dipoles (2 × ± z -direction from the SMBH). This very simpletopology is made necessary by the limitation of spatial resolu-tion. More sophisticated choices, e.g., considering helical mag-netic fields (Ruszkowski et al. 2007) require higher spatial reso-lution, i.e. ≤ . E B , jet = (cid:15) B , jet E jet , with (cid:15) B , jet = .
1. How-ever, this does not directly correspond to the magnetisation ofjets since the thermal energy is released isotropically around theSMBH, while the magnetic energy is directed . For this reason,the typical magnetic energy of our jets (see e.g. Tab. 1) is a factor ∼
10 larger than the kinetic and thermal energy assigned to cellsat the poles of SMBHs, where jets are launched, making theminitially magnetically dominated (e.g. Xu et al. 2009; Massagliaet al. 2019).The scenarios we explored in this work are: – Run0 : No feedback from SMBH is activated. This simulationis essentially a standard non-radiative MHD simulation, inwhich the ICM magnetic fields are the result of compressionand, partly, of small-scale dynamo amplification of the initial B = . . – Run1 : In this case the SMBH at the centre of the most mas-sive halo in the high-resolution region is activated at redshift z =
1, with an initial mass of M BH = M (cid:12) and remains ac-tive for ≈
200 Myr, releasing its thermal and magnetic feed-back into the ICM. We set α Bondi = to boost the accretionrate. The initial magnetic seed field is the same as in Run0.Hence, the magnetisation of the ICM is the combined resultof seed field amplification and AGN feedback. – Run2 : as in Run0 until z = .
5, then a M BH = M (cid:12) SMBHis placed at the centre of the most massive halo in the high-resolution region (which is inside an already formed galaxygroup). The SMBH releases thermal/magnetic feedback for ≈
10 Myr. As in Run1, the magnetisation of the ICM stemsfrom the combination of seed field amplification and AGNfeedback. In this case, we use α Bondi =
10 to boost the (unre-solved) mass accretion rate onto the SMBH. The initial mag-netic seed field is the same as in Run0.Table 1 gives more details of the list of parameters describingthe jet launching in all runs. These are measured within the jetvolume, after a time t jet . Possible caveats of our method are dis-cussed in Sec.4. CRaTer -simulations of Lagrangian tracer particles
We use the Lagrangian code
CRaTer to follow the spatialevolution of the cosmic-ray electrons in our runs. In previouswork,
CRaTer has been used to study cosmic rays and turbulence Note that we oriented the poles always along the same coordinateaxis. In comparison with our most resolved simulations in Vazza et al.(2018) and Domínguez-Fernández et al. (2019), the fraction of volumewhere the Alfvén scale is well-resolved is smaller, owing to the lessaggressive AMR strategy employed here. Hence, small-scale dynamoamplification can develop only in a smaller volume in galaxy clusters (e.g. Wittor et al. 2016, 2017b,a). For detailsof the implementation, we point to these references.In post-processing, we injected ∼ · particles in all runs,starting from the highest-density peaks of the simulation at z = z = . m trac = · M (cid:12) , and evolved them using all snap-shots of the simulation ( ∼
120 for Run1 and ∼
100 for Run2).The various grid quantities, e.g. density, velocity, are assignedto the tracers using a cloud-in-cell (CIC) interpolation method.Further details on the full procedure implemented to advect trac-ers in our simulations, in order to properly sample the advectionof gas matter performed by the
ENZO (Eulerian) calculation, werefer to Wittor (2017) and Wittor et al. (2016) .Fig. 1 shows the spatial evolution of all tracers in runs Run0and Run1, where tracers were initialised at z =
1. The mapsare projected along the jet axis to show the dispersal of tracersin directions perpendicular to the jets. Although the large-scaledynamics of the cluster is the same in both runs, the impulsivefeedback from the radio galaxy, even with a single event, appearsto increase the filling factor of tracer particles as a function oftime (see a longer discussion in Sec. 3.3).While all tracers (also in Fig. 8) were injected at the highestdensity peaks, only ∼ tracers were used to compute the radiospectra. Those tracers were selected based on the values of mag-netic field strength ( ≥ µ G) and the magnitude of the velocitywith respect to the SMBH ( | v | ≥
600 km / s) at the moment of thefirst injection of AGN jets, so that only the gas directly entrainedby jets is initially enriched with cosmic rays.After their injection, we used a temperature-jump basedshock finder to detect when tracer particles traverse a shock (Sec.2.2 in Wittor et al. 2017b). For each detected shock they computethe Mach number according to: M = (cid:115) T new T old ρ new ρ old + . (2)Our tracer particles also keep track of the local fluiddivergence, ∇ · v , and of the fluid vorticity, ∇ × v , whichserve as proxies for the local turbulence experienced by thetracer particles. Hence, we can estimate the effects of com-pression/rarefaction as well as turbulent re-acceleration onthe electron energy distribution (Sec. 2.3). In particular, weused the gas vorticity to estimate the solenoidal turbulenceexperienced by particles, σ v = |∇ × v | l scale , where for simplicitywe used the same fixed reference scale of l scale =
27 kpc usedto compute vorticity via finite differences (i.e. 3 cells on thehigh-resolution mesh). The vorticity is used to compute theturbulent re-acceleration model outlined below, as well as the(solenoidal) turbulent kinetic energy flux. In our model forturbulent re-acceleration (Sec. 2.3) the kinetic energy flux isthe key input. We assume that this is constant across the Kol-mogorov turbulent cascade, provided the turbulence injectionscale is always ≥ l scale . The assumption of a turbulent cascadeclose to a Kolmogorov model, and of an injection scale ≥ Article number, page 4 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 2.
Evolution from z = .
49 to z = .
39 of the radio emission (50 MHz) from jets in our Run2. The red contours (spaced in square rootof the emission) are drawn to better visualise the location of emission peaks in the jet structure. The colorbar gives the emitted power per pixel( = .
86 kpc comoving).
Fig. 3.
Evolution from z = .
99 to z = .
89 of the radio emission (50 MHz) from jets in our Run1. The red contours (spaced in square rootof the emission) are drawn to better visualise the location of emission peaks in the jet structure. The colorbar gives the emitted power per pixel( = .
86 kpc comoving).
We solve the time-dependent diffusion-loss equation of relativis-tic electrons represented by tracer particles, using the standardChang & Cooper (1970) finite-difference scheme, implementedin the programming language Julia (https://julialang.org). Wetypically used N b = p min ≤ p ≤ p max momentum range (with P = γ m e v and p = P / ( m e c )is the normalised momentum of electrons and p min = p max = · (hence d p = N ( p ), can be computed separately for each tracer parti-cle: ∂ N ∂ t = ∂∂ p (cid:34) N (cid:32)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p τ rad (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) + (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p τ c (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) + p τ adv − (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p τ acc (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:33)(cid:35) + ∂ ∂ p (cid:16) D pp N (cid:17) , (3) where D pp is the particle diffusion coefficient in momentumspace. We can thus define p = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p τ rad (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) + (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p τ c (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) + p τ adv − (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p τ acc (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) . (4)The loss timescales for the radiative, Coulomb and expansion(compression) processes are give by the following formulae,adapted from Brunetti & Jones (2014): τ rad = γ/ (cid:20)(cid:16) B . µ G (cid:17) + (1 + z ) (cid:21) , (5) τ c = (cid:40) n / − γ/ (cid:32) . +
175 ln (cid:32) γ/ / − (cid:33)(cid:33)(cid:41) − (6) Article number, page 5 of 21 &A proofs: manuscript no. main and τ adv =
951 Myr ∇ · v / − , (7)where n measured in [cm − ], B in [ µ G] and ∇ · v in [1 / s ]. We ne-glect bremsstrahlung losses since their timescale is significantlylarger than the ones of all other loss channels for the ICM phys-ical condition considered here.Concerning energy gains for relativistic electrons, we con-sidered the contribution from Fermi-I type acceleration fromshocks (i.e. via Diffusive Shock Acceleration) and from Fermi-II type acceleration (i.e. from adiabatic-stochastic-acceleration).The shock kinetic energy flux that is converted into the acceler-ation of cosmic rays (e.g. Ryu et al. 2003), is given by Ψ CR = ξ e η ( M ) ρ u V s dx t , (8)where ρ u is the pre-shock gas density, V s is the shock veloc-ity and the combination ξ e η ( M ) gives the cosmic-ray accel-eration efficiency, which comprises a prescription for the en-ergy going into cosmic rays, η ( M ), and the ratio of accelera-tion rates of electrons to that of protons, ξ e . For convenienceof implementation, we use the polynomial approximation givenby Kang & Jones (2007) for η ( M ), which includes the effectsof finite Alfvén wave drift and wave energy dissipation in theshock precursor. dx t is the surface associated to each tracers,which we can compute considering that dx t = dx / n tracers is theinitial volume associated to every tracer at the epoch of theirinjection ( n tracer being the number of tracers in every cell) and dx t ( z ) = dx t · ρ t /ρ ( z ) is the relative change of the volume as-sociated to each tracer as a function of z , based on the ratio be-tween the density at injection, ρ and the density of cells whereeach tracer sits as a function of redshift, ρ ( z ).For the electron-to-proton ratio, ξ e ∼ − would be com-monly used to model strong supernova remnant shocks (e.g.Uchiyama et al. 2007). However, the exact ratio is extremelyuncertain for weak shocks, and several sophisticated particle-in-cell simulations have been performed over the last years toinvestigate this (e.g. Riquelme & Spitkovsky 2011; Guo et al.2014a,b; Xu et al. 2020). In order to compare to previous workon fossil electrons, we strictly follow the previous work byPinzke et al. (2013) and link the injection of fresh electronsto the possible one protons injected by DSA, by requiring anequal number density of cosmic-ray electrons and protons abovea fixed injection momentum, which yields ξ e = ( m p / m e ) (1 − δ inj ) / .This approach gives ξ e ∼ − for an injection spectral indexof δ inj ≈ .
3, in line with the injection spectral index of localGalactic supernova remnants.For simplicity, we also neglect a possible dependence of theshock acceleration efficiency on the local magnetic field topol-ogy, which could change the cosmic ray content of galaxy clus-ters (e.g. Wittor et al. 2020; Banfi et al. 2020; Guo et al. 2014a;Xie et al. 2020). We inject relativistic electrons with a momen-tum distribution that follows a power-law (e.g. Kardashev 1962;Sarazin 1999): Q inj ( p ) = K inj p − δ inj (cid:32) − pp cut (cid:33) δ inj − , (9)in which the initial slope of the input momentum spectrum, δ inj ,follows from the standard Diffuse Shock Acceleration predic-tion, δ inj = M + / ( M − p cut is the cut-off momentum, which is the defined for every shocked tracer as the maximummomentum, beyond which the radiative cooling timescale getsshorter than the acceleration timescale, τ DSA : τ DSA = D ( E ) V s · r ( r + r − , (10)where r is the shock compression factor and D ( E ) is the electrondiffusion coefficient as a function of energy (e.g. Gabici & Blasi2003). The specific energy-dependent value of D ( E ) is poorlyconstrained and it depends on the local turbulent conditions ofthe plasma undergoing shocks. It is however critical to set themaximum energy that can reached by shock acceleration (e.g.Kang et al. 2012). Nevertheless, this is not an issue for our work,as all plausible choices of D ( E ) in Eq. 10 give an accelerationtimescale which is many orders of magnitude smaller than thetypical cooling time of radio emitting electrons, whose momen-tum distribution at injection can be assumed to follow a powerlaw within out momentum range of interest.This also motivates the fact that we can model shock injec-tion by DSA by adding the newly created population of particlesacross timesteps (see Eq. 17 below), without integrating a sourceterm as needed for the much slower re-acceleration by turbulence(see below).The expression for the normalisation factor, K , can thus bederived by equating the cosmic ray energy flux crossing eachtracer volume element, and the product between the total energyof cosmic rays advected with a post-shock velocity ( v d ): Ψ CR dx t = v d E CR , (11)where v d is the downstream (post-shock) velocity E CR = (cid:90) p cut p inj Q inj ( p ) T ( p ) d p , (12)with Q inj ( p ) defined as above and T ( p ) = ( (cid:112) + p − m e c .The integration yields (e.g. Pinzke et al. 2013) E CR = K inj m e c δ inj − · (cid:34) B x (cid:32) δ inj − , − δ inj (cid:33) + p − δ inj cut (cid:32) (cid:113) + p − (cid:33) , (cid:35) (13)where B x ( a , b ) is the incomplete Bessel function and x = / (1 + p ). Beside the direct injection of relativistic electronsby shocks, we also include re -acceleration by shocks waves (e.g.Markevitch et al. 2005; Kang & Ryu 2011; Kang et al. 2012).According to DSA, the input particle spectrum, N ( x ), becomes N ( p ) = ( δ inj + · p − δ inj (cid:90) pp min N ( x ) x δ + dx , (14)where δ inj is the local slope in momentum space and δ is theslope corresponding to the new shock, according to DSA.Turbulence is thought to play a major role in producing giantradio halos in galaxy clusters, in connection with mergers (e.g.Brunetti & Jones 2014, and references therein). Some residuallevel of turbulence is expected to be present at all times (e.g.Vazza et al. 2011; Angelinelli et al. 2020), caused by minormergers (e.g. Norman & Bryan 1999; Cassano & Brunetti 2005; Article number, page 6 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Iapichino & Niemeyer 2008; Wittor et al. 2017b), AGN feed-back (e.g. Gaspari et al. 2011; Bourne & Sijacki 2020) and coolcore sloshing (e.g. ZuHone et al. 2010). The solenoidal compo-nent of turbulence is prominent in clusters (e.g. Miniati 2014;Vazza et al. 2017). According to Brunetti & Lazarian (2016),solenoidal turbulence can re-accelerate particles via stochasticinteraction with diffusing magnetic field lines in super-Alfvenicturbulence developed in the ICM. In addition to radio halos, thismechanism has also been used in models of Gamma-Ray Bursts(Xu & Zhang 2017) and radio bridges between galaxy clusters(Brunetti & Vazza 2020; Bonafede et al. 2021). The implied ac-celeration timescale is thus: τ ASA = p D pp = . × Myr L / (0 .
5) B (cid:112) n / − δ V / , (15)where L is measured in Mpc and δ V is the local gas turbulentvelocity within the scale L , measured in cm/s. Based on this for-mula, it can be seen that for typical ICM conditions, only for(solenoidal) turbulent rms velocities of δ V ≥ rmkm / s canthis mechanism produce acceleration on timescales of ≤ Gyr(e.g. Sec.3.3.1).The interaction between particles and turbulence is a stochas-tic process that can be thought of as the combination of twoeffects: First, electrons are systematically accelerated on atimescale given by Eq. 15 and, second, their energy is changedstochastically (Eq. 3). Starting from a monoenergetic initial dis-tribution of electrons, these two effects result in the progressivewidening of the momentum distribution of electrons, with a me-dian value determined by the systematic acceleration (capturedby our solver).Finally, the injection of relativistic electrons from radiogalaxies assumes that a new population of electrons is generatedat each tracer located in a jet launching cell, as K r (cid:90) p max p min p − δ r d p = φ e m trac , (16)where we assume an injection spectrum of δ r = φ e so that the number of all injected relativistic electrons is 10 − ofthe total number of thermal electrons for each tracer.In the following, we neglect the stochastic term, D pp , in Eq.3, the numerical solution is obtained with the Chang & Cooper(1970) finite difference scheme: N ( p , t + dt ) = N ( p , t ) / dt + N ( p + d p , t + dt ) p / dt + p / d p + Q inj , (17)where in the splitting-scheme for finite differences we assumed N ( p + d p / = N ( p + d p ) and N ( p − d p / = N ( p ), and Q inj accounts for the injection by radio galaxies or shocks, when-ever present, regarded as an instantaneous process owing to thetimescales much shorter than the time step of our integration (asclear from Eq. 10).With this simple approach to solve Fokker-Planck equations,we can only model the systematic acceleration by turbulent re-acceleration. In order to model the stochastic acceleration, whichcan lead to higher electron energies, one has to resort to morecomplex solvers (e.g. Donnert & Brunetti 2014, for an applica-tion to numerical simulations).In order to quantify the effects of different mechanisms onthe observed properties of radio emission, we will compare theoutcome of three models applied to our simulated electron evo-lution: – a model in which relativistic electrons get injected by ra-dio jets (assuming an initially an unbroken power-law energydistribution with slope δ r = − .
0) and are subsequently onlysubject to energy losses as they are advected in the ICM; – a model in which, beside the injection of relativistic elec-trons by jets and their continuous cooling during advection,new relativistic electrons can be injected at shocks waves(via the Q inj term in Eq. 17), as well as by considering there-acceleration by weak shocks ( M ≤ – a model in which, beside the two acceleration mechanismsconsidered above and the continuous cooling, electrons canalso be re-accelerated by the turbulent re-acceleration viasecond-order acceleration (Eq. 15).By comparing the radio properties of cosmic-ray electrons inthese three models, we will infer the most likely distribution offossil electrons, injected by radio galaxies.The main routines to solve for the evolution of relativisticelectrons are written and parallelized in the Julia (v1.4.0) lan-guage. This allows us to simultaneously evolve three differentrealisations of the momentum spectra over N b = ∼ tracers with ∼
100 s / step, using 16 threads on 8 Intel I9cores.
3. Results
Using cosmological, magnetohydrodynamical simulations, weaddressed the following questions: – What is the impact of jets from radio galaxies on the thermo-dynamic, kinematic and magnetic properties of the ICM? – How much do these effects influence the distribution of rela-tivistic electrons in the ICM? – How do the properties of relativistic electrons change de-pending on loss and (re)acceleration mechanisms arisingfrom the interaction with the ICM?In the following sections, we will discuss these question inthe light of our results in turn.
Radio galaxies have been historically separated into two classes,based on the radio surface brightness of their cocoons (Fa-naroff & Riley 1974). While Fanaroff & Riley type II (FRII)sources show well-defined jet termination shocks around edge-brightened cocoons, type I sources (FRI) have peaks of sur-face brightness closer to the core. Both class of FR sourcesare thought to initially expand supersonically, with the jet-environment interaction likely to have a decisive role in deter-mining when FRII sources evolve into FRI systems, followingthe disruption of their jets (e.g. Kaiser & Alexander 1997; Turner& Shabala 2015). Moreover, many radio sources associated withlow-z galaxies and detected by recent surveys are compact andsmall in size, and can be regarded as one of the extremes (FR0)of a continuous population of radio galaxies, with a broad dis-tribution of sizes and luminosities (e.g. Capetti et al. 2019, seehowever Hardcastle & Croston 2020 for a different take on FR0sources). A serial version of the code, containing all most important fea-tures and a sample sequence of tracer particle data is publiclyavailable here: https://github.com/FrancoVazza/JULIA/tree/master/CR_solver_pub . Article number, page 7 of 21 &A proofs: manuscript no. main
Fig. 4.
Evolution of the X-ray emission and of the radio emission (log contours) for the Run1 model, in which a radio galaxy is activated at z = .
5. The epochs in the panels are z = . z = . z = . z = . z = .
31 and z = .
16. No observational cut is applied to the maps,and the radio data have been smoothed to a × δ x = .
86 kpc).
Our simulated radio sources appear to have broad morpho-logical similarity with FRII-type sources, but in both cases showa dynamic evolution. Figures 2-3 show the evolution of radioemission (at 50 MHz) in the first evolutionary stages ( ≤ ∼ erg / s / Hz 100 Myr after their injection,and it fades to ≤ − erg / s / Hz 200 Myr after injection,due to the combination of radiative cooling and adiabatic lossesexperienced by lobes.Jets in Run2 are quickly distorted and bent by the interactionwith the surrounding turbulent ICM, while they remain more col-limated in Run1, owing to the lower and narrower amplitude ofthe pressure profile of the gas atmosphere in the cluster at z ∼ ∼ −
200 Myr sincethe injection of jets, both sources resemble FRII-type galaxies.However, an analysis of the distance of the radio emission peakfrom the jet core (e.g. Vardoulaki et al. 2019; Mingo et al. 2019)suggests that Run1 source is an FRI galaxy in a fraction of itsevolution. However, in the absence of a dominant relativisticcomponents in our simulated jets, our model cannot reproducehot spots of real radio sources.
We start by focusing on the dynamics of the ICM as a results ofradio galaxies launched at two different times. In Fig. 4 and Fig. 5, we show the evolution of the X-ray emis-sion (colors) overlayed with the projected distribution of radioemission at 50 MHzfrom our tracers (contour) for the two runs.In both cases, no observational cuts were applied in drawing theradio contours. Altogether, we show six snapshots to display themost salient features of AGN feedback: – i) the release of powerful shocks ahead of jets (with M ∼ . M ∼ . – ii) the formation of pairs of cavities, partially void of X-rayemitting gas in the region inflated by the two radio jets (sec-ond and third panel of each figure); – iii) in the Run2 case, the development of two ∼ −
250 kpclong jets, which remain straight at least ∼
400 Myr sincetheir ejection; in the Run1 case, jets are dissolved already ∼
200 Myr after their ejection; – iv) the subsequent evolution of jet-inflated radio lobes, whichexpand laterally (i.e. perpendicular to the jet axis) up to ∼ −
500 kpc; – v) the progressive mixing of the radio-emitting plasma withthe ICM, which covers an increasingly larger number of cellsin the innermost cluster regions. While in Run2 the distribu-tion of electrons has a large covering factor in the innermostcluster regions down to z = . Article number, page 8 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 5.
Evolution of the X-ray emission and of the radio emission (log contours) for the Run2 model, in which a radio galaxy is activated at z = .
0. The epochs in the panels are z = . z = . z = . z = . z = .
37 and z = .
19. No observational cut is applied to the maps,and the radio data have been smoothed to a × δ x = .
86 kpc). overall properties of the ICM is far from dramatic in these runs.The difference to Run0 models are limited to quite localised ex-cesses of gas heating due to the thermal feedback by AGN-drivenshock waves, and to slightly more pronounced excesses of mag-netic fields due to radio jets. However, while the differences areclear close to the epoch of jet launching ( z ≥ . ∼ −
10) remains visible in the radialprofile of the cluster even at later redshift, within ≤
250 kpc.Very similar trends are found when comparing the radial profilesof Run1 and Run0 (not shown).As a caveat, we should stress that our limited (and typicallydecreasing with radius) spatial resolution prevents the Reynoldsnumber from getting large enough for an efficient amplificationof magnetic fields seeded by radio galaxies, which can proceedonly when the spatial resolution of the simulation is enoughto resolve the Alfven length-scale (e.g. Donnert et al. 2018;Brunetti & Vazza 2020, and references therein). Hence, thesesimulations likely provide a lower limit to the ICM magneticfield strength after the jet activity has ceased, and future higherresoluton (and more expensive) simulations will be needed toassess the maximal possible amplification here. However, recentsimulations by Ehlert et al. (2020) support the fact that magneticfields injected by jets only affect the close proximity of AGN,while the cluster-wide magnetic field distribution is dominatedby large-scale turbulence.
We now focus on the properties of the ICM as viewed by thetracer particles as they move through the ICM. Fig. 7 gives theevolution of the distance travelled by tracers, in which we com- pared the median distance covered by tracer particles initiallyplaced in jets (or in corresponding cells, in Run0 models). Trac-ers ejected by radio jets always travel a larger distance comparedto their corresponding particles in Run0 models. This followsfrom the combined effect of outflows, as well as of large-scalemixing by turbulent motions, which are typically found to be sig-nificant at large radii (e.g. Angelinelli et al. 2020). The effect ismore significant in Run1, where the jet power is large enough topropel a fraction of the lobe to beyond the cluster virial radius,and in general is sufficient to spread the jets material over theentire cluster volume. In Run2, the jets are not powerful enoughto escape through the denser and hotter ICM, and the jet mate-rial settles to a similar distance as the Run0 case over time. Inboth jet simulations, being dominated by the fast bipolar trans-port, the average spatial separation of tracers is faster than aRichardson ( ∝ t / ) diffusion, which applies to the pair disper-sion statistics of passive tracers in Kolmogorov-type turbulence(e.g. Vazza et al. 2010). The latter trend is instead overall mea-sured initially in the pair dispersion statistics of tracers in Run0.Later on ( ≥ − ∝ t / ), approximatelyafter tracers have travelled a distance larger than the typical sizeof turbulent eddies. Interestingly, tracers ejected at z = Article number, page 9 of 21 &A proofs: manuscript no. main
Fig. 6.
Radial (mass weighted) profiles of gas density, gas temperature and proper magnetic filed strength for runs Run0 and Run2 at three differentepochs. the quantities for tracers found within the same radius from thecluster centre ( ≤
500 kpc) at all times. The changes of the am-bient medium that tracer particles are subjected to, are typicallyviolent and involve velocity fluctuations with amplitude several ∼
100 km / s, even during the later stages of lobe evolution, seee.g. Fig. 9 for Run2.Nevertheless, the activity of radio galaxies in both runs addsadditional fluctuations to the dynamical evolution of Lagrangiantracers, especially in the first hundreds of Myr since the AGNbursts. In particular, the magnetisation, temperature and vortic-ity of particles tracing jets are found to be significantly higher(up to an order of magnitude) up to ∼ −
400 Myr afterthe release of jets. Interestingly, the magnetic field strengths car-ried by injected tracers show the largest difference compared toRun0. In the case of Run1, the injected particles carry a signif-icantly higher magnetic field strength even ∼ ∼
400 Myr since their injection. Also in Run2, par-ticles released by jets have a significantly higher magnetic fieldthan their counterparts in Run0, up to ∼ ∼ − In our runs, relativistic electrons are always first seeded in theICM by radio jets, and are later subject to energy losses and/orre-acceleration processes, depending on the local gas conditionsthey encounter during their propagation.Most electrons released by radio galaxies suffer from signif-icant adiabatic losses as the radio lobes expand into the ICM,and hence their radio synchrotron power decreases quickly.However, a fraction of the electrons can undergo further(re)acceleration by shocks and turbulence.Fig. 10 shows the electron acceleration and losses in typicalconditions of a dynamical ICM. We give the median values ofthe timescales for energy losses τ loss , acceleration τ acc and advec-tion τ adv , using p = electrons as a reference. The timescales Article number, page 10 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 7.
Evolution of the median distance (in comoving Mpc) covered bytracers since their initial injection in our runs, as a function of time, forRun2 vs Run0 (top, injection at z = .
5) and for Run1 vs Run0 (bottom,injection at z = D ∝ t / , ∝ t and ∝ t / trends, to guide the eye. for shock acceleration ( τ DSA ) are much smaller because shocksonly cover a tiny fraction of the cluster volume. Hence, we showthe turbulent re-acceleration timescale, τ AS A (Sec. 2.3). Like-wise, the timescale for Coulomb losses ( τ c ) is much larger thanthe radiative timescale for IC and synchrotron losses ( τ rad ) for p = , hence only the second is displayed in the plot.We should point out some interesting trends: the radiativeloss time for electrons steadily increases with time, following thedecrease of the ∝ (1 + z ) term related to Inverse Compton losses(Eq. 5). The latter dominates the radiative loss terms consideringthat the typical magnetic field strength probed by most tracers incells gets smaller than the CMB-equivalent field strength most ofthe time. The contributions of adiabatic compression/expansionoften switche roles in their impact on the evolution of parti-cles, i.e. from enhancing to decreasing the particle energies, ontimescales of ≥ −
200 Myr. In both runs, the particles arefirst subjected to strong adiabatic losses, followed after 5 × yr by adiabatic gains due to strong compression by the externalICM. Later on, particles follow more erratic sequences of ex-pansion/compression as shocks waves cross the ICM. Moreover,particles are dispersed in a less dense ICM, particularly in Run1.The acceleration timescale related to the turbulence also signif-icantly varies during the evolution and depending on the localICM conditions; in general turbulent re-acceleration is less effi-cient to balance cooling losses for p = electrons in the Run1case, at least for half of the simulated evolution, because parti- cles diffuse to a larger cluster radius, where the relative impactof solenoidal turbulence is generally smaller (e.g. Vazza et al.2017), with the exception of a boost of turbulent accelerationin at least one peripheral sector of the group ( ∼ −
13 Gyr)following the last major merger experienced by the cluster (seeFig.1).Conversely, electrons ejected by the radio galaxy in Run2remains more confined in the innermost ICM, where solenoidalturbulence is dominant, and the median τ AS A is measured to besmaller (or of the same order) of the median synchrotron cool-ing time for p = electrons for most of the investigatedevolution. The complex evolution of the amplitude of differentloss/acceleration terms once more stresses the importance of fullnumerical simulations for a detailed simulation of the full evolu-tion of radio emitting electrons in a realistic environment.The panels in Fig. 11 shown an example of the spectral en-ergy evolution of 100 electrons injected by the radio galaxy inRun2, in which particles are subject to several possible combi-nations of re-acceleration events (or absence of thereof): onlylosses, losses and injection and re-acceleration by shocks, lossesand only shock re-acceleration, losses and only turbulent re-acceleration, or all mechanisms together, from top to bottom,respectively.In the presence of cooling and adiabatic losses alone, thisfamily of particles steadily accumulates at a momentum of p ∼ by the end of the run, consistent with the literature (Brunettiet al. 2001, 2007; Pinzke et al. 2013). The re-acceleration byshocks and turbulence instead largely increases the energy bud-get of particles in the p ≤ range, while the direct injectionof "fresh" electrons by shocks (here seeded by a M ∼ . ∼ r ≤
300 kpc with a solid line, 300 kpc < r ≤
600 kpcwith a dashed line and r >
600 kpc with a dotted line). Weobserve a similar evolution at all radii of a complex spectrum(albeit smoothed by the combination of the many more spectraconsidered here), with tails of p ≥ electrons. This is a conse-quence of the crossing of shock waves, and the steady build-up ofa fossil electron reservoir for p ≤ due to the combination ofshocks and turbulent re-acceleration. The p ∼ bump causedby turbulent re-acceleration is more evident in Run2, where theelectrons remain more confined in the cluster innermost regions.In those core regions, ∼ −
50% of the volume is typically filledby solenoidal turbulent motions (Vazza et al. 2017), while only ≤ .
1% of cells in clusters are shocked (e.g. Vazza et al. 2009).However, we should point out that our model cannot account forthe ∼ Article number, page 11 of 21 &A proofs: manuscript no. main
Fig. 8.
Evolution of the median magnetic field strength, temperature and absolute value of vorticity for tracers in our runs (the error bars showthe 16-84th percentiles of the distributions at each redshift) comparing Run2 vs Run0 (left panels, injection at z = .
5) and Run1 vs Run0 (right,injection at z = ≤ . Cassano et al. 2010), hence shorter acceleration times can be ex-pected for more massive objects.On the other hand, few electrons can exist with momentalarger than a few ∼ if they are only subject to radiativecooling and adiabatic expansion. As seen already above, re-accelerated electrons can reach a higher energy in Run2, whilein Run1 the bulk of the energy distribution of electrons spreadsout to larger cluster radii due to their earlier ejection from theirsource.Such differences in energy distribution obviously have a di-rect impact on the observable radio emission which depends onthe slope of available CR electrons. The distribution of the radio spectral index of the detectable emission in the two runs at dif-ferent epochs is shown in Fig. 14. For clarity, only the extremecases with all re-acceleration terms included, or only with losses,are shown. The radio-weighted distributions are generated onlyfor pixels above the putative detection threshold of a LOFARLBA observation. As expected, the radio spectral index distribu-tion of models where re-acceleration terms are included showsthe marked tendency of having flatter radio spectra, comparedto the pure cooling scenario, at each epoch and more increas-ingly in more evolved epochs. For most of their evolution after z ≤ . Article number, page 12 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 9.
Average velocity module for tracers in the jet at three different epochs ( z = . z = . z = . spectrum α ≥ .
25 (measured from 50 to 138 MHz), while theresidual emission in case only cooling and adiabatic losses areincluded is typically a factor ∆ α ∼ . ≤
140 MHz), and almost invisible at highfrequencies ( ≥ . ≥ σ detection level of LOFAR-LBA ∼ σ ≈ . / beam for a 30” ×
30” resolution beam,e.g. Bonafede et al. 2021), placing in all cases our source at a ref-erence distance of z = .
1. Figure 16 also shows the emission at1.4 GHz, limited to the last snapshot of the above panel, only forthe detectable emission above a realistic approximate sensitivityof a ∼ σ ≈ µ Jy / beam fora 5” ×
5” resolution beam, e.g. Rajpurohit et al. 2020). Clearly,only when all sources of re-acceleration from the ICM are in-cluded, the diffuse emission from fossil electrons injected byradio galaxies several billions years ago can become again de-tectable in the radio band. In the absence of turbulent and shockre-acceleration, the brightest emission patches of fossil electronscan hardly be seen as connected on ∼ Mpc scales at high fre-quencies,Hence, the late radio morphologies ( ≥
4. Caveats
This paper is only a first step towards a more realistic modelof the complex interplay between radio galaxies and their envi-ronment, and it is admittedly very limited in several importantaspects.First, our runs do not include radiative gas cooling, hencethere is no self-consistent coupling between the feeding of gasundergoing overcooling and the duty cycle of SMBH jets. Thiswould introduce a number of numerical and physical compli-cations in the modelling of the ICM, which are not entirelysolved in cosmological simulations (e.g. Dubois & Teyssier2008; Collins et al. 2010; Ruszkowski et al. 2011; Tremmel et al.2017) but are nevertheless crucial for any realistic modelling ofICM atmospheres (e.g. Voit et al. 2020). While this choice wasmotivated to best focus on the dynamical effect of the feedbackfrom radio galaxies (or absence thereof) on the amplification of
Article number, page 13 of 21 &A proofs: manuscript no. main
Fig. 10.
Evolution of the median timescales for acceleration or losses(see Sec. 2.3) experienced by electrons (here we consider γ = as areference) released by radio galaxies in Run1 (top) and Run2 (bottom).Timescales for loss terms and shown as negative, while timescales forgain terms are shown as positive (to ease the comparison with radiativelosses, the latter timescale is also mirrored in the upper panel). Thetimescales for shock acceleration ( τ DSA ) and for Coulomb losses ( τ C )are not visible as they are much smaller and much bigger, respectively,than the time range covered by vertical axis in the plots. magnetic fields in the ICM and on the evolution of relativis-tic electrons, removed by the additional effects of compressionfrom gas cooling, future developments will have to couple thelaunching of radio jets with the duty cycle of AGN feedback(e.g. Bourne & Sijacki 2020; Thomas et al. 2020). While ourchoice of having single and well defined episodes of jets releaseby radio galaxies is convenient to study the dynamics of ejectedelectrons in a clean way, recent results from radio surveys haveinstead suggested that at least the most massive radio galaxiescan be "on" for most of their lifetime (e.g. Sabater et al. 2019),unlike the discontinuous activity scenario explored in this work.Only through a more consistent coupling of cooling/feedback inour simulations we can approach the complexity of real duty cy-cles of radio galaxies. The extension of this study with clustersin a different dynamical state will also allow us to explore differ-ent combinations of shocks and turbulent re-acceleration events,naturally triggered by different merger configurations.Secondly, our simulations cannot reproduce the actual inter-action between the relativistic (either electron or proton domi-nated) content of jets and the thermal ICM, as our jets are filledentirely with magnetic fields and hot thermal gas, as commonlydone by simulations of this kind (e.g. Xu et al. 2009; Gaspari Fig. 11.
Evolution (from z = . z = . ∼ Fig. 12.
Evolution of the spectral energy distributions of electrons emitted by the radio galaxy in Run2, for three different acceleration/coolingmodels (colors) and for three different radial ranges: r ≤
300 kpc (solid), 300 kpc < r ≤
600 kpc (dashed) and r >
600 kpc (dotted). et al. 2011; Teyssier et al. 2011). In order to directly test theimpact of a significant relativistic energy content in protons orelectrons, future development will have to resort to two-fluidmodelling of the ICM (e.g. Ehlert et al. 2018; Yang et al. 2019),similar (e.g. Vazza et al. 2013).The use of passive tracer particles to track the propagation ofrelativistic electrons introduces an element of discreteness in ourmodelling, in the sense that compared to a fluid approach (i.e. viatwo-fluid modelling, e.g. Girichidis et al. 2020, Ogrodnik et al.2020) it cannot fully map the spatial diffusion of jet ejecta overtime. In this respect, the small-scale morphology of our radiojets may change if the spatial diffusion of relativistic electrons istaken into account, even if the global dynamics of lobes is notexpected to change significantly.Our Fokker-Planck model for particles (Sec. 2.3) is only con-cerned with specific choices for the acceleration of relativisticelectrons via shock waves. It also includes a simplified treat-ment of the turbulent re-acceleration, in which only a systematicacceleration effect on the global distribution can be captured.Significant effects on the highest energy tail of the electron dis-tributions can be expected for different choices in accelerationscenarios (e.g. Wittor et al. 2020), and will be explored in futurework with more sophisticated numerical schemes (e.g. Donnert& Brunetti 2014).Finally, our simulations only evolved a P M = R M / R e = ν/η ≈
5. Conclusions
In this paper, we simulated the evolution of relativistic electronsinjected in the ICM by active radio galaxies. Using a Fokker-Planck method on Lagrangian tracer particles, we tracked thecooling and re-acceleration of electrons as a function of time.We find that relativistic electrons from radio galaxies efficientlyspread across the intracluster volume, as a result of cluster-wideturbulence and the momentum of powerful jets. The magneticoutput of jets dominates close to the radio galaxies ( ≤
500 kpc),yet the large-scale distribution of magnetic fields seems to bedominated by the growth of the galaxy cluster. Hence, the in-put from jets is subdominant with respect to the overall clustermagnetic field - with the exception of the cavities inflated byjets, which preserve a higher and distinct magnetic field up to ∼ Article number, page 15 of 21 &A proofs: manuscript no. main
Fig. 13.
Evolution of the spectral energy distributions of electrons emitted by the radio galaxy in Run1, for three different acceleration/coolingmodels (colors) and for three different radial ranges: r ≤
300 kpc (solid), 300 kpc < r ≤
600 kpc (dashed) and r >
600 kpc (dotted). on the details of (re)acceleration events caused by shock wavesand turbulence.Spatially and spectrally resolved radio observations (espe-cially in the low-frequency regime explored by LOFAR, MWA,and in the future by the SKA-LOW) will be able to retrace thesequence of radio-mode feedback and probe a number of plasmaprocesses.In conclusion, electrons and magnetic fields seeded by radiogalaxies in galaxy clusters are important for modelling the non-thermal emission on the largest scales in the Universe. This hasimportant applications to the study of feedback physics and ofdiffuse radio emissions in the intracluster medium.
Acknowledgments
References
Abdulla, Z., Carlstrom, J. E., Mantz, A. B., et al. 2019, ApJ, 871, 195Ackermann, M., Ajello, M., Albert, A., et al. 2014, ApJ, 787, 18Angelinelli, M., Vazza, F., Giocoli, C., et al. 2020, MNRAS, 495, 864Banfi, S., Vazza, F., & Wittor, D. 2020, arXiv e-prints, arXiv:2006.10063Beresnyak, A. & Miniati, F. 2016, ApJ, 817, 127Berezinsky, V. S., Blasi, P., & Ptuskin, V. S. 1997, ApJ, 487, 529Bîrzan, L., Rafferty, D. A., Brüggen, M., et al. 2020, MNRAS, 496, 2613Bîrzan, L., Rafferty, D. A., Nulsen, P. E. J., et al. 2012, MNRAS, 427, 3468Bodo, G., Rossi, P., Massaglia, S., et al. 1998, A&A, 333, 1117Bonafede, A., Brunetti, G., Vazza, F., et al. 2021, ApJ, 907, 32Bonafede, A., Intema, H. T., Brüggen, M., et al. 2014, ApJ, 785, 1Booth, C. M. & Schaye, J. 2009, MNRAS, 398, 53Borse, N., Acharya, S., Vaidya, B., et al. 2020, arXiv e-prints, arXiv:2009.13540Botteon, A., Brunetti, G., Ryu, D., & Roh, S. 2020a, A&A, 634, A64Botteon, A., Brunetti, G., van Weeren, R. J., et al. 2020b, ApJ, 897, 93Bourne, M. A. & Sijacki, D. 2020, arXiv e-prints, arXiv:2008.12784Brighenti, F. & Mathews, W. G. 2000, ApJ, 535, 650Brüggen, M. & Kaiser, C. R. 2002, Nature, 418, 301Brüggen, M., Reiprich, T. H., Bulbul, E., et al. 2020, arXiv e-prints,arXiv:2012.08775Brunetti, G. & Jones, T. W. 2014, International Journal of Modern Physics D,23, 1430007
Article number, page 16 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 14.
Evolution of the (radio power weighted) distribution of radio emission spectra between 50 and 138 MHz, for four epochs in the Run2 andonly including radio detectable pixels at 50MHz (assuming a LOFAR LBA observation). Only the cooling and cooling+shocks+turbulence casesare shown for clarity.
Brunetti, G. & Lazarian, A. 2011, MNRAS, 410, 127Brunetti, G. & Lazarian, A. 2016, MNRAS, 458, 2584Brunetti, G., Setti, G., Feretti, L., & Giovannini, G. 2001, MNRAS, 320, 365Brunetti, G. & Vazza, F. 2020, Phys. Rev. Lett., 124, 051101Brunetti, G., Venturi, T., Dallacasa, D., et al. 2007, ApJ, 670, L5Bryan, G. L., Norman, M. L., O’Shea, B. W., et al. 2014, ApJS, 211, 19Bykov, A. M., Vazza, F., Kropotina, J. A., Levenfish, K. P., & Paerels, F. B. S.2019, Space Sci. Rev., 215, 14Candelaresi, S. & Del Sordo, F. 2020, ApJ, 896, 86Capetti, A., Baldi, R. D., Brienza, M., Morganti, R., & Giovannini, G. 2019,A&A, 631, A176Cassano, R. & Brunetti, G. 2005, MNRAS, 357, 1313Cassano, R., Ettori, S., Giacintucci, S., et al. 2010, ApJ, 721, L82Chang, J. S. & Cooper, G. 1970, Journal of Computational Physics, 6, 1Churazov, E., Brüggen, M., Kaiser, C. R., Böhringer, H., & Forman, W. 2001,ApJ, 554, 261Collins, D. C., Xu, H., Norman, M. L., Li, H., & Li, S. 2010, ApJS, 186, 308Croston, J. H., Ineson, J., & Hardcastle, M. J. 2018, MNRAS, 476, 1614de Gasperin, F., Intema, H. T., Shimwell, T. W., et al. 2017, Science Advances,3, e1701634de Gasperin, F., Orrú, E., Murgia, M., et al. 2012, A&A, 547, A56Domínguez-Fernández, P., Vazza, F., Brüggen, M., & Brunetti, G. 2019, MN-RAS, 486, 623Donnert, J. & Brunetti, G. 2014, MNRAS, 443, 3564Donnert, J., Vazza, F., Brüggen, M., & ZuHone, J. 2018, ArXiv e-prints[ arXiv:1810.09783 ]Dubois, Y., Devriendt, J., Slyz, A., & Teyssier, R. 2010, MNRAS, 409, 985Dubois, Y. & Teyssier, R. 2008, A&A, 482, L13Ehlert, K., Weinberger, R., Pfrommer, C., Pakmor, R., & Springel, V. 2018, MN-RAS, 481, 2878Ehlert, K., Weinberger, R., Pfrommer, C., & Springel, V. 2020, arXiv e-prints,arXiv:2011.13964Fabjan, D., Borgani, S., Tornatore, L., et al. 2010, MNRAS, 401, 1670Fanaroff, B. L. & Riley, J. M. 1974, MNRAS, 167, 31P Furlanetto, S. R. & Loeb, A. 2001, ApJ, 556, 619Gabici, S. & Blasi, P. 2003, ApJ, 583, 695Gaspari, M., Brighenti, F., D’Ercole, A., & Melioli, C. 2011, MNRAS, 415, 1549Gaspari, M., Ruszkowski, M., & Sharma, P. 2012, ApJ, 746, 94Gendron-Marsolais, M.-L., Hull, C. L. H., Perley, R., et al. 2021, arXiv e-prints,arXiv:2101.05305Girichidis, P., Pfrommer, C., Hanasz, M., & Naab, T. 2020, MNRAS, 491, 993Guo, X., Sironi, L., & Narayan, R. 2014a, ApJ, 794, 153Guo, X., Sironi, L., & Narayan, R. 2014b, ApJ, 797, 47Hahn, O. & Abel, T. 2011, MNRAS, 415, 2101Hardcastle, M. J. & Croston, J. H. 2020, New A Rev., 88, 101539Hardcastle, M. J. & Krause, M. G. H. 2014, MNRAS, 443, 1482Heinz, S., Brüggen, M., Young, A., & Levesque, E. 2006, MNRAS, 373, L65Iapichino, L. & Niemeyer, J. C. 2008, MNRAS, 388, 1089Ignesti, A., Shimwell, T., Brunetti, G., et al. 2020, A&A, 643, A172Jaffe, W. J. & Perola, G. C. 1973, A&A, 26, 423Jones, T. W., Nolting, C., O’Neill, B. J., & Mendygral, P. J. 2017, Physics ofPlasmas, 24, 041402Kaiser, C. R. & Alexander, P. 1997, MNRAS, 286, 215Kang, H. 2018, Journal of Korean Astronomical Society, 51, 185Kang, H. & Jones, T. W. 2007, Astroparticle Physics, 28, 232Kang, H. & Ryu, D. 2011, ApJ, 734, 18Kang, H., Ryu, D., & Jones, T. W. 2012, ApJ, 756, 97Kardashev, N. S. 1962, Soviet Ast., 6, 317Kim, J.-h., Wise, J. H., Alvarez, M. A., & Abel, T. 2011, ApJ, 738, 54Kronberg, P. P., Lesch, H., & Hopp, U. 1999, ApJ, 511, 56Lau, E. T., Kravtsov, A. V., & Nagai, D. 2009, ApJ, 705, 1129Locatelli, N., Vazza, F., & Domínguez-Fernández, P. 2018, Galaxies, 6, 128Mandal, S., Intema, H. T., van Weeren, R. J., et al. 2020, A&A, 634, A4Markevitch, M., Govoni, F., Brunetti, G., & Jerius, D. 2005, ApJ, 627, 733Massaglia, S., Bodo, G., Rossi, P., Capetti, S., & Mignone, A. 2016, A&A, 596,A12Massaglia, S., Bodo, G., Rossi, P., Capetti, S., & Mignone, A. 2019, A&A, 621,A132
Article number, page 17 of 21 &A proofs: manuscript no. main
Fig. 15.
Evolution of the radio emission for Run2, for three epochs and the three re-acceleration models investigated in the paper. Only the emissionabove the putative 3 σ detectability using LOFAR-HBA (140 MHz) is shown. Fig. 16.
Radio emission for Run2 at 1.4 GHz, for the three re-acceleration models investigated in the paper. Only the emission above the putative3 σ detectability using JVLA is shown. Mathews, W. G. & Brighenti, F. 2007, ApJ, 660, 1137McCarthy, I. G., Schaye, J., Ponman, T. J., et al. 2010, MNRAS, 406, 822McCourt, M., Parrish, I. J., Sharma, P., & Quataert, E. 2011, MNRAS, 413, 1295Mendygral, P. J., Jones, T. W., & Dolag, K. 2012, ApJ, 750, 166Mignone, A., Massaglia, S., & Bodo, G. 2005, Space Sci. Rev., 121, 21Mignone, A., Rossi, P., Bodo, G., Ferrari, A., & Massaglia, S. 2010, MNRAS,402, 7Mingo, B., Croston, J. H., Hardcastle, M. J., et al. 2019, MNRAS, 488, 2701Miniati, F. 2014, ApJ, 782, 21Nolting, C., Jones, T. W., O’Neill, B. J., & Mendygral, P. J. 2019a, ApJ, 876, 154Nolting, C., Jones, T. W., O’Neill, B. J., & Mendygral, P. J. 2019b, ApJ, 885, 80Norman, M. L. & Bryan, G. L. 1999, in Lecture Notes in Physics, Berlin SpringerVerlag, Vol. 530, The Radio Galaxy Messier 87, ed. H.-J. Röser & K. Meisen-heimer, 106–+ Norman, M. L., Winkler, K. H. A., Smarr, L., & Smith, M. D. 1982, A&A, 113,285Ogrodnik, M., Hanasz, M., & Wólta´nski, D. 2020, arXiv e-prints,arXiv:2009.06941O’Neill, S. M. & Jones, T. W. 2010, ApJ, 710, 180Owen, F. N., Rudnick, L., Eilek, J., et al. 2014, ApJ, 794, 24Perucho, M. & Martí, J. M. 2007, MNRAS, 382, 526Pinzke, A., Oh, S. P., & Pfrommer, C. 2013, MNRAS, 435, 1061Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2016, A&A, 594, A13Puchwein, E., Sijacki, D., & Springel, V. 2008, ApJ, 687, L53Rajpurohit, K., Hoeft, M., Vazza, F., et al. 2020, A&A, 636, A30Ramatsoku, M., Murgia, M., Vacca, V., et al. 2020, A&A, 636, L1Rasia, E., Borgani, S., Murante, G., et al. 2015, ApJ, 813, L17
Article number, page 18 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. 17.
Radio spectral index maps for the 140-50 MHz range, for three epochs and two lines of sight, for the Run2 simulation. Only regionsdetectable with LOFAR-HBA (as in Fig.15) are shown.
Ricarte, A., Tremmel, M., Natarajan, P., & Quinn, T. 2019, MNRAS, 489, 802Rincon, F., Califano, F., Schekochihin, A. A., & Valentini, F. 2016, Proceedingsof the National Academy of Science, 113, 3950Riquelme, M. A. & Spitkovsky, A. 2011, ApJ, 733, 63Ruszkowski, M., Enßlin, T. A., Brüggen, M., Heinz, S., & Pfrommer, C. 2007,MNRAS, 378, 662Ruszkowski, M., Lee, D., Brüggen, M., Parrish, I., & Oh, S. P. 2011, ApJ, 740,81Ryu, D., Kang, H., Hallman, E., & Jones, T. W. 2003, ApJ, 593, 599Sabater, J., Best, P. N., Hardcastle, M. J., et al. 2019, A&A, 622, A17Sarazin, C. L. 1999, ApJ, 520, 529Schekochihin, A. A., Cowley, S. C., Taylor, S. F., Maron, J. L., & McWilliams,J. C. 2004, ApJ, 612, 276Sijacki, D., Springel, V., Di Matteo, T., & Hernquist, L. 2007, MNRAS, 380, 877Stuardi, C., Bonafede, A., Wittor, D., et al. 2019, MNRAS, 489, 3905Teyssier, R., Moore, B., Martizzi, D., Dubois, Y., & Mayer, L. 2011, MNRAS,414, 195Thomas, N., Dave, R., Jarvis, M. J., & Angles-Alcazar, D. 2020, arXiv e-prints,arXiv:2010.11225Tremblay, G. R., Oonk, J. B. R., Combes, F., et al. 2016, Nature, 534, 218Tremmel, M., Karcher, M., Governato, F., et al. 2017, MNRAS, 470, 1121Turner, R. J. & Shabala, S. S. 2015, ApJ, 806, 59Uchiyama, Y., Aharonian, F. A., Tanaka, T., Takahashi, T., & Maeda, Y. 2007,Nature, 449, 576Vachaspati, T. 2020, arXiv e-prints, arXiv:2010.10525van Weeren, R. J., Andrade-Santos, F., Dawson, W. A., et al. 2017, Nature As-tronomy, 1, 0005van Weeren, R. J., de Gasperin, F., Akamatsu, H., et al. 2019, Space Sci. Rev.,215, 16Vardoulaki, E., Jiménez Andrade, E. F., Karim, A., et al. 2019, A&A, 627, A142Vazza, F., Brueggen, M., Gheller, C., et al. 2017, Classical and Quantum GravityVazza, F. & Brüggen, M. 2014, MNRAS, 437, 2291Vazza, F., Brüggen, M., & Gheller, C. 2013, MNRAS, 428, 2366Vazza, F., Brunetti, G., Brüggen, M., & Bonafede, A. 2018, MNRAS, 474, 1672Vazza, F., Brunetti, G., & Gheller, C. 2009, MNRAS, 395, 1333Vazza, F., Brunetti, G., Gheller, C., Brunino, R., & Brüggen, M. 2011, A&A,529, A17+Vazza, F., Eckert, D., Brüggen, M., & Huber, B. 2015, MNRAS, 451, 2198 Vazza, F., Gheller, C., & Brunetti, G. 2010, A&A, 513, A32+Vazza, F., Jones, T. W., Brüggen, M., et al. 2017, MNRAS, 464, 210Voit, G. M., Bryan, G. L., Prasad, D., et al. 2020, ApJ, 899, 70Völk, H. J. & Atoyan, A. M. 1999, Astroparticle Physics, 11, 73Wilber, A., Brüggen, M., Bonafede, A., et al. 2018, MNRAS, 473, 3536Wittor, D. 2017, PhD thesis, Universität Hamburg, Von-Melle-Park 3, 20146HamburgWittor, D. & Gaspari, M. 2020, MNRAS, 498, 4983Wittor, D., Jones, T., Vazza, F., & Brüggen, M. 2017a, MNRAS, 471, 3212Wittor, D., Vazza, F., & Brüggen, M. 2016, Galaxies, 4, 71Wittor, D., Vazza, F., & Brüggen, M. 2017b, MNRAS, 464, 4448Wittor, D., Vazza, F., Ryu, D., & Kang, H. 2020, MNRAS, 495, L112Xie, C., van Weeren, R. J., Lovisari, L., et al. 2020, arXiv e-prints,arXiv:2001.04725Xu, H., Li, H., Collins, D. C., Li, S., & Norman, M. L. 2009, ApJ, 698, L14Xu, H., Li, H., Collins, D. C., Li, S., & Norman, M. L. 2011, ApJ, 739, 77Xu, R., Spitkovsky, A., & Caprioli, D. 2020, ApJ, 897, L41Xu, S. & Zhang, B. 2017, ApJ, 846, L28Yang, H. Y. K., Gaspari, M., & Marlow, C. 2019, ApJ, 871, 6Yang, H.-Y. K., Sutter, P. M., & Ricker, P. M. 2012, MNRAS, 427, 1614ZuHone, J. A., Markevitch, M., & Johnson, R. E. 2010, ApJ, 717, 908
Article number, page 19 of 21 &A proofs: manuscript no. main
Appendix A: Visual impression of the two runs.
Figure A.1 and Figure A.2 give the visual impression ofthe refined evolution of Run1 and Run2, using RGB colorcoding of different information: radio emission (pink col-ors), X-ray emission (green) and gas temperature (blue),in all cases applying logarithmic stretching and self-renormalisation of all emissions in order to enhances thevisibility of structures. The full movies of the two simu-lations can be found at https://vimeo.com/490397871 andhttps://vimeo.com/490399056 .
Article number, page 20 of 21. Vazza, D. Wittor, G. Brunetti & M. Brüggen: Simulated Radio Galaxies in ENZO
Fig. A.1.
Rendering of the evolution of radio emission (pink colors), X-ray emission (green) and gas temperature (blue) for roughly equally spacedtimesteps from z = .
49 to z = . Fig. A.2.
Rendering of the evolution of radio emission (pink colors), X-ray emission (green) and gas temperature (blue) for roughly equally spacedtimesteps from z = .
99 to z = ..