Driving of ICM Turbulence by AGN

# Inefficient Driving of Bulk Turbulence by Active Galactic Nuclei in a Hydrodynamic Model of the Intracluster Medium

Christopher S. Reynolds11affiliation: Department of Astronomy, University of Maryland, College Park, MD 20742-2421, USA; chris@astro.umd.edu 22affiliation: Joint Space-Science Institute (JSI), College Park, MD 20742-2421, USA , Steven A. Balbus33affiliation: Department of Physics, University of Oxford, Oxford OX1 3NP, UK , and Alexander A. Schekochihin44affiliation: Rudolf Peierls Centre for Theoretical Physics, University of Oxford, Oxford OX1 3NP, UK 55affiliation: Merton College, Oxford OX1 4JD, UK
###### Abstract

Central jetted active galactic nuclei (AGN) appear to heat the core regions of the intracluster medium (ICM) in cooling-core galaxy clusters and groups, thereby preventing a cooling catastrophe. However, the physical mechanism(s) by which the directed flow of kinetic energy is thermalized throughout the ICM core remains unclear. We examine one widely discussed mechanism whereby the AGN induces subsonic turbulence in the ambient medium, the dissipation of which provides the ICM heat source. Through controlled inviscid 3-d hydrodynamic simulations, we verify that explosive AGN-like events can launch gravity waves (-modes) into the ambient ICM which in turn decay to volume-filling turbulence. In our model, however, this process is found to be inefficient, with less than 1% of the energy injected by the AGN activity actually ending up in the turbulence of the ambient ICM. This efficiency is an order of magnitude or more too small to explain the observations of AGN-feedback in galaxy clusters and groups with short central cooling times. Atmospheres in which the -modes are strongly trapped/confined have an even lower efficiency since, in these models, excitation of turbulence relies on the -modes’ ability to escape from the center of the cluster into the bulk ICM. Our results suggest that, if AGN-induced turbulence is indeed the mechanism by which the AGN heats the ICM core, its driving may rely on physics beyond that captured in our ideal hydrodynamic model.

galaxies: clusters: intracluster medium — hydrodynamics — turbulence
slugcomment:

## 1. Introduction

It is now widely accepted that active galactic nuclei (AGN) regulate the growth of the most massive galaxies (see review by Fabian, 2012, and references therein). Cooling-core clusters of galaxies provide a remarkably accessible view of these AGN feedback processes at work. Most of the baryonic matter in galaxy clusters resides in an atmosphere of hot (temperatures ) and tenuous (particle-number densities ) intracluster medium (ICM). Without any energy input, the observed atmospheres would radiatively cool on a timescale of or less leading to a cooling catastrophe, levels of star formation far greater than those observed, and central galaxies that are much more massive than found in nature (Fabian, 1994; McNamara & Nulsen, 2007). Instead, real clusters appear to be in a state of approximate thermal equilibrium, with some heat source appearing to balance the radiative cooling (Peterson & Fabian, 2006). Given the failings of models based on thermal conduction (Binney & Cowie, 1981; Ruszkowski & Begelman, 2002; Kim & Narayan, 2003; Conroy & Ostriker, 2008; Bogdanović et al., 2009), energy injection by an AGN hosted by the central galaxy appears to be the most viable source of this heat (Churazov et al., 2000; Reynolds, Heinz & Begelman, 2002; Churazov et al., 2002; Omma et al., 2004).

Indeed, the vast majority of cool-core clusters have central galaxies that host radio-loud AGN (Burns, 1990). Furthermore, when sufficiently deep X-ray images are available, we invariably see evidence of interaction between the jets of these AGN and the ICM. The most common signature of interaction are bi-polar ICM cavities that have been blown by jet-inflated cocoons (Fabian et al., 2000; Churazov et al., 2001; Heinz et al., 2002; Bîrzan et al., 2004). In a hydrodynamic model, blowing well-defined low-density cavities before they are shredded by Rayleigh-Taylor instabilities requires supersonic initial expansion and hence we expect the presence of shocks within the ICM. While strong shocks are never observed (implying that the strong shock phase may be short lived), there is evidence for weak shocks in a number of sources (e.g., Jones et al., 2002; Million et al., 2010; Randall et al., 2015), supporting the basic picture of supersonic jet-driven cavity expansion.

The basic energetics of the AGN heating of ICM cores also make sense. In the Perseus cluster, for example, the X-ray luminosity emerging from within the cooling radius (usually defined as the radius within which the cooling time falls below 3 Gyr) is . To create the inner, most well-defined pair of cavities requires in total approximately , assuming adiabatic inflation. However, they are bounded by weak shocks which add an additional to the energy budget (Graham, Fabian & Sanders, 2008). Noting that these power estimates suffer from order-unity uncertainty owing to projection effects and uncertainties in the timescale of cavity inflation, we see that the AGN can clearly inject sufficient energy to heat the ICM core, provided that the process is reasonably efficient. More generally, studies of samples of groups and clusters show that cavities are ubiquitous in systems with short central cooling times, and that between of energy per cavity (where is the pressure of the surrounding ICM and is the cavity volume) needs to be thermalized in order to offset the radiative cooling (Bîrzan et al., 2004; Panagoulia et al., 2014).

However, going beyond these global energy arguments, the actual mechanism by which the AGN heating of the ICM occurs is not well understood. Energy injection could be just that—direct injection of hot bubbles that mix with the cluster gas (Heinz & Churazov, 2005; Hillel & Soker, 2015). However, at least in the most powerful observed outbursts in intermediate temperature clusters, any such mixing must be slow (i.e., comparable to the cooling timescale) or else the cooling core would have been destroyed by the injected enthalpy (Pfrommer, Chang & Broderick, 2012; McNamara & Nulsen, 2012). Potentially of equal or greater importance is heating of the ICM by cosmic-ray (CR) protons as they diffuse/stream through the ICM (e.g., Loewenstein, Zweibel & Begelman, 1991; Guo & Oh, 2008; Pfrommer, 2013). The most likely source of CRs relevant for heating the cluster core is the jet-inflated cocoon/bubble itself and, hence, this is closely related to the bubble mixing issue already discussed.

Alternatively, the AGN energy injection can take the form of mechanical agitation that propagates into the ICM atmosphere initially in the form of waves, with the wave energy ultimately being thermalised (Churazov et al., 2004). In this case, the heating is a secondary process. There are multiple paths by which the thermalization might occur, depending upon the nature of the perturbation itself. Sound waves propagate rapidly through the core and tend to steepen into weak shocks. Buoyancy, or “-modes” are much less compressive and propagate much less rapidly and extensively (indeed, as discussed in Section 2, they can be trapped within a finite radius). These waves might be themalized directly if conduction or viscosity are sufficiently large (Ruszkowski, Brüggen & Begelman, 2004; Fabian et al., 2005; Heinz & Churazov, 2005), or indirectly on small scales via a turbulent cascade (Fujita, 2005; Schekochihin & Cowley, 2006; Enßlin & Vogt, 2006; Rebusco et al., 2006; Kunz et al., 2011; Zhuravleva et al., 2014a). The turbulent path could be particularly important when direct wave dissipation is inefficient. None of these pathways to thermalization of mechanical energy has been systematically explored to see whether they are dynamically viable.

It is this important gap that our current work seeks to fill —or begin to fill—by isolating and studying the process of turbulent dissipation of the g-modes generated by stirring in a cluster-like environment. The turbulent heating picture recently received a major boost from a detailed analysis of X-ray surface brightness (SB) fluctuations in the Perseus and Virgo clusters as seen in deep Chandra images. Interpreting the SB fluctuations as being due to turbulence111Of course, if there are any additional SB fluctuations arising from clumping/multiphase structure within the ICM, this analysis method would over-estimate the turbulent heating., Zhuravleva et al. (2014a) determine the velocity power spectrum of the turbulence and proceed to estimate the volumetric turbulent dissipation as a function of radius. They show that, within each radial bin that they probe, the derived turbulent dissipation approximately matches the observed radiative cooling. The turbulence itself is envisaged to be due to the decay of -modes driven into the bulk of the ICM through the action of AGN (during the inflation and subsequent buoyant rise of the jet-blown cavities; Omma et al., 2004). Of course, the many other pathways by which wave thermalization might occur, including the role of magnetic fields, are a priori no less important, but lie outside the scope of the current work.

In this paper, we investigate the theoretical foundations of the AGN-driven turbulent heating picture. We construct a toy-model atmosphere in which we then explore the driving of the hydrodynamic -modes by explosive events (representing the AGN outbursts) and the subsequent decay of those modes into turbulence. Our principal finding is that the AGN-like energy injection can indeed drive volume-filling turbulence throughout the ICM atmosphere but, within our particular model system, the transfer of energy from the AGN into the turbulent cascade is rather inefficient. In our fiducial case, less than 1% of the injected energy ends up as the kinetic energy of turbulence in the ambient ICM. The majority of the injected energy is put into compressible modes (weak shocks and sound waves) that radiate out of the top boundary of our atmosphere, or are rapidly thermalized in the strongly turbulent cavity-debris zone. Furthermore, we find that atmospheres in which the -modes are strongly trapped are even more inefficient because, in our models, driving turbulence in the bulk of the ICM is contingent on the -modes being able to escape outwards from the center of the cluster. If the AGN-driven-turbulence picture gains strong support from future observational studies (for example, from direct X-ray spectroscopic measurements of large turbulent velocities by Astro-H), our work suggests that physics beyond ideal hydrodynamics is required to understand the driving of the turbulence.

This paper is organized as follows. In Section 2, we review some necessary theoretical concepts including the properties of -modes in ICM atmospheres. Our toy model is described in Section 3. Section 4 describes the results for our fiducial parameter set with only weak -mode trapping; the strong-trapping cases are treated in Section 5. The implications of our findings, as well as the limitations and shortcomings of our approach, are discussed in Section 6. Finally, Section 7 draws together our conclusions.

## 2. Theoretical Preliminaries

Before launching into a discussion of our model atmospheres and their dynamics, we first review in brief some relevant theoretical background.

We are concerned with the adiabatic hydrodynamic disturbances of an atmosphere about some hydrostatic equilibrium. In linear theory, there are two families of modes corresponding to whether pressure or gravity acts as the restoring force. The pressure modes are the familiar sound waves which propagate at the sound speed, where is the thermal pressure, the density, and the ratio of specific heat capacities. The pressure modes radiate out of any given region of interest on the sound crossing time.

The gravity modes, or -modes, are of central importance to our discussion. To lowest order in the (small) Mach number of the velocity perturbations, they are incompressible, , where is the velocity field. This implies that the excitation of -modes leads to the generation of vorticity, . The evolution of the vorticity in the case of an ideal inviscid fluid subject to only conservative external forces is

 DωDt=(ω⋅∇)v−ω(∇⋅%$v$)+1ρ2∇ρ×∇P, (1)

where is the convective derivative. To linear order in the small velocities this gives

 ∂ω∂t=1ρ2∇ρ×∇P. (2)

Hence, vorticity production is associated with departures between surfaces of constant density and those of constant pressure. In the case of linear -modes in an atmosphere, these departures are a consequence of the fact that surfaces have first-order (in the Mach number of the velocity perturbations) corrugations relative to equipotentials whereas surfaces only have second-order corrugations. Given that the pressure gradient is then essentially in the vertical direction (as defined by the local gravitational field), the term and hence the generated vorticity will lie principally in the horizontal plane. We shall use this fact later in the paper to aid in the identification of -modes.

Gravity modes have the important property that they can be trapped. To see this, we note that the dispersion relation of a local -mode with a space-time dependence is

 Ω2=k2hk2N2. (3)

Here, , is the magnitude of the wavenumber projected onto the horizontal plane, is the Brunt-Väisälä frequency given by

 N2=g(z)γ∂∂zln(Pρ−γ), (4)

and the gravitational field is defined to be ( is the unit vector in the -direction). Since , -modes have a maximum possible frequency of , at which point . If the Brunt-Väisälä frequency is a decreasing function of height , -modes of a given frequency will be trapped/confined below the height at which . We shall see that -mode trapping may have a significant effect on the ability of localized AGN driving to generate turbulence in the ambient medium.

## 3. The Model and Computational Setup

Our goal is to examine the basic physical ingredients of the AGN-driven turbulent heating scenario rather than construct a realistic global model of a turbulently heated, radiatively-cooling ICM atmosphere. Thus, we abstract the problem to one of a plane-parallel atmosphere that, at the initial time, is in hydrostatic equilibrium with a time-invariant gravitational field . Then, working within the framework of ideal-gas hydrodynamics, we mimic AGN jet driving by detonating thermal bombs close to the base of the atmosphere at random time intervals. At some given point in time, we stop the driving and study the decay of the motions (turbulent or otherwise) within the atmosphere.

We adopt a Cartesian coordinate system . Neglecting self-gravity of the ICM plasma, we assume a (fixed) gravitational field of the form

 g=−g01+z/z0^z. (5)

At time , the atmosphere is isothermal and in hydrostatic equilibrium,

 ∇P=−ρg, (6)

where the pressure and the density are related by the (initially uniform) isothermal sound speed in the usual manner, . Thus, with our chosen form for the gravitational field, the initial density profile satisfying equation (6) is

 ρ=ρ0(1+zz0)−z0/h, (7)

where is the pressure scale height at the base () of the atmosphere. Equation (7) tends to the usual exponential atmosphere as . With , where is the internal energy, we find that the Brunt-Väisälä frequency at the initial time is

 N=g0(γ−1)1/2ciso,0γ1/2(1+zz0)−1. (8)

Hence, the Brunt-Väisälä  frequency decreases more rapidly, leading to stronger trapping of -modes, if we consider smaller values of .

We use the PLUTO code (Mignone et al., 2007) to model the dynamics of this atmosphere subject to AGN driving. PLUTO evolves the equations of compressible ideal hydrodynamics in conservative form,

 ∂ρ∂t+∇⋅(ρv) = 0, (9) ∂∂t(ρv)+∇⋅(ρvv+PI) = −ρ∇Φ, (10) ∂∂t(E+ρΦ)+∇⋅[(E+P+ρΦ)v] = 0,

where is the fluid velocity, is the unit rank-two tensor, is the total energy density of the fluid,

 E=u+12ρ|v|2, (12)

and is the gravitational potential, which for our problem is, from equation (6),

 Φ=g0z0ln(1+zz0). (13)

Note that we evolve the total energy equation in conservative form — apart from losses at any open boundaries and energy that is explicitly injected into the system, energy is conserved to machine accuracy and hence the dissipation of turbulence by grid viscosity will be captured as heating.

Driving occurs in the form of thermal bombs close to the base of the atmosphere that occur randomly in time according to Poisson statistics with a mean recurrence time . When such an AGN driving event is triggered, we instantaneously increase the thermal pressure in a small sphere of radius centered at by an amount . We stop the driving and allow the atmospheric dynamics to decay after . For all simulations in this paper, we use

 tinj = 10h/ciso,0, (14) tstop = 200h/ciso,0, (15) rbub = 0.1h, (16) zbub = 0.2h, (17) ΔPbub = 5ρ0c2iso,0. (18)

The lateral location ( and ) of each detonation is chosen from a uniform random distribution.

Our computational domain is the Cartesian box defined by , and which we cover with computational zones. We use PLUTO in its dimensionally-unsplit mode using the hllc Riemann solver, linear spatial interpolation, and RK2 time-stepping. We employ periodic boundary conditions in the and directions and reflecting boundary conditions at the bottom of the atmosphere (). At the top of the atmosphere (), the boundary conditions on and are given by extrapolating an isothermal hydrostatic equilibrium solution into a ghost zone, whereas the velocity is set by a standard zero-gradient outflow () condition.

We record the full state of the atmosphere every . In addition to the standard hydrodynamic variables, we also keep track of two additional passive scalar tracers, and , which evolve according to

 ∂∂t(ρμi)+∇⋅(ρμiv)=0. (19)

The first tracer is injected along with the additional heat into the thermal bombs and hence gives us a view of the region that has been directly “polluted” by the debris of the AGN-driven cavities — in real systems, this would be the region containing the relativistic particles that started out in the AGN jet/cocoon if we were to ignore the diffusive transport of these particles. The second tracer is laid down at time into thin horizontal layers at , . The evolution of this tracer allows us to gain a quasi-Lagrangian viewpoint of perturbations and disturbances within the atmosphere, at least until layers are mixed together.

We consider three model atmospheres distinguished by the value of . Our fiducial case has . In this model, -modes are only weakly trapped. We also consider cases of strong () and very strong () trapping. For all three models, we choose . Each model was run on 1024-cores of the Deepthought2 cluster at the University of Maryland, College Park, with a run time of approximately 30,000 core hours per model.

## 4. Results for the Weak-Trapping Case

### 4.1. Basic Dynamics and Evolution

We begin by discussing the evolution of the weak-trapping case (). Following each injection event, the strong overpressure within the heated sphere drives the supersonic inflation of a low-density bubble, launching a shock that starts strong and rapidly weakens (an example is shown in the left panel of Fig. 1). The rapid expansion of the bubble stops when it reaches approximate pressure balance with its surroundings, after which point it begins to rise buoyantly. As is well known (e.g., Brüggen & Kaiser, 2002; Friedman, Heinz & Churazov, 2012), the buoyant rise of the bubble induces Kelvin-Helmholtz and Richtmyer-Meshkov instabilities, which transform the bubble into a vortex ring followed by a turbulent wake (right panel of Fig. 1). As the bubble-ring rises, it mixes with its lower-entropy surroundings and eventually stops rising when its entropy has been reduced to the background level. The repeated generation of these buoyant bubbles quickly forms a turbulent layer of bubble-debris. The mixing process for the later generations of bubbles is more effective because they are rising into an already turbulent region. In our fiducial model, the mixing and equilibration occur below and the bubble-debris (as defined by tracer ) always remains confined to this region. From hereon, we shall refer to this region as the bubble-debris layer. From early times, the bubble-debris layer can be seen oscillating at approximately the Brunt-Väisälä frequency.

Above , the atmosphere remains pristine in the sense that it is unpolluted by the bubble debris. However, it is still agitated by the weak shocks and sound waves from the AGN driving as well as by the buoyant oscillations of the bubble-debris zone below. A central question of this study is whether this pristine region is driven into a state of turbulence. We assess this by examining the vorticity, . As discussed in Section 2, we expect -modes in the atmosphere to generate - and -components of vorticity, but, at linear order in the velocity perturbations, not the -component . Small-scale vorticity in all directions will be produced, however, by non-linear interactions characterizing turbulence. Thus, the production of small-scale and volume-filling is a useful probe of the development of turbulence. The values of on a slice in the -plane at representative snapshots in time are shown in Fig. 2. At very early times (Fig. 2a), we see the generation of -vorticity during the rise of an individual bubble. Later (Fig. 2b), there is significant vorticity generation by the turbulence within the bubble-debris layer but very little in the ambient medium (apart from a few streaks of vorticity high in the atmosphere formed when weak shocks intersect). Then, between times (Fig. 2c–e), there is a volume-filling generation of small scale in the ambient medium. This is a strong indication for the generation of volume-filling turbulence in the ambient medium.

### 4.2. Decay of g-modes and Transition to Turbulence

The transition to turbulence can be seen by studying the time-evolution of the passive tracer . In Fig. 3, we consider a “core sample” drilled though the middle () of the simulated atmosphere, and show the variation of the passive tracer as a function of height and time. The tracer layers in the lower regions of the atmosphere () are rapidly mixed by the action of the buoyant bubbling. We also see that the highest parts of the atmosphere (those with initial positions above ) are rapidly pushed towards the top boundary of the simulation domain due to the heating and slight resultant expansion of the bubble-debris. At intermediate heights, however, we clearly see the onset of oscillatory motions. Determining the frequency of these oscillations by eye, we find good agreement (at the 10–20% level) with the Brunt-Väisälä frequency as given by equation (4). Examining maps of the three components of the vorticity (e.g. Fig. 4), we find that these oscillations are associated with the production of vorticity in the -plane but not in the -direction. Thus, we are secure in identifying them as -modes.

From Fig. 3, we clearly see that the -modes undergo a transition to turbulence, at which point the layers start to mix with the surrounding ICM. Another view of this process is presented in Fig. 5, which shows the tracer across -slices in snapshots of time that span the transition to turbulence. At the beginning of this sequence (; Fig. 5a), the layers are still well-defined and display a distinct sine-wave profile — it is clear that the dominant -modes have the largest non-trivial horizontal wave-numbers that will fit within our domain. As time progresses, the -modes “break” (see, e.g., the layer at in Fig. 5c) and turbulence is driven. From approximately onwards, turbulence as signaled by mixing of layers and the presence of small-scale is present at all altitudes in the atmosphere.

Thus, at least in our fiducial atmosphere, we have verified two of the central tenets of the turbulent heating scenario — that the AGN activity can indeed drive turbulence in the ambient ICM and that it is the non-linear decay of -modes that gives rise to this turbulence. Having established these basic properties of our model atmosphere, we turn to a quantitative discussion of its energetics.

### 4.3. Comparison of Compressible and Incompressible Motions: Inefficiency of AGN-Driven Turbulence

The dynamics of our atmosphere involves both compressible (weak shocks and sound waves) and incompressible (-modes and turbulence) motions. To enable a quantitative study of the energetics, we need to decompose the velocity field into its compressible and incompressible components,

 v=vc+vi, (20)

where

 vc = ∇ψ, (21) vi = ∇×a, (22)

and we determine the scalar field by solving the Poisson equation,

 ∇2ψ=∇⋅v. (23)

(The Cartesian nature of our simulation and the periodicity in the and directions greatly simplify the decomposition by allowing us to utilize standard Fast Fourier Transform methods.)

Once this decomposition is performed, we can examine the kinetic energy of the compressible and incompressible components for both the bubble-debris layer and the ambient atmosphere:

 Kbub,i/c = ∫V12ρ|v2i/c|dV(μ1≥ϵ), (24) Kamb,i/c = ∫V12ρ|v2i/c|dV(μ1<ϵ), (25)

where is taken to be a small number () that, when compared with the tracer , operationally divides volume into the bubble-debris layer and the ambient/unpolluted medium222We note that, due to the small but finite numerical diffusion in PLUTO, the whole computational domain becomes contaminated by very low levels () of the tracer . However, the boundary of the true bubble-debris layer is well delineated and our discussion is insensitive to the precise choice of provided it is in the range ..

Figure 6 shows the compressible and incompressible components of the kinetic energy for the bubble-debris layer and ambient medium in our fiducial model. The evolution divides into early times (during the active AGN driving phase) and late times (when the driving is shut off). During the early/active phase, the kinetic energy of the bubble-debris layer is dominated by incompressible motions. Combining this with the vorticity information discussed above, it is clear that this region is turbulent from an early stage. When we compare the detailed time-evolution of the kinetic energy of incompressible motions in the bubble-debris layer with the timing of the injection events (Fig. 6b), we find that each injection deposits approximately units of kinetic energy into the incompressible motions of the bubble ash; this amounts to between 3–10% of the total injected energy. The injected kinetic energy then decays on a timescale () that is comparable to the Brunt-Väisälä period, suggesting that the freshly-injected turbulent motions in the bubble-debris layer are thermalized on a buoyancy timescale. Compressible modes are always subdominant within the bubble-debris layer because they escape this relatively narrow region on a sound crossing time.

In the ambient material during this early/active phase, the situation is essentially reversed. The kinetic energy is dominated by the compressible modes, with each injection depositing approximately units of energy into these modes. The incompressible modes are sub-dominant during this phase, with a kinetic energy that grows steadily with time until approximately , at which time it appears to saturate. We note that this saturation time coincides with the transition to turbulence discussed in Section 4.2. At early times, there is direct contact between freshly injected bubbles and the ambient material, hence we see a direct correspondence between injection events and spikes in the energy of these ambient incompressible modes. The rapid decay of these spikes likely reflects the mixing of the bubbles and ICM (and hence a reassignment of fluid elements from ambient gas to the bubble-debris layer). At later times, a complete bubble-debris layer forms that “shields” the ambient ICM from the early stages of evolution of a given bubble, thus ending that correspondence.

At late times, after the cessation of AGN driving, the kinetic energy in the compressible modes of both the bubble-debris layer and the ambient ICM quickly diminishes as the sound waves escape the upper (open) boundary (Fig. 6a). The whole atmosphere rapidly becomes dominated by incompressible modes, with approximately twice as much kinetic energy residing in the bubble-debris. The absolute magnitude of the kinetic energies reveal a result of central importance — the transfer of energy from the AGN injection through to the turbulence of the ambient medium is very small, significantly less than 1%.

### 4.4. Power Spectra of the Turbulent Ambient ICM

Before leaving the fiducial case, we briefly discuss the energy power spectrum of the turbulence that we have noted. In our quick treatment here, we only examine the dependence of these power spectra on the horizontal wavenumber — we focus our attention on the power-spectrum of incompressible motions as a function of within horizontal cross-sections of the domain, i.e. at some given . We do not discuss the dependence on the vertical wavenumber due to the heterogeneity and lack of periodicity in the -direction However, motivated by the expectation that the turbulence may be stratified (Zhuravleva et al., 2014a, b), we examine separately the energy power spectrum of horizontal and vertical incompressible motions,

 Eh(k) = πkh[|~vix(kh)|2+|~viy(kh)|2], (26) Ez(k) = πkh|~viz(kh)|2, (27)

respectively, where

 ~viα(kh)=∫viα(x)eik⋅xd3x,α∈{x,y,z}, (28)

and

 kh=kx^x+ky^y. (29)

Figure 7 shows the late-time spectra spectra and for a cross-sectional slice of the atmosphere at . This height is chosen as a representative location where the -mode decay clearly drives turbulence. Both spectra have the same shape and approximate normalization, with a steep decay on scales smaller that . This scale corresponds to only computational zones and so we attribute this steepening to the effects of numerical dissipation. On scales between , both power-spectra have approximately the Kolmogorov slope of . The emergence of a broad power-law spectrum confirms the visual impression that a turbulent state is achieved, with very similar magnitudes of horizontal and vertical motions suggesting approximate isotropy. For the largest scales, both power-spectra are steeper than the Kolmogorov slope which we attribute to the continuing presence of long-wavelength -modes.

## 5. Behavior of the Atmosphere with g-mode trapping

We now turn to a discussion of atmospheres with strong -mode trapping, i.e., our cases. One can imagine competing effects of mode trapping. On the one hand, trapping prevents the -modes from propagating out of the bubble-debris zone in the first place, diminishing the possibility of AGN-driven turbulence. On the other hand, those modes that can propagate into the ambient medium will be confined in the lower regions of the atmosphere and hence may be locally enhanced and promote turbulence. We will see that, at least in our model atmospheres, the first of these factors is dominant and that trapping acts to suppress the AGN driving of turbulence.

Figure 8 shows the time evolution of a core-sample through the tracer layers as a function of time for these two cases. Firstly we discuss the case. As in the fiducial case, we see strong mixing within the bubble-debris layer which now reaches up to by the end of the simulation. The increased vertical extent of the bubble-debris layer is due to flatter entropy profile — the mixed debris of the bubbles rise higher before equalizing to background entropy level. The tracer layers reveal oscillations in the ambient medium starting in the lower layers and, at later times, extending to the higher layers. As in the fiducial case, these oscillations are found to have approximately the Brunt-Väisälä frequency and are associated with the production of vorticity in the -plane, hence we associate them with -modes. In contrast to the fiducial case, there is only a very modest level of breaking of these modes. Only very weak turbulence is generated, with most of the energy remaining in the -modes by the end of the simulation.

Many of these trends hold for the atmosphere with the strongest trapping, . The bubble-debris layer now extends to due to an even flatter entropy profile, and at late times (once driving stops) the ambient upper atmosphere supports low-frequency -modes that remain essentially linear through the extent of the simulation. No turbulence is generated within the ambient atmosphere.

The role of -mode trapping is revealed most directly through an analysis of the incompressible component of the kinetic energy in the ambient matter, as shown in Fig. 9. During the early/driven phase, these three runs show very similar levels of incompressible dynamics. However, once driving ceases, it is apparent that the atmosphere becomes progressively more quiescent as the trapping increases (i.e., decreases).

Finally, by way of a technical comment, we note that, during the driving phase () the atmosphere is seen to expand and, at all heights, to undergo high-frequency oscillations (Fig. 8, right panel). These oscillations are much faster than the Brunt-Väisälä frequency (with a period of approximately compared with , evaluating at ), have no associated vorticity generation, and are coherent across heights. We identify this as the fundamental tone of a global acoustic mode resonating within the “cavity” defined by our reflecting boundary condition at and velocity outflow (zero-gradient) boundary condition at . This, of course, flags these oscillations as being artificial and an artifact of our finite computational domain. They do, however, appear to dissipate quickly once the driving stops, leaving the -modes which are of primary interest to us.

The fact that these global acoustic oscillations are not seen in the other cases () can be understood in terms of the density profile characterizing each run. The (initial) density ratio between the bottom and top boundaries is

 ρtopρbottom=(1+8hz0)−z0/h, (30)

and, while the density profile evolves somewhat during the driving phase, it retains a strong memory of this initial state. Evaluating the ratioes for the three cases gives respectively. For sufficiently extreme density ratios, the atmosphere at the top of the box is so rarified that the dynamics of the lower atmosphere is insensitive to the existence of the upper computational boundary (formally, the eigenmode characterizing global acoustic modes becomes insensitive to the upper boundary condition). Our results demonstrate that the and cases are in this regime. For the case, on the other hand, the upper atmosphere maintains sufficient density to enable an acoustic resonance to be established.

## 6. Discussion

The goal of this work is to examine, or at least start the examination of, the basic fluid dynamics underlying the AGN driving of ICM turbulence. Working within the framework of ideal hydrodynamics, we find that AGN activity (which we model as explosive events at the base of the atmosphere) can indeed drive -modes into the atmosphere, which subsequently decay into volume-filling turbulence. Within our model, however, this process is inefficient — less than 1% of the injected energy eventually ends up as the kinetic energy of turbulence in the ambient ICM, whereas the vast majority of the energy is lost in compressible modes or rapidly thermalized within the bubble-debris layer. Strong trapping of the -modes further decreases this efficiency and can even render them so weak as to not generate turbulence at the Reynolds numbers afforded by the resolution of our simulations.

Taken at face value, this inefficiency of turbulent dissipation is inconsistent with observations. The energetics of AGN in cool core clusters are usually estimated through a determination of the product , where is the (measured) pressure of the ICM surrounding a cavity of (measured) volume . On purely energetic grounds, the AGN needs to transfer of heat per cavity to the cluster core in order to balance radiative cooling (Bîrzan et al., 2004; Panagoulia et al., 2014). In our model, we find that the initial expansion phase of the bubbles (which take the bubbles from an overpressure by a factor of 6 down to pressure equilibrium with the surrounding medium) only increases their volume by a factor of 2.5–3, with approximately 70% of the energy being used to do work on and directly shock-heat the surrounding shell of gas. Thus, examining just the energy content of a cavity will under-estimate the injected power by a factor of . However, even adopting a generous value of 1% for the efficiency with which the AGN transfers its energy to turbulence, we still infer that the AGN transfers only per cavity to the ICM turbulence.

Of course, our model is highly idealized in order to facilitate a detailed study of the hydrodynamics, and we must critically assess the degree to which the idealizations affect our results.

Firstly, real clusters are spherical, not plane-parallel, and geometry is always important for the detailed dynamics. For example, in a spherical cluster with a realistic radial profile of pressure, acoustic modes generated at the center will weaken (i.e., will decrease) as the wave propagates outwards and carries an energy flux that diminishes as . In our plane-parallel atmospheres with periodic boundary conditions in the - and -directions, on the other hand, the upper atmosphere is seeing quasi-planar acoustic modes which can re-steepen and become shocks. Indeed, examining the entropy of our tracer layers in the fiducial case reveals heating of the uppermost regions of the atmosphere that we attribute precisely to the action of steepened acoustic modes.

Secondly, we have adopted the framework of inviscid ideal hydrodynamics. One consequence of this is that we may underestimate how long the structural integrity of the buoyantly rising bubbles is preserved (Reynolds et al., 2005), which may, in turn, diminish the transference of energy from the bubble-induced upwell into the -modes. Magnetic fields can play an important role here. While the magnetic field in the bulk of the ICM is expected to be weak (with magnetic pressures that are approximately 1% of the gas pressure), the draping of field lines across the rising-cap of a buoyant bubble may significantly reduce the disruption of the bubble and subsequent mixing with the ambient ICM (Ruszkowski et al., 2007; Dursi & Pfrommer, 2008), potentially impacting the efficiency with which the AGN drives ICM turbulence. Even if magnetic forces are weak enough so as to be unimportant, the magnetic field strongly affects transport processes, introducing a very strong anisotropy into the viscous stress tensor and thermal conductivity. In general terms, the anisotropic viscosity and conduction can affect the convective stability of the ICM (Balbus, 2000; Quataert, 2008; Balbus & Reynolds, 2010; Kunz et al., 2011, 2012) and both the conduction and the viscosity can dramatically enhance the dissipation of compressible modes. Of more direct interest to the results of this paper, the combined action of the anisotropic conduction and the radiative cooling can drive -modes to be overstable (Balbus & Reynolds, 2010); this maybe be one possible route for increasing the efficiency of energy transfer from the AGN to the ambient ICM.

## 7. Conclusions

We have constructed model atmospheres tailored to test the physical assumptions underlying the AGN-induced turbulent heating model of cooling-core clusters. We find that AGN-like disturbances can indeed induce volume-filling turbulence in the ambient ICM via the decay of -modes. In our ideal hydrodynamic model, however, we find that the transfer of energy from the AGN-injection process into the turbulence of the ambient gas is quite inefficient (less than 1% for our fiducial models). For the models assessed here, the trapping of -modes acts to decrease the efficiency further because the excitation of turbulence in our model relies on -mode propagation from the cluster center where they are launched by explosive events to the bulk of the ICM core. We stress, however, that this is just the first step in the theoretical assessment of the AGN-driven turbulent heating model and certainly not the last word. Next steps must examine models that include the action of radiative cooling, magnetic fields, (anisotropic) thermal conduction and (Braginskii) viscosity. While these are ingredients that, from first principles, should be present, we advocate a systematic and controlled study of models of increasing complexity. Only then can the essential underlying physical mechanisms be revealed.

CSR is grateful for financial support from the Simons Foundation (through a Simons Fellowship in Theoretical Physics), a Sackler Fellowship (hosted by the Institute of Astronomy, Cambridge), and the National Science Foundation (under grant AST1333514). SAB acknowledges support from the Hintze Family Charitable Foundation, and from the Royal Society in the form of a Wolfson Research Merit Award. CSR thanks the Department of Physics at the University of Oxford and the Institute of Astronomy at the University of Cambridge for their hospitality during extended visits while this work was being performed. The simulations presented in this paper were performed on the Deepthought2 cluster, maintained and supported by the Division of Information Technology at the University of Maryland College Park.

## References

• Balbus (2000) Balbus, S. A., 2000, ApJ, 534, 420
• Balbus & Reynolds (2010) Balbus, S. A., & Reynolds, C. S., 2010, ApJ, 720, L97
• Binney & Cowie (1981) Binney, J., & Cowie, L. L., 1981, ApJ, 247, 464
• Bîrzan et al. (2004) Bîrzan, L., Rafferty, D. A., McNamara, B. R., Wise, M. W., & Nulsen, P. E. J., 2004, ApJ, 607, 800
• Bogdanović et al. (2009) Bogdanović, T., Reynolds, C. S., Balbus, S. A., & Parrish, I. J., 2009, ApJ, 704, 211
• Brüggen & Kaiser (2002) Brüggen, M., & Kaiser, C. R., 2002, Nature, 418, 301
• Burns (1990) Burns, J. O., 1990, AJ, 99, 14
• Churazov et al. (2001) Churazov, E., Brüggen, M., Kaiser, C. R., Böhringer, H., & Forman, W., 2001, ApJ, 554, 261
• Churazov et al. (2000) Churazov, E., Forman, W., Jones, C., & Böhringer, H., 2000, A&A, 356, 788
• Churazov et al. (2004) Churazov, E., Forman, W., Jones, C., Sunyaev, R., & Böhringer, H., 2004, MNRAS, 347, 29
• Churazov et al. (2002) Churazov, E., Sunyaev, R., Forman, W., & Böhringer, H., 2002, MNRAS, 332, 729
• Conroy & Ostriker (2008) Conroy, C., & Ostriker, J. P., 2008, ApJ, 681, 151
• Dursi & Pfrommer (2008) Dursi, L. J., & Pfrommer, C., 2008, ApJ, 677, 993
• Enßlin & Vogt (2006) Enßlin, T. A., & Vogt, C., 2006, A&A, 453, 447
• Fabian (1994) Fabian, A. C., 1994, ARA&A, 32, 277
• Fabian (2012) Fabian, A. C., 2012, ARA&A, 50, 455
• Fabian et al. (2005) Fabian, A. C., Reynolds, C. S., Taylor, G. B., & Dunn, R. J. H., 2005, MNRAS, 363, 891
• Fabian et al. (2000) Fabian, A. C., et al., 2000, MNRAS, 318, L65
• Friedman, Heinz & Churazov (2012) Friedman, S. H., Heinz, S., & Churazov, E., 2012, ApJ, 746, 112
• Fujita (2005) Fujita, Y., 2005, ApJ, 631, L17
• Graham, Fabian & Sanders (2008) Graham, J., Fabian, A. C., & Sanders, J. S., 2008, MNRAS, 386, 278
• Guo & Oh (2008) Guo, F., & Oh, S. P., 2008, MNRAS, 384, 251
• Heinz et al. (2002) Heinz, S., Choi, Y.-Y., Reynolds, C. S., & Begelman, M. C., 2002, ApJ, 569, L79
• Heinz & Churazov (2005) Heinz, S., & Churazov, E., 2005, ApJ, 634, L141
• Hillel & Soker (2015) Hillel, S., & Soker, N., 2015, ArXiv e-prints
• Jones et al. (2002) Jones, C., Forman, W., Vikhlinin, A., Markevitch, M., David, L., Warmflash, A., Murray, S., & Nulsen, P. E. J., 2002, ApJ, 567, L115
• Kim & Narayan (2003) Kim, W.-T., & Narayan, R., 2003, ApJ, 596, 889
• Kunz et al. (2012) Kunz, M. W., Bogdanović, T., Reynolds, C. S., & Stone, J. M., 2012, ApJ, 754, 122
• Kunz et al. (2011) Kunz, M. W., Schekochihin, A. A., Cowley, S. C., Binney, J. J., & Sanders, J. S., 2011, MNRAS, 410, 2446
• Loewenstein, Zweibel & Begelman (1991) Loewenstein, M., Zweibel, E. G., & Begelman, M. C., 1991, ApJ, 377, 392
• McNamara & Nulsen (2007) McNamara, B. R., & Nulsen, P. E. J., 2007, ARA&A, 45, 117
• McNamara & Nulsen (2012) McNamara, B. R., & Nulsen, P. E. J., 2012, New Journal of Physics, 14, 055023
• Mignone et al. (2007) Mignone, A., Bodo, G., Massaglia, S., Matsakos, T., Tesileanu, O., Zanni, C., & Ferrari, A., 2007, ApJS, 170, 228
• Million et al. (2010) Million, E. T., Werner, N., Simionescu, A., Allen, S. W., Nulsen, P. E. J., Fabian, A. C., Böhringer, H., & Sanders, J. S., 2010, MNRAS, 407, 2046
• Omma et al. (2004) Omma, H., Binney, J., Bryan, G., & Slyz, A., 2004, MNRAS, 348, 1105
• Panagoulia et al. (2014) Panagoulia, E. K., Fabian, A. C., Sanders, J. S., & Hlavacek-Larrondo, J., 2014, MNRAS, 444, 1236
• Peterson & Fabian (2006) Peterson, J. R., & Fabian, A. C., 2006, Phys. Rep., 427, 1
• Pfrommer (2013) Pfrommer, C., 2013, ApJ, 779, 10
• Pfrommer, Chang & Broderick (2012) Pfrommer, C., Chang, P., & Broderick, A. E., 2012, ApJ, 752, 24
• Quataert (2008) Quataert, E., 2008, ApJ, 673, 758
• Randall et al. (2015) Randall, S. W., et al., 2015, ArXiv e-prints
• Rebusco et al. (2006) Rebusco, P., Churazov, E., Böhringer, H., & Forman, W., 2006, MNRAS, 372, 1840
• Reynolds, Heinz & Begelman (2002) Reynolds, C. S., Heinz, S., & Begelman, M. C., 2002, MNRAS, 332, 271
• Reynolds et al. (2005) Reynolds, C. S., McKernan, B., Fabian, A. C., Stone, J. M., & Vernaleo, J. C., 2005, MNRAS, 357, 242
• Ruszkowski & Begelman (2002) Ruszkowski, M., & Begelman, M. C., 2002, ApJ, 581, 223
• Ruszkowski, Brüggen & Begelman (2004) Ruszkowski, M., Brüggen, M., & Begelman, M. C., 2004, ApJ, 611, 158
• Ruszkowski et al. (2007) Ruszkowski, M., Enßlin, T. A., Brüggen, M., Heinz, S., & Pfrommer, C., 2007, MNRAS, 378, 662
• Schekochihin & Cowley (2006) Schekochihin, A. A., & Cowley, S. C., 2006, Physics of Plasmas, 13, 056501
• Zhuravleva et al. (2014a) Zhuravleva, I., et al., 2014a, Nature, 515, 85
• Zhuravleva et al. (2014b) Zhuravleva, I., et al., 2014b, ApJ, 788, L13
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters