<figure><img src="image_1.png" /></figure>

HAL
open science

# AGN jet-driven stochastic cold accretion in cluster cores

Deovrat Prasad, Prateek Sharma, Arif Babul

► To cite this version:

Deovrat Prasad, Prateek Sharma, Arif Babul. AGN jet-driven stochastic cold accretion in cluster cores. Monthly Notices of the Royal Astronomical Society, 2017, 471 (2), pp.1531-1542. 10.1093/mn-ras/stx1698. hal-01584682

HAL Id: hal-01584682

https://hal.science/hal-01584682v1

Submitted on 26 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# AGN jet-driven stochastic cold accretion in cluster cores

Deovrat Prasad,<fnref n="1" />* Prateek Sharma<fnref n="1" />* and Arif Babul<fnref n="2" />,<fnref n="3" />,<fnref n="4" />*

Accepted 2017 July 4. Received 2017 July 3; in original form 2016 November 4

## ABSTRACT

Several arguments suggest that stochastic condensation of cold gas and its accretion on to the central supermassive black hole (SMBH) is essential for active galactic nuclei (AGNs) feedback to work in the most massive galaxies that lie at the centres of galaxy clusters. Our 3-D hydrodynamic AGN jet-ICM (intracluster medium) simulations, looking at the detailed angular momentum distribution of cold gas and its time variability for the first time, show that the angular momentum of the cold gas crossing ≲1 kpc is essentially isotropic. With almost equal mass in clockwise and counterclockwise orientations, we expect a cancellation of the angular momentum on roughly the dynamical time. This means that a compact accretion flow with a short viscous time ought to form, through which enough accretion power can be channeled into jet mechanical energy sufficiently quickly to prevent a cooling flow. The inherent stochasticity, expected in feedback cycles driven by cold gas condensation, gives rise to a large variation in the cold gas mass at the centres of galaxy clusters, for similar cluster and SMBH masses, in agreement with the observations. Such correlations are expected to be much tighter for the smoother hot/Bondi accretion. The weak correlation between cavity power and Bondi power obtained from our simulations also matches observations.

**Key words:** galaxies: clusters: intracluster medium – galaxies: haloes – galaxies: jets.

## 1 INTRODUCTION

There is a general consensus that the intracluster medium (ICM) in low-entropy galaxy cluster cores is able to remain in rough global thermal equilibrium because the powerful jets from the central active galactic nucleus (AGN) inject sufficient energy to compensate for radiative losses (for a review, see McNamara & Nulsen 2007, 2012; Fabian 2012; Soker 2016). There are, however, a number of key issues associated with this 'radio feedback schema that have yet to be firmly pinned down, one of which is: how is the central supermassive black hole (SMBH) able to transport sufficient energy to the cluster core (~10s of kpc) at a time-scale shorter than the core cooling time (≈0.5 Gyr). There are two scenarios proposed for AGN feedback: (1) hot/Bondi accretion that invokes spherical accretion with the mass flow rate set by the hot gas density and temperature on the cluster core scales (Bondi 1952); and (2) the competing cold mode accretion model, in which a fraction of the gas in the cluster cores condenses into cold clouds and filaments (when the ratio of the cooling time to the free-fall time [$t_{\text{cool}}/t_{\text{ff}}$] falls below a threshold value in the range 5–20; e.g. McCourt et al. 2012;

Voit et al. 2015; but see Hogan et al. 2017a), which then decouple from surrounding hot ICM, rain down upon the central SMBH and power AGN feedback strong enough to prevent catastrophic cooling in the core.

Of the two, the Bondi accretion model, though more frequently invoked in the literature, has a number of shortcomings. For one, Bondi accretion is only applicable if the accretion flow is spherically symmetric and single phase, both conditions that are at odds with the multiphase gas with angular momentum observed in several cool cores (e.g. see Macchetto et al. 1997; McDonald et al. 2010; McNamara, Rohanizadegan & Nulsen 2011; Tremblay et al. 2016). Moreover, for a polytropic hot flow $p\alpha\rho^\gamma$, the mass accretion rate is given by (Shu 1992; Frank, King & Raine 2002)

$$ 
\begin{aligned}
\dot{M}_B &= 4\pi\lambda(\gamma)G^2 M_{BH}^2 \frac{\rho_{\infty}}{c_{s,\infty}^3} \\
&\approx 0.0015 M_{\odot} \text{yr}^{-1} M_{BH,9}^2 n_{\infty,0.1} T_{\infty,\text{keV}}^{-3/2}
\end{aligned}
\quad (1) $$

where $c_{s,\infty} = [\gamma k_B T_\infty / \mu m_p]^{1/2}$ is the sound speed in the ambient medium far from the black hole; the mean particle mass $\mu = 0.62$; $M_{BH,9}$ is the black hole mass in units of $10^9 M_\odot$; $n_{\infty,0.1}$ and $T_{\infty,\text{keV}}$ are the ambient particle number density and temperature scaled to 0.1 cm$^{-3}$ and 1 keV, respectively; and since $\lambda(\gamma)$ varies weakly with $\gamma$ ($\lambda$ varies from 1/4 for $\gamma = 5/3-1.1$ for $\gamma = 1$), we adopt $\gamma = 5/3$

<footnote n="1">¹Department of Physics & Joint Astronomy Programme, Indian Institute of Science, Bangalore 560012, Karnataka, India</footnote>

<footnote n="2">²Department of Physics and Astronomy, University of Victoria, Victoria, BC V8P 1A1, Canada</footnote>

<footnote n="3">³Institute of Computational Science, Centre for Theoretical Astrophysics and Cosmology, University of Zurich, Winterthurerstrasse 190, 8057, Zurich, Switzerland</footnote>

<footnote n="4">⁴Institut d'Astrophysique de Paris-98 bis boulevard Arago-75014, Paris, France</footnote>

<footnote>* E-mail: deovrat@physics.iisc.ernet.in (DP); prateek@physics.iisc.ernet.in (PS); babul@uvic.ca (AB)</footnote>
for convenience. In three nearby systems where both the ambient conditions governing the Bondi flow rate and the mass accretion rate on to the black hole can be deduced, the former is at least two orders of magnitude too large. These nearby examples are Sgr A*, the Galactic Centre BH [Bondi rate of $10^{-6}-10^{-5} M_{\odot}yr^{-1}$ (Baganoff et al. 2003) versus the accretion rate at $\sim 10R_{Sch}$ of $\lesssim 10^{-7} M_{\odot}yr^{-1}$ (Marrone et al. 2007 and references therein)]; NGC 3115 [Bondi rate of $2.2 \times 10^{-2} M_{\odot}yr^{-1}$ versus the accretion rate on to the BH at least two orders of magnitude smaller (Wong et al. 2011, 2014)]; and M87 [Bondi rate of 0.1–0.2 $M_{\odot}yr^{-1}$ versus the accretion rate on to the BH of $\lesssim 10^{-3} M_{\odot}yr^{-1}$ (Russell et al. 2015)].

Nemmen & Tchekhovskoy (2015) find that the median efficiency (defined as the ratio of the observed jet power and $\dot{M}_{BH}c^2$; $\dot{M}_{BH}$ is the SMBH accretion rate) required to supply the cavity power of nearby radio galaxies is 300 per cent, assuming $\dot{M}_{BH} \sim 0.01 \dot{M}_{B}$ typical of similar nearby hot accretion flows. This implies that hot accretion, with reasonable assumptions, is insufficient for powering most of the cavities observed in cluster cores. At the extreme end of the spectrum, the estimated jet/cavity power in some of the galaxy clusters (e.g. MS0735, Cen A) is found to be much larger than $0.1 \dot{M}_{B}c^2$; i.e. the Bondi mass flow rate is insufficient to power the observed X-ray cavities even if all the hot gas in the Bondi flow is accreted by the SMBH (McNamara et al. 2011). And finally, the accretion rate in the Bondi regime takes a long time to adjust to the conditions $(n, T)$ changing in the cluster cores ($\sim 1$ kpc) at the core cooling time-scale approximately few 100 Myr (e.g. Soker, Sternberg & Pizzolato 2009).

The cold mode feedback model sidesteps many of the challenges associated with the hot/Bondi model. Very briefly (we refer the readers to Sharma et al. 2012a; Voit et al. 2015; Prasad, Sharma & Babul 2015 [hereafter PSB15 for a more expansive discussion), the model is based on the realization that even if the ICM in the cluster core exists in rough global thermal balance, if the ratio of the cooling time and the gravitational free-fall time ($t_{cool}/t_{ff}$) falls below a critical value ($\sim 10$ or so; Gaspari, Ruszkowski & Sharma 2012; Sharma et al. 2012a; Li et al. 2015; Prasad, Sharma & Babul 2015; Voit et al. 2015; Choudhury & Sharma 2016), the gas will become susceptible to local thermal instabilities and fragment, leading to the formation of a multiphase medium consisting of cold dense clouds condensing from the hot diffuse ICM itself. The clouds then fall towards central AGN, resulting in increased SMBH accretion and feedback that, in turn, quenches runaway cooling in the cluster core. The main question with this model is: will the infalling cold gas, which has the non-zero angular momentum, end up forming a viscosity-mediated standard accretion flow (SAF) in which the gas will flow inwards on a (long) viscous accretion time-scale, or will the gas circularize sufficiently close to the SMBH and accrete on a time-scale shorter than the core cooling time.

A particularly promising accretion mechanism, overcoming the angular momentum barrier by stochastically feeding cold gas clumps and fueling AGN activity, was invoked by Pizzolato & Soker (2005), Hopkins & Hernquist (2006), Pizzolato & Soker (2010), and Nayakshin & King (2007). The AGNs are stochastically fed by cloud clumps that cancel each other's angular momentum through inelastic collisions, leading to the formation of a small transient disc accreting on to the black hole (Nayakshin & King 2007). This turbulent accretion mechanism of cold gas clouds on to SMBH was later called cold chaotic accretion (Nayakshin, Power & King 2012). Using isothermal simulations of idealized turbulent accretion flow over central 100 pc, Hobbs et al. (2011) show that dense gas accretes ballistically on to the central SMBH

at a rate few orders of magnitude larger than the case without turbulence.

Gaspari, Ruszkowski & Oh (2013) carried out simulations of cluster cores (although only run for 40 Myr, a duration much shorter than the core cooling time) with idealized turbulent driving (which is fixed in time), and showed that the condensation of cold clumps from the hot ICM can boost the mass accretion rate on to the SMBH by $\sim 100$ compared to the Bondi accretion rate. These simulations were later generalized to include an initially rotating ICM (Gaspari, Brighenti & Temi 2015) and cooling down to much lower temperatures allowing for three phases (hot, atomic and molecular; Gaspari, Temi & Brighenti 2017), although still using idealized turbulence and running for much less than a core cooling time. The process of colliding cold gas clouds losing the angular momentum and boosting SMBH accretion has been called chaotic cold accretion. Other terminology has also been used in the past (ballistic, stochastic, forced accretion, etc.). In this paper we use a related but somewhat different term stochastic cold accretion (SCA) to highlight that a turbulent system such as the ICM should be described statistically (Monin & Yaglom 1971).

While AGN jet-ICM simulations have been run successfully on cosmological time-scales for some time now (e.g. Gaspari et al. 2012; Li et al. 2015; Prasad et al. 2015), we investigate the detailed angular momentum distribution of cold gas and its implications on SMBH accretion in realistic jet-ICM simulations over cosmological time-scales for the first time. Here is a brief outline of the paper. Section 2 briefly presents our numerical setup. In Section 3, we calculate the circularization radius and viscous time of the cold gas based on the standard accretion physics. We show that the angular momentum distribution of cold gas crossing $\sim 1$ kpc is close to isotropic, which will lead to angular momentum cancellation and the formation of an efficient accretion flow. We use the mass distribution of angular momentum to estimate the accretion rate on to the SMBH. In Section 4, we compare some of the results from our simulations against key observations. In Section 5, we discuss the implications of our results, highlighting the need for better understanding of accretion. We summarize our paper in Section 6.

## 2 NUMERICAL SIMULATIONS

In this paper we discuss three 3-D AGN jet-ICM simulations carried out in spherical $(r, \theta, \phi)$ coordinates. The details of the numerical set up are given in PSB15. The three runs are: (i) an initially hydrostatic ICM in a fixed NFW potential (this is the fiducial run in PSB15 with the halo mass of $7 \times 10^{14} M_{\odot}$); (ii) the same run but with the inner and outer radii reduced by half; (iii) and a NFW+BCG (the latter is modelled as a singular isothermal potential with circular velocity $V_c = 350 \text{ km s}^{-1}$) potential run. The NFW+BCG run is carried out with PLUTO (Mignone et al. 2007) code, whereas the NFW runs are done with ZEUS-MP (Hayes et al. 2006). The primary purpose of the NFW+BCG runs is to facilitate a detailed analysis of how the addition of the BCG potential impacts the various aspects of the cold mode accretion model, a topic of considerable interest (cf. Hogan et al. 2017a,b). We will present this analysis in a follow-up paper. Most of the results and analysis in this paper are based on the NFW runs although we occasionally call upon the NFW+BCG simulations to highlight that our key results are robust.

The initial entropy profiles for both NFW and NFW+BCG runs are identical (see equation 7 in PSB15). The initial ICM density for the NFW+BCG run is almost twice of the NFW runs (which have identical initial density profiles). The different runs have different inner ($r_{in}$) and outer ($r_{out}$) radial boundaries (see Table 1). The
<table><caption>**Table 1.** List of runs.</caption><thead><tr><th>Run</th><th>r<sub>in</sub><br>(kpc)</th><th>r<sub>out</sub><br>(kpc)</th><th>ε<sup>a</sup></th><th>Run Time<br>(Gyr)</th><th>M&#x0307;<sub>in,cold/hot</sub><br>(M<sub>☉</sub> yr<sup>-1</sup>)</th><th>M<sub>cold</sub><br>(10<sup>11</sup> M<sub>☉</sub>)</th><th>M&#x0307;<sub>B</sub><br>(M<sub>☉</sub> yr<sup>-1</sup>)<br>PL, β</th><th>M&#x0307;<sup>b</sup><sub>BH,SAF</sub><br>(M<sub>☉</sub> yr<sup>-1</sup>)</th><th>M&#x0307;<sup>b</sup><sub>BH,SCA</sub><br>(M<sub>☉</sub> yr<sup>-1</sup>)</th></tr></thead><tbody><tr><td>NFW<sup>c</sup></td><td>1</td><td>200</td><td>6 × 10<sup>-5</sup></td><td>4</td><td>18.1/3.5</td><td>1</td><td>0.2, 0.007</td><td>0.015</td><td>0.2</td></tr><tr><td>NFW</td><td>0.5</td><td>100</td><td>6 × 10<sup>-5</sup></td><td>5</td><td>10.7/0.9</td><td>5</td><td>0.54, 0.01</td><td>0.1</td><td>1.27</td></tr><tr><td>NFW+BCG</td><td>0.5</td><td>500</td><td>5 × 10<sup>-4</sup></td><td>3</td><td>6.4/2.8</td><td>0.1</td><td>0.008, 0.003</td><td>0.001</td><td>0.03</td></tr></tbody></table>

Notes. The resolution of all runs, done in spherical ($r_{min} \le r \le r_{max}$, $0 \le \theta \le \pi$, $0 \le \phi \le 2\pi$) coordinates, is 256 × 128 × 32. Logarithmic grid is used in the $r$ direction, and a uniform one in others. $\dot{M}_{in,cold/hot}$ is the average cold/hot mass flow rate across $r_{in}$; $\dot{M}_{B}$ is the estimate of average Bondi accretion rate using power-law (PL) and isothermal beta model ($\beta$) extrapolations; $M_{cold}$ is the total cold ($T < 0.1$ keV) gas mass in the simulation domain by the end. $\dot{M}_{in,cold/hot}$, $\dot{M}_{BH,SAF}$, $\dot{M}_{BH,SCA}$ and $\dot{M}_{B}$ are averaged from 1 Gyr till the end of the run.

<sup>a</sup>Jet efficiency relative to $\dot{M}_{in}$ (see equation 6 in PSB15).

<sup>b</sup>$\dot{M}_{BH,SAF/SCA}$ is the SMBH accretion rate calculated using equation (4)/(8).

<sup>c</sup>The fiducial run.

density and pressure at the outer radial boundary are fixed to their initial values. The total mass (including both hot and cold phases) crossing the inner boundary is calculated ($\dot{M}_{in}$), and a fraction $\epsilon$ of $\dot{M}_{in}c^2$ is put back into the ICM in the form of kinetic power ($P_{jet}$) of AGN jets.

All our runs behave in a very similar manner. There are cooling and heating cycles, and the feedback-driven heating events are less frequent for a higher feedback efficiency. Since the feedback efficiency, $\epsilon$, is 10 times higher for the NFW+BCG run, the amount of cold gas and accretion rate through the inner boundary are suppressed relative to NFW runs. The NFW run with $r_{in} = 0.5$ kpc uses the same efficiency parameter ($\epsilon$) as the fiducial run with $r_{in} = 1$ kpc. However, because of a smaller inner radius and non-zero angular momentum of the condensing gas, the mass accretion rate across the inner boundary is lower as some of the gas circulates at $0.5$ kpc $< r < 1$ kpc. Consequently, the feedback energy input is smaller and the run shows a higher accumulation of cold gas in the central regions (see $M_{cold}$ in Table 1). Density and temperature profiles, and other X-ray and jet properties, are discussed in detail in PSB15, and will not be repeated in this paper. Here we focus on the angular momentum distribution of cold gas in our simulations and the plausibility of stochastic accretion on to SMBHs (Section 3), and comparison of simulation results with observations (Section 4).

## 3 ANGULAR MOMENTUM OF COLD GAS AND SMBH ACCRETION

Our simulations, like Gaspari et al. (2012), Li et al. (2015), and Yang & Reynolds (2016), use the total (dominated by cold gas) accretion rate at $\sim 1$ kpc to estimate the feedback energy deposited by jets. One of the major open problems with cold mode feedback is: how does the cold gas originating at $\gtrsim 1$ kpc lose the angular momentum and get rapidly accreted on to the SMBH. This is needed for feedback to prevent a cooling flow. Although it is numerically formidable to resolve both the cluster core and accretion flow on to the SMBH, we attempt to estimate the accretion rate on to the SMBH based on the angular momentum distribution of cold gas within our simulation domain. While the large angular momentum gas forms a massive disc at $\gtrsim 1$ kpc, the smaller angular momentum gas crossing the inner boundary has stochastic angular momentum changing over short time-scales. Rapid angular momentum cancellation in this stochastic cold gas can allow substantial cold gas to be channeled to the SMBH sufficiently fast to prevent a cooling flow (e.g. see Pizzolato & Soker 2010).

### 3.1 Standard accretion estimates

In this section we review the standard accretion physics that we apply later to estimate the accretion rate on to the SMBH (Section 3.3.1). For a given specific angular momentum ($l$), the circularization radius of the gas is given by

$$R_{\text{circ}} \equiv \frac{l^2}{GM_{\text{enc}}} \approx 0.24 \text{ kpc } M_{\text{enc},9} l_{28}^2, \quad (2)$$

where $M_{enc}$ is the mass enclosed (including only SMBH + dark matter contributions; i.e. ignoring gas mass) within the circularization radius. The viscous accretion time at the circularization radius can be expressed in terms of the specific angular momentum,

$$ \begin{aligned} t_{\text{visc}} &\approx \frac{1}{\alpha \Omega (H/R)^2} \approx \frac{l^3}{\alpha G^2 M_{\text{enc}}^2 (H/R)^2} \\ &\approx 1.8 \text{ Gyr } \alpha_{0.1}^{-1} \left(\frac{H}{R}\right)_{0.1}^{-2} M_{\text{enc},9}^{-2} l_{28}^3, \end{aligned} \quad (3) $$

where $\Omega \equiv (GM_{enc}/R^3)^{1/2}$ is the angular velocity, $\alpha$ is the viscosity parameter (Shakura & Sunyaev 1973) and $H/R$ is the ratio of the disc height and radius. We scale the numerical values using $\alpha = 0.1$. This is consistent with results from MHD simulations of magnetized accretion flows, which find that the effective value of $\alpha$ is 0.1 over the bulk of the flow (Hawley & Krolik 2001; Babul, Sharma & Reynolds 2013). Although the realistic $H/R$ for standard AGN discs is $\sim 10^{-3}$ (e.g. see equation 6 in Babul et al. 2013), we scale our results to $H/R = 0.1$ because the standard choice gives a very small SMBH accretion rate. Moreover, stochastic accretion simulations show a thick cold flow ($H/R \gtrsim 0.1$; see Section 3.2) supported by turbulent motions (see Fig. 16 in Hobbs et al. 2011; Gaspari et al. 2013). Incidentally, radiatively inefficient accretion flows (RIAFs) also have similar thick discs (Das & Sharma 2013; Yuan & Narayan 2014).

Fig. 1 shows the viscous time (top panel) and circularization radius (bottom panel), as a function of the specific angular momentum, for the two gravitational potentials (NFW and NFW+BCG) used in our simulations. Note that the inner radius of our computational domain is well beyond the Bondi radius, and the SMBH contribution to gravity is negligible everywhere in the computational domain ($r_{in} \le r \le r_{out}$). The slope of $t_{visc}$ and $R_{circ}$ as a function of $l$ changes at the radius at which the enclosed mass is dominated by the extended potential rather than the SMBH. The arrows in Fig. 1 mark important quantities such as Bondi radius, inner radius of the computational domain and a viscous time of 0.2 Gyr. Accretion flow with a viscous time more than the core cooling time (here taken to be 0.2 Gyr, a conservative estimate) would not be able to respond faster than
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Circularization radius (equation 2; bottom panel) and viscous time (equation 3; top panel) as a function of specific angular momentum ($l$) for NFW+SMBH and NFW+BCG+SMBH potentials. As expected, the SMBH dominates the enclosed mass within the Bondi radius. Magenta arrows represent the Bondi radius, the corresponding $l$ (equation 2) and $l_{\text{visc}}$ (equation 3). Black arrows represent quantities corresponding to the inner simulation radius of 0.5 kpc (dot-dashed line for NFW and dashed line for NFW+BCG) and 1 kpc (solid line). Brown arrows show the viscous time of 0.2 Gyr. These scales are useful to estimate SMBH accretion rates in Section 3.3.</figcaption></figure>

the cooling rate. A faster response of the SMBH accretion rate is required to prevent a cooling flow.

## 3.2 Stochastic cold accretion

In this section we discuss the angular momentum distribution of cold gas in our numerical simulations. We show that the time-averaged angular momentum distribution of cold gas crossing our inner simulation boundary is stochastic, with almost equal mass going around in different directions. This implies that in the region of circularization, angular momentum cancellation will take place almost at the local dynamical time. Moreover, the cold gas with $l < 10^{28} \text{ cm}^2 \text{ s}^{-1}$ changes its mean angular momentum over a time-scale shorter than the core cooling time; i.e. cold gas with small angular momentum is able to respond fast enough to close the feedback loop.

### 3.2.1 Time-averaged angular momentum probability density function

To start with, we emphasize that when considering gas whose specific angular momentum ($l$) distribution is stochastic, it is essential to explicitly account for the vector nature of the specific angular momentum, $l = r \times v$ ($= -rv_{\phi}\hat{\theta} + rv_{\phi}\hat{\phi} = v_{\phi}R\hat{z} - r[v_{\phi}\cos\theta\cos\phi + v_{\phi}\sin\phi]\hat{x} + r[-v_{\phi}\cos\theta\sin\phi + v_{\phi}\cos\phi]\hat{y}$, where $R/r$ is the cylindrical/spherical radius). In standard accretion disc literature, $l$ is generally taken to be $l_z$ since the gas is considered to be coplanar. However, in our simulations the low angular momentum cold gas distribution shows almost equal importance of all the components of $l$.

To examine this, we consider the quantity $(dM/dl)\Delta l$, the incremental mass in the computational domain with the angular mo-

mentum between $l$ and $l + \Delta l$ (we plot the mass distribution in a logarithmic bin $dM/d\log_{10}l = [l/\log_{10}e]dM/dl$, as one can simply read off the mass from the figure). We can choose $l$ to be $|l|$, or $l_x, l_y$ or $l_z$. The left-hand panel of Fig. 2 shows the time-averaged cold mass distribution as a function of $|l_z|$ (both clockwise and counterclockwise rotations are shown) for our three runs. The right-hand panel shows the time-averaged mass distribution as a function of the total specific angular momentum $|l|$. Comparing the $|l|$ and $|l_z|$ distributions, we find a dramatic lack of cold gas with small $|l|$ but a non-negligible mass of cold gas with very low $|l_z|$. This gas with low $|l_z|$ has large $|l_x|$ and $|l_y|$ components, and therefore is not expected to circularize in the $x-y$ plane.

Focusing on gas with $l < l_{\text{in}}$ (corresponding to the circularization radius in equation (2) equal to the inner radial boundary $r_{\text{in}}$), where $l_{\text{in}} \sim 10^{29} \text{ cm}^2 \text{ s}^{-1}$, we find that in all our simulations the low angular momentum gas ($l < l_{\text{in}}$) has roughly similar distribution in $l_z$ and $|l|$, with $dM/dl_z \propto l_z^0$ and $dM/d|l| \propto |l|^{-1}$. In fact, both the clockwise and counterclockwise components of angular momentum distribution have very similar profiles for $l \lesssim l_{\text{in}}$. Since the mass probability density function (pdf) in Fig. 2 does not peak as abruptly for simulations with lower $l_{\text{in}}$ (and smaller $r_{\text{in}}$), we assert that the low-$l$ mass pdfs are robust (in the sense that the low-$l$ pdf will be similar if $r_{\text{in}}$ is reduced further).

Also, as apparent in Fig. 2, there is an excess mass at specific angular momenta larger than $l_{\text{in}}$. This matter with large angular momentum (primarily $l_z$) circularizes within the computational domain (see the cold torus in Fig. 3 of PSB15) and does not cross the inner boundary, unlike the low angular momentum cold gas.

### 3.2.2 Angular momentum pdf of mass crossing $r_{\text{in}}$

Apart from the angular momentum pdf of cold gas mass within the computational domain, it is instructive to consider the pdf of the cold gas crossing the inner boundary, which is more relevant for SMBH accretion estimates. Fig. 3 shows the time-averaged angular momentum pdf ($d\dot{M}/d\log_{10}l$) with respect to $l_x, l_y, l_z$ and $|l|$ of cold gas accreting through the inner boundary of our fiducial run. The pdfs relative to all three components ($l_x, l_y, l_z$) are similar, including the clockwise and counterclockwise components, suggesting that the cold gas crossing the inner boundary is roughly isotropic. The $d\dot{M}/d\log_{10}l$ pdf with respect to $|l|$ is truncated at higher $|l|$, cutting off sharply at $|l| \approx 4 \times 10^{28} \text{ cm}^2 \text{ s}^{-1}$, which corresponds to the specific angular momentum at the inner boundary ($l_{\text{in}}$; see Fig. 2); the gas with larger angular momentum cannot overcome the centrifugal barrier and fall in through the inner boundary. The accretion rate pdfs in Fig. 3 and the mass pdfs in Fig. 2 have similar slopes at low $l$ with respect to $l_x, l_y, l_z$ and $|l|$.

The left-hand panel of Fig. 4 represents the rotational direction ($\cos\theta = l_z/|l|$) of cold gas at all times in the fiducial NFW run. The right-hand panel shows the corresponding rms $\cos\theta = \langle l_z^2/|l|^2 \rangle^{1/2}$ ($\langle \cdot \rangle$ stands for mass-weighted averaging). The blue solid (black dashed) line is the time-averaged, mass-weighted orientation of cold gas within (crossing) the computational domain. The grey-shaded region shows the $1\sigma$ scatter of cold gas within the computational domain while the yellow-shaded region represents the $1\sigma$ scatter of the cold gas crossing the inner boundary. The cold gas with

<footnote n="1">¹ While we have not investigated the origin of these power laws, we note that the power law exponents are analogous to the Maxwellian distribution at low velocities (compared to thermal); namely the 1-D distribution $f \propto v_x^0$ (with respect to the velocity along a fixed axis) is shallower compared to the distribution with respect to the absolute value of the total velocity $f \propto v^2$.</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Time-averaged (between 2 and 3 Gyr) mass distribution of specific angular momentum (left-hand panel is with respect to $l_z$ and the right-hand panel with respect to $|\mathbf{l}|$) in our simulations. To calculate probability density functions (pdfs), equal-sized bins are used in $\log_{10}|\mathbf{l}|$ space. Note that the slopes of the mass distribution for $l_z$ and $|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{s}^{-1}$ in all three simulations are very similar. At smaller angular momenta ($|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{s}^{-1}$), the mass distributions scale roughly as $d\mathcal{M}/dl_z \propto l_z^0$ and $d\mathcal{M}/d|\mathbf{l}| \propto |\mathbf{l}|$. Notice the sudden rise in cold gas mass for specific angular momentum larger than the value corresponding to the inner simulation domain ($l_{in}$; marked by black arrows). This cold gas with large angular momentum corresponds to the massive rotationally supported torus seen in our simulations. Note that in later figures we show averages from 1 Gyr till the end of the simulation to get better statistics of time-averaged quantities (see Table 1). We confirm that the slopes of the pdfs do not vary with the duration of averaging.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Time-averaged (between 1 and 4 Gyr) cold $\dot{M}$ (a mass inflow rate of cold gas across the inner simulation boundary) distribution (probability density function [pdf]) of $l_x$, $l_y$, $l_z$ and $|\mathbf{l}|$ for the fiducial NFW run with $r_{\text{in}} = 1$ kpc. The slopes are similar to the mass distributions in Fig. 2. The mass accretion rate is multiplied by $t_{\text{ff,in}}$, the free-fall time-scale at the inner boundary (15 Myr for the fiducial run). Not only are the $l_x$, $l_y$ and $l_z$ pdfs similar, the clockwise and counterclockwise pdfs are also very similar. There is a slight bias in the $l_z$ distribution towards larger values of $l_z$ in the clockwise sense, which is a signature of the massive rotating torus seen more prominently for the same run in Fig. 2.</figcaption></figure>

$|\mathbf{l}| \lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$ has equal scatter around the value expected for the mean and rms angle between $\mathbf{l}$ and the $z$-axis for isotropic distribution ($\langle \cos \theta \rangle = 0$ and $\langle \cos^2 \theta \rangle^{1/2} = 1/\sqrt{3} \approx 0.58$), implying that the time-averaged angular momentum distribution of the cold gas with small angular momentum within the simulation domain and crossing the inner boundary is almost isotropic. The larger angular momentum gas, on the other hand, has a clockwise bias as also seen in Figs 2 and 3.

### 3.2.3 Time dependence of angular momentum pdf

Now that we have shown that the time-averaged angular momentum distribution of cold gas with $|\mathbf{l}| \lesssim 10^{28} \text{ erg s}^{-1}$ is roughly isotropic, we want to study its time dependence. In particular, stochastic accretion of cold gas is a viable solution of the cooling flow problem only if the low angular momentum cold gas changes direction on a time-scale shorter than the core cooling time (taken to be 0.2 Gyr as a conservative estimate; for a longer cooling time, more cold gas can power SMBH accretion). If the angular momentum variability time-scale is shorter, the cold gas can undergo angular momentum cancellation (and hence accretion) faster than the core cooling time.

Fig. 5 shows the orientation of the average mass-weighted angular momentum distribution of cold gas crossing the inner radial boundary ($\langle l_z/|\mathbf{l}| \rangle$) and rms $\langle l_z^2/|\mathbf{l}|^2 \rangle^{1/2}$, same quantities as in Fig. 4) as a function of time. The orientation of the total cold gas and $|\mathbf{l}| < 10^{29} \text{ cm}^2 \text{ s}^{-1}$ and $|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{ s}^{-1}$ cold gas are shown. The $l_z/|\mathbf{l}|$ ratios of the total cold gas and with $|\mathbf{l}| < 10^{29} \text{ cm}^2 \text{ s}^{-1}$ show that this gas only changes its orientation (as measured by $\cos \theta$ crossing zero) till only 1.5 Gyr, and attains a fixed sense of rotation thereafter. Even though this cold gas (dominated by large angular momentum gas) shows small fluctuations on short time-scales ($<30$ Myr), it still maintains its orientation for the rest of the
<figure><img src="image_7.png" /><figcaption>**Figure 4.** The orientation of cold gas from 1 to 4 Gyr in the fiducial simulation characterized by $\theta$, the angle between the z-axis and the specific angular momentum vector ($\mathbf{l}$). The left panel shows the distribution of $\cos\theta = \langle l_z/|\mathbf{l}| \rangle$ and the right-hand panel shows rms $\cos\theta = (\langle l_z^2/|\mathbf{l}|^2 \rangle)^{1/2}$ ($\langle\rangle$ represents mass-weighted averaging for each $|\mathbf{l}|$ bin). For perfectly isotropic distribution, the mean value of $\cos\theta$ should be 0 and that of rms $\cos\theta$ should be $1/\sqrt{3} = 0.58$ (these values are marked by the horizontal solid red lines). The blue solid line (black dashed line) is the mass-weighted, time-averaged orientation as a function of $|\mathbf{l}|$ for the cold gas within the domain (crossing the inner boundary). The shaded yellow (grey) region represents the $1\sigma$ spread around the mean of cold gas within the computational domain (crossing the inner boundary). Unlike the cold gas within the simulation box, there is no cold gas with large angular momentum crossing the inner boundary (also seen by comparing Figs 2 and 3). Note that the time-averaged angular momentum distribution becomes isotropic for $|\mathbf{l}| \lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Average mass-weighted $\langle l_z/|\mathbf{l}| \rangle$ (bottom panel) and rms $\langle l_z^2/|\mathbf{l}|^2 \rangle^{1/2}$ (top panel) of cold gas crossing the inner boundary as a function of time for the fiducial NFW run. The dashed red lines are for all the cold gas crossing the inner boundary. The purple dot-dashed lines with stars show the average orientation of cold gas with $|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{ s}^{-1}$ and solid blue lines are for $|\mathbf{l}| < 10^{29} \text{ cm}^2 \text{ s}^{-1}$. The orientation of cold gas with $|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{ s}^{-1}$ changes on a short time-scale compared to the other two cases.</figcaption></figure>

simulation time. Time of 1.5 Gyr coincides with the formation of the massive cold torus in the fiducial simulation (see Fig. 3 in PSB15). The rms $l_z^2/|\mathbf{l}|^2$ value for this cold gas remains well above 0.58, suggesting that most of it settles down in a disc-like structure. Note that both $\langle l_z/|\mathbf{l}| \rangle$ and rms $\langle l_z^2/|\mathbf{l}|^2 \rangle^{1/2}$ (see the purple dashed lines with markers in Fig. 5) are much more variable for low angular

momentum cold gas ($|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{ s}^{-1}$) as compared to the total cold gas. The mean angle between $\mathbf{l}$ and the z-axis ($\langle\cos\theta\rangle$) fluctuates around zero and the rms $\langle\cos^2\theta\rangle$ fluctuates around the isotropic value (0.58) on a time-scale shorter than 0.2 Gyr. Therefore, only the low angular momentum gas ($|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{ s}^{-1}$) is expected to participate in feedback heating of the ICM.

Results from our AGN jet-ICM simulations (Fig. 5) show that cold gas accretion at $\lesssim 1$ kpc is stochastic, with the angular momentum of cold gas with $|\mathbf{l}| < 10^{28} \text{ cm}^2 \text{ s}^{-1}$ changing on time-scales shorter than the core cooling time. The mixing and inelastic cloud-cloud collisions of these randomly oriented cold gas clouds near the circularization radius will result in angular momentum cancellation. After losing angular momentum, these clouds will fall ballistically to the centre, where they ought to establish a compact turbulent accretion flow with $H/R \gtrsim 0.1$ (see Hobbs et al. 2011) and viscous time shorter than 0.2 Gyr, the cooling time of the cluster core. From our simulation results, a specific angular momentum $\lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$ can be taken as the limit for SCA (see Figs 4 and 5).

### 3.3 SMBH accretion estimates

In this section we estimate the mass accretion rate on to the central SMBH, first using the standard steady accretion flow model and then using the more plausible stochastic cold gas accretion. For the standard flow model estimate, we assume that cold gas forms a thick extended disc before accreting on to the central SMBH on the local viscous time at the circularization radius. For SCA, the angular momentum of the randomly oriented gas clouds cancel each other and the clouds settle down in a compact accretion flow very close to central SMBH. In this scenario accretion happens at the dynamical time (roughly the time-scale for angular momentum cancellation) because instead of relying on viscosity to transport out angular momentum, accretion is powered by angular momentum cancellation among infalling clouds with random orientations.

Typical efficiency required with respect to the mass accretion rate $\dot{M}_{\text{in}}$ at $\sim 1$ kpc to suppress the cooling rate to at least 10 per cent
of the pure cooling flow value, $\dot{M}_{cf}$, is ~$10^{-4}$. Our cluster mass is $M_{200} = 7 \times 10^{14} M_{\odot}$; this $\epsilon \equiv P_{jet}/\dot{M}_{in}c^2$ is expected to be smaller for lower mass haloes (see Fig. 8 in PSB15), where $P_{jet}$ is the jet power. This means that assuming an SMBH mechanical efficiency of 10 per cent ($P_{jet} = 0.1 M_{BH}c^2$), $\dot{M}_{BH} \gtrsim 10^{-3} \dot{M}_{1kpc}$ is required to sufficiently suppress a cooling flow; i.e. at least 0.1 per cent of the cold gas infalling at 1 kpc must be accreted by the SMBH on a time-scale shorter than the core cooling time. Thus, for our fiducial run (for which $\dot{M}_{in} \approx 20 M_{\odot} \text{yr}^{-1}$; see Table 1) average $\dot{M}_{BH} \gtrsim 0.02 M_{\odot} \text{yr}^{-1}$ is required to suppress a cooling flow.

### 3.3.1 Viscosity-mediated SAF estimate

The time-scale on which the gas crossing the inner boundary reaches the SMBH is determined by the nature of the accretion flow. If the accreting cold gas has zero angular momentum, it falls at the free-fall rate, but in the presence of angular momentum, the time for matter to cross the inner boundary and reach the SMBH is the sum of the local free-fall time at the inner radius ($t_{in}$) and the viscous time-scale at the circularization radius ($t_{visc}$). Using these considerations, we can estimate the viscosity-mediated SAF rate on to the SMBH, using the mass distribution of the cold gas crossing the inner boundary ($d\dot{M}/dl$; here $l = |l|$), as

$$ \dot{M}_{\text{BH,SAF}} \approx \int_{l_{\text{ISCO}}}^{l_{\text{max}}} dl \frac{d\dot{M}}{dl} \frac{t_{\text{ff,in}}}{[t_{\text{ff,in}} + t_{\text{visc}}(l)]} \quad (4) $$

where $l_{max}$ is the maximum angular momentum of the cold gas that can contribute to accretion and $l_{ISCO}$ is the angular momentum corresponding to the innermost stable circular orbit (ISCO), which for a Schwarzschild black hole is, $r_{ISCO} = 6GM_{BH}/c^2 \approx 2.9 \times 10^{-4} M_{BH,9}$ pc. An appropriate upper angular momentum cut-off ($l_{max}$) is the angular momentum for which the viscous accretion time equals 0.2 Gyr, the core cooling time (we show later in this section that $\dot{M}_{BH,SAF}$ is rather insensitive to $l_{max}$). From equation (3), the specific angular momentum for a viscous time of 0.2 Gyr is $0.48 \times 10^{28} \text{cm}^2 \text{s}^{-1}$. Gas with a larger angular momentum accretes on a time-scale longer than the core cooling time and hence cannot stop a cooling flow. Moreover, matter with angular momentum smaller than $l_{ISCO} \equiv \sqrt{12}GM_{BH}/c \approx 1.5 \times 10^{25} M_{BH,9} \text{cm}^2 \text{s}^{-1}$ (assuming a non-rotating BH) falls directly into the SMBH with no way to extract its gravitational potential energy.

The low angular momentum distribution of the gas crossing the inner simulation boundary (with respect to $|l|$; see Fig. 3) can be roughly approximated as $d\dot{M}(l)/dl_{28} \sim 0.56 M_{\odot} \text{yr}^{-1} l_{28}$ ($l_{28} = l/10^{28} \text{cm}^2 \text{s}^{-1}$). Equation (4), with this approximation, can be written as

$$ \dot{M}_{\text{BH,SAF}} \approx 0.56 M_{\odot} \text{yr}^{-1} \int_{l_{\text{ISCO},28}}^{l_{\text{max},28}} dl_{28} \frac{l_{28}}{[1 + t_{\text{visc}}(l_{28})/t_{\text{ff,in}}]} \quad (5) $$

From equation (3), we can write

$$ \frac{t_{\text{visc}}(l_{28})}{t_{\text{ff,in}}} = 120 \left( \frac{t_{\text{ff,in}}}{15\text{Myr}} \right)^{-1} \alpha_{0.1}^{-1} \left( \frac{H}{R} \right)^{-2} M_{\text{BH},9}^{-2} l_{28}^{3} \quad (6) $$

Taking $l_{max,28} = 0.48$, which corresponds to a viscous time equal to the core cooling time (equation 3) for the fiducial parameters,<fnref n="2" />

$$ \dot{M}_{\text{BH,SAF}} = 0.56 \int_{1.5 \times 10^{-3}}^{0.48} dl_{28} \frac{l_{28}}{(1 + 120 l_{28}^3)} \approx 0.018 M_{\odot} \text{yr}^{-1} \quad (7) $$

This value is uncomfortably close to the minimum value required to prevent a cooling flow ($\epsilon \dot{M}_{in}/0.1 \sim 0.02 M_{\odot} \text{yr}^{-1}$). Moreover, we have used a large $H/R = 0.1$ instead of the standard thin AGN disc value of $H/R \sim 10^{-3}$. Using $H/R = 10^{-3}$ gives $l_{max,28} = 0.022$ and $\dot{M}_{\text{BH,SAF}} = 1.3 \times 10^{-4} M_{\odot} \text{yr}^{-1}$, way too small compared to the required accretion rate. The value of the $\alpha$ viscosity parameter is uncertain, with the variability observations giving an order of magnitude larger value ($\gtrsim 0.1$) compared to MHD simulations without net flux (see King, Pringle & Livio 2007, and references therein). A smaller $\alpha$ makes the case for SCA even stronger. We list the actual value of the integral in equation (4) (rather than making a power-law approximation) in Table 1 for $H/R = 0.1$.

### 3.3.2 SCA estimate

Accretion for the stochastic cold gas (Figs 4 and 5 show that $|l| < 10^{28} \text{cm}^2 \text{s}^{-1}$ cold gas is stochastic in space and time over 0.2 Gyr) is expected to proceed at almost the local free-fall time till the gas settles in a thick disc close to the ISCO where the viscous time is very short compared to the free fall at the inner boundary of our simulations, $t_{ff,in}$. From equation (3), the viscous accretion time for $l = 2 \times 10^{27} \text{cm}^2 \text{s}^{-1}$ is 14 Myr. We assume that the stochastic cold gas with initial $l \lesssim 10^{28} \text{cm}^2 \text{s}^{-1}$ cancels its angular momentum to $\lesssim 2 \times 10^{27} \text{cm}^2 \text{s}^{-1}$ on a time-scale shorter than 15 Myr (the free-fall time at the inner boundary). The SCA rate, therefore, can be approximated as

$$ \dot{M}_{\text{BH,SCA}} \approx \int_{l_{\text{ISCO}}}^{l_{\text{max}}} dl \frac{d\dot{M}}{dl} \quad (8) $$

where $l = |l|$. This estimate assumes that most of the mass accreting on to the SMBH originates as cold gas beyond ~1 kpc. Not much gas condenses within 1 kpc because at most times $t_{\text{cool}}/t_{\text{ff}} > 10$ at these radii. Moreover, the available hot gas mass (out of which cold gas can condense) decreases inward. We take $l_{max} = 10^{28} \text{cm}^2 \text{s}^{-1}$ as the upper cut-off of the integral in equation (8) because cold gas is stochastic only for angular momentum lower than this.

The integral in equation (8) can be approximated as $\dot{M}_{\text{BH,SCA}} \approx 0.2 M_{\odot} \text{yr}^{-1} l_{\text{max,28}}^2$ (see Table 1 for the exact value) for our fiducial run for which $d\dot{M}/dl_{28} \approx 0.56 M_{\odot} \text{yr}^{-1} l_{28}$. Note that the SCA estimate is independent of disc parameters ($\alpha$ and $H/R$) as we assume that accretion in the SCA regime happens at the local free-fall time (which is shorter than the free-fall time at the inner boundary). Unlike viscosity-mediated SAF integral (equation 4) which is rather insensitive to $l_{max}$, SCA integral (equation 8) increases quadratically with $l_{max}$. Note that the SMBH mass accretion rate in the SCA regime is comfortably larger than the required $\dot{M}_{\text{BH}} \sim 0.02 M_{\odot} \text{yr}^{-1}$ to quench a cooling flow. In fact, we can also comfortably accommodate the reduction of $\dot{M}_{\text{BH}}$ in the RIAF regime due to outflows.

Fig. 6 shows the jet power, and the mass accretion rates estimated from SAF (equation (4) with $H/R = 0.1$) and SCA (equation 8) as a function of time. As expected, the mass accretion rate is higher on average as compared to viscosity-mediated SAF. Also $\dot{M}_{\text{SCA}}$, in comparison to $\dot{M}_{\text{SAF}}$, is variable on a shorter time. Notice the spikes in $\dot{M}_{\text{SCA}}$ in the cooling phases of the core and the lack of stochastic accretion when jets overheat the core. Similar trends are also seen for $\dot{M}_{\text{in}}$; i.e. all accretion rates are higher in the cooling phases when jet power is small and suppressed in the overheating phases. We discuss the application of SCA estimate based on the angular momentum distribution of cold gas as the basis of a more robust feedback prescription in Section 3.4.

<footnote n="2">² Equation (7) is not sensitive to the upper cut-off as $\dot{M}_{\text{BH,SAF}} \propto l_{28}^{-1}$.</footnote>
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Mass accretion rate estimates as a function of time based on SAF (SAF in black thick solid line; equation 4) and stochastic cold accretion (SCA in green dashed line; equation 8) for our fiducial run. Also shown in magenta solid lines is the expected mass accretion rate on to the SMBH assuming that the SMBH accretion rate is $\epsilon \dot{M}_{\text{in}}/0.1 = 6 \times 10^{-4} \dot{M}_{\text{in}}$. Red dot-dashed line is the jet power (see Section 4.3 on how we calculate the jet/cavity power) as a function of time. All quantities show large variations because of heating/cooling cycles driven by condensation/overheating. The accretion rate at $r_{\text{in}}$ is less time variable compared to $\dot{M}_{\text{SCA}}$ and $\dot{M}_{\text{SAF}}$ but overall trends are similar.</figcaption></figure>

### 3.4 Improving cold feedback prescription

Table 1 lists important time-averaged quantities from our numerical simulations. The accretion rate estimates based on the SAF (equation (4) with $H/R = 0.1$) and SCA (equation (8) with $l_{\text{max}} = 4.8 \times 10^{27} \text{ cm}^2 \text{ s}^{-1}$) are included, as is the average rate of mass crossing the inner boundary. The NFW run with $r_{\text{in}} = 0.5$ kpc but the same $\epsilon$ as the fiducial NFW run (with $r_{\text{in}} = 1$ kpc) shows that the average mass accretion (hot+cold) through the inner boundary is reduced by a factor of $\approx 2$ (11.6 $M_{\odot} \text{ yr}^{-1}$ as opposed to 21.6 $M_{\odot} \text{ yr}^{-1}$ for the fiducial run) but the cold gas mass accumulating within the computational domain (dominated by the high angular momentum gas; e.g. see Fig. 2) by the end of the simulation is about five times higher. In fact, the average cold gas mass accumulation rate for the $r_{\text{in}} = 0.5$ kpc run is $5 \times 10^{11} M_{\odot}/5 \text{ Gyr} \sim 100 M_{\odot} \text{ yr}^{-1}$, only a factor of 2 lower than the pure cooling flow rate. A lot of cold gas crossing 1 kpc circularizes before crossing the inner boundary at 0.5 kpc (compare the angular momentum pdfs and $l_{\text{in}}$ for different runs in the right-hand panel of Fig. 2), and the mass accretion rate at $r_{\text{in}} = 0.5$ kpc ($\dot{M}_{\text{in}}$) is much smaller than the average mass cooling/deposition rate. Importantly, a smaller $\dot{M}_{\text{in}}$ for a fixed $\epsilon$ gives a much lower feedback heating, allowing a large amount of cold gas with sufficient angular momentum to accumulate within the computational domain. The NFW+BCG simulation also uses $r_{\text{in}} = 0.5$ kpc but with a larger feedback efficiency ($\epsilon = 5 \times 10^{-4}$), and therefore the mass of the accumulated cold gas is reasonable. These numbers illustrate the weakness of the current feedback models in which the suitable $\epsilon$ depends, rather sensitively, on the radius at which $\dot{M}_{\text{in}}$ is measured.

Instead of measuring $\dot{M}_{\text{in}}$ at an arbitrary inner radius and tuning $\epsilon$ to get a reasonable match to observations as we are presently doing, our simulations suggest a more reliable way to estimate the SMBH accretion rate based on the mass distribution of low angular momentum cold gas (equation (8) with $|l| \lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$). Table 1 shows that $\dot{M}_{\text{SCA}}$ for the NFW run with $r_{\text{in}} = 0.5$ kpc is about

6.4 times larger than the NFW run with $r_{\text{in}} = 1$ kpc but $\dot{M}_{\text{in}}$ is about a factor of 2 smaller for the former. This implies that $\dot{M}_{\text{BH}}$ estimate based on the angular momentum distribution of cold gas is more robust because a larger mass accretion rate (as obtained from SCA) will prevent excessive cooling and mass deposition seen in the NFW run with $r_{\text{in}} = 0.5$ kpc. Therefore, we anticipate similar outcomes (similar cold gas mass, small scale accretion rate, etc.) for feedback simulations with identical $\epsilon$ as if the SMBH accretion rate is estimated using the angular momentum distribution of the cold gas. In future, we plan to carry out such simulations.

While the results in this section strongly suggest the importance of stochastic cold feedback in maintaining thermal equilibrium in cluster cores, our simulations have several limitations like the absence of star formation, magnetic fields and anisotropic thermal conduction. In the absence of stellar feedback and gas consumption by star formation, we see a build-up of a massive cold gas torus. The next section explores the implications of cold gas depletion due to star formation, via post-processing using a simple model. In the following section, we also compare our cold mode feedback results with Bondi/hot model and with observations.

## 4 COLD VERSUS HOT FEEDBACK: SIMULATIONS CONFRONT OBSERVATIONS

In this section we compare the results of cold feedback simulations with salient observations.

### 4.1 Gas depletion due to star formation

Unlike Li et al. (2015), our numerical simulations do not allow for the depletion of cold gas due to star formation. This leads to a secular build-up of cold gas to unrealistically large values. We can account for the influence of star formation on cold gas depletion using a simple post-processing model. The cold gas mass in the presence of cold gas depletion due to star formation is given by

$$ \frac{dM_{\text{cold}}}{dt} = \dot{M}_{\text{cond}} - \frac{M_{\text{cold}}}{\tau}, \qquad (9) $$

where $M_{\text{cold}}$ is the total cold mass in the computational domain, $\dot{M}_{\text{cond}}$ is the rate of condensation of cold gas from the hot ICM (i.e. the derivative of the $\tau = \infty$ line in Fig. 7) and $\tau$ is the cold gas depletion time-scale (due to star formation). Therefore, the star formation rate $dM_{\text{cond}}/dt = M_{\text{cold}}/\tau$. While in reality, star formation and associated feedback will impact the surrounding gas distribution, we do not take this into account as stellar feedback is subdominant on global scales in cluster cores (although it may help drive local turbulence in the cold gas; e.g. Hobbs et al. 2011). For a constant $\tau$, equation (9) can be solved analytically to give

$$ M_{\text{cold}}(t) = \int_{0}^{t} e^{-(t-t')/\tau} \dot{M}_{\text{cond}}(t')dt'. \qquad (10) $$

Fig. 7 shows the cold gas mass as a function of time for different cold gas mass depletion time-scales ($\tau$). For a gas depletion time-scale $\lesssim 0.5$ Gyr, the cold gas mass is always $\lesssim 5 \times 10^{10} M_{\odot}$ and most of the gas is channeled into stars. For $\tau \lesssim 0.2$ Gyr, there are times when there is no cold gas present in the cluster core. Cold gas mass in our simulations is on the higher side and would lead to an average star formation rate of $\approx 10^{11} M_{\odot}/(5 \text{ Gyr}) = 20 M_{\odot} \text{ yr}^{-1}$. This rate is comparable to that seen in some of the BCGs in cool core clusters (e.g. Bildfell et al. 2008; Mittal, Whelan & Combes 2015; Loubser et al. 2016) but is on the higher side of the observed distribution. This is understandable because in PSB15 we deliberately chose the
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Cold gas mass as a function of time in the fiducial NFW run for different mass depletion time-scales due to star formation. Without depletion of cold gas due to star formation ($\tau = \infty$; this is the same as the cold gas mass as a function of time in Fig. 9 of PSB15), cold gas secularly builds up with time. Depletion of cold gas due to star formation leads to cycles in cold gas mass as a function of time.</figcaption></figure>

accretion efficiency of $\epsilon = 6 \times 10^{-5}$, the minimum value required to bring down the cooling rate by a factor of 0.1 with respect to cooling flow in the absence of AGN feedback. This suggests that the feedback efficiency (with respect to the accretion rate at 1 kpc) in typical clusters is between $10^{-4}$ and $10^{-3}$. Accretion efficiency much above $10^{-3}$ (as in Gaspari et al. 2012; Li & Bryan 2014) leads to strong feedback that wipes out cold gas at radii >10 kpc at late times, and maintains the cluster in the hot state for most of the time.

## 4.2 Cold gas mass, SMBH and halo mass

The mass accretion rate on to the SMBH in hot/Bondi accretion is strongly dependent on the mass of the SMBH (equation 1). In the standard Bondi scenario, one does not normally consider the gas to cool to very low temperatures (e.g. Quataert & Narayan 2000). Cold gas in cluster cores can be naturally understood as a result of local thermal instability in the core satisfying rough global thermal balance. Under the assumption that the cold gas mass in the Bondi scenario is proportional to the SMBH accretion rate, the cold gas mass in cluster cores is expected to be $\propto M_{BH}^2 n T_{keV}^{-3/2}$ (i.e. strongly dependent on SMBH and halo properties). Fig. 8, adapted from McNamara et al. (2011), shows the cold gas mass as a function of the SMBH mass for different groups and clusters. The scatter in cold gas mass is very large, given the rather uniform core entropies observed in cool cores, unlike the expectation from the Bondi accretion rate. Also the observations do not show a strong correlation (anticorrelation) between the SMBH mass (cluster temperature; indicated by the colours of markers) and the mass of molecular gas. In fact, the largest molecular gas mass occurs in the hottest cluster, unlike what is expected from $T_{keV}^{-3/2}$ scaling of $M_B$ (equation 1; this assumes that the temperature of the hot gas at the Bondi radius scales with the core X-ray temperature; X-ray spectra do not show X-ray emitting gas below 0.3 the cluster temperature; e.g. Peterson et al. 2003; see also fig. 1 a of Hogan et al. 2017a). Similarly, in Fig. 8 the data point with the smallest SMBH mass has a large cold gas mass. There is a large scatter in cold gas mass, irrespective of

<figure><img src="image_11.png" /><figcaption>**Figure 8.** Molecular gas mass in centres of cool groups and clusters, obtained using CO observations, as a function of the SMBH mass. The filled circles represent detections and the triangles are upper limits. The cluster temperature is represented by colour (obtained from Cavagnolo et al. 2009). There is a large scatter in the molecular gas mass for a given SMBH mass. Three vertical lines on extreme right denote the range of cold gas mass in our three simulations (measured from 1 Gyr till the end of the run, without accounting for cold gas depletion). A vertical line and an upper limit in extreme left are obtained from the simulations of Li et al. (2015) with their feedback efficiency of 0.01 and 0.001 (bottom panels of their Figs 2 and 7). Our model for cold gas depletion (see equation 9 and Fig. 7) can also reduce the cold gas mass (of course, at the expense of creating more stars).</figcaption></figure>

the SMBH mass and the cluster temperature, a hallmark of cold mode feedback that is stochastic and shows cooling and heating cycles (Li et al. 2015; Prasad et al. 2015).

Fig. 9 in PSB15 and figs 2 and 7 in Li et al. (2015) show that the cold gas mass, for the same SMBH and halo, is highly time variable and shows a broad distribution. In fact, at times, one may not see any cold gas because star formation and stellar feedback can, in principle, consume all the cold gas before it is replenished by ICM cooling (this is seen in fig. 2 of Li et al. 2015, who explicitly model star formation and stellar feedback in their simulations, as well as in our Fig. 7, where we attempt to infer the impact of star formation, but not stellar feedback, on the cold gas mass in our simulation during post-processing). The range of cold gas mass seen in our numerical simulations is indicated by the three vertical intervals on the right-hand side in Fig. 8. Even without star formation and stellar feedback, the lower end of the cold gas mass distribution in our simulations overlaps with the observed range of cold gas mass in cluster cores. Typically, however, the cold gas mass is too high by up to two orders of magnitude. Allowing for gas depletion by star formation on time-scales of 0.2–0.5 Gyr during post-processing (Fig. 7) brings our gas mass into agreement with the observations. Nonetheless, this is an issue that requires further study, and we are in the process of explicitly including a realistic treatment of star formation and stellar feedback, as well as allowing for the effects of quasar-mode feedback that may occur when the AGN occasionally switches to that mode, in our simulations.

## 4.3 Bondi power and cavity power

Allen et al. (2006) noted a tight correlation between the Bondi power and the cavity power for nearby elliptical galaxies and argued
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Bondi power (estimated using two different extrapolations of the spherically averaged hot gas [$T > 0.1$ keV] density profiles to the Bondi radius) as a function of cavity power from our fiducial NFW run. Simulation data are sampled every 10 Myr and there is a large scatter for both extrapolations. Also shown as green triangles is the instantaneous accretion power ($\epsilon \dot{M}_{in}c^2$) supplied by the SMBH. At most times Bondi power is within $10^{43}-10^{45}$ erg s$^{-1}$, but cavity power has a much larger variation. A similar behaviour is seen in our other runs listed in Table 1. Data points from Russell et al. (2013) and the best-fit from Allen et al. (2006) are also shown.</figcaption></figure>

for the validity of Bondi accretion. Russell et al. (2013) recently refined the determination of Bondi and cavities powers from Allen et al.'s sample and included more data points, and found that the correlation is much weaker. Fig. 9 shows the Bondi-power-cavity-power relationship obtained for our fiducial NFW run where, as we have discussed, the jet is in fact triggered by cold gas condensation, not Bondi accretion. Bondi power $P_B = 0.1 \dot{M}_{BH}c^2$, where $\dot{M}_B$ is given by equation (1), is calculated using the hot gas density profile in our domain, extrapolated to the Bondi radius,

$$R_B \equiv 2GM_{BH}/c_{s,\infty}^2 \approx 35 \text{ pc } M_{BH,9}T_{\infty,\text{keV}}^{-1}. \quad (11)$$

Emulating what Russell et al. (2013) do with the observational data, we also use two different best-fitting models for the hot gas density applied from $r_{in}$ to 10 kpc: a power-law [PL] model (red squares in Fig. 9); and a $\beta$ model (blue stars) with the hot gas density going like $n_0[1 + (r/r_0)^2]^{-3\beta/2}$. A power-law profile gives a larger Bondi accretion rate (see equation 1) because the density extrapolated to the Bondi radius is higher. The temperature at the Bondi radius (needed to calculate $c_{s,\infty}$ in equation 11) is assumed to be the same as that at $r_{in}$; results are unchanged with other reasonable extrapolations of the temperature. Also shown in Fig. 9 is the instantaneous value of the feedback power, where $\epsilon \dot{M}_{in}c^2$ (green triangles), $\dot{M}_{in}$ is the total [hot+cold] accretion rate measured at the inner radius $r_{in}$.

**Estimating cavity power:** The cavity power in Fig. 9 is calculated in the usual way. To estimate the jet/cavity energy, we compute $4\int pdV$ over the volume of low-density regions (density lower than 0.2 times the initial minimum density; results are insensitive to the exact value of the threshold as long as it is smaller than the initial minimum density), corresponding to the enthalpy of a relativistic bubble. To convert this into power, we divide the energy by an estimate of the jet sound crossing time using the outermost radius of the bubble material ($r/[600 \text{ km s}^{-1}]$; $600 \text{ km s}^{-1}$ is the sound

speed for the 1.5 keV cluster core plasma; see the bottom-right panel in fig. 2 of PSB15). This definition is slightly different from what was used in PSB15, but the results are independent of the details of jet power and Bondi estimates. As discussed in PSB15, there is a large scatter between cavity power and the instantaneous accretion power because of hysteresis in cold mode feedback.

More importantly, the scatter and the lack of obvious correlation in the observational data from Russell et al. (2013) in Fig. 9 are similar to those for the Bondi power-cavity power scatter seen in our fiducial simulation. Notice that the variation in Bondi power is smaller than the variation in cavity power in the simulations. This is because the hot gas properties in the core (on which the Bondi estimate depends) do not vary as much in time as the cold gas accretion rate. Moreover, cooling and accretion events occur closely spaced in time and our feedback prescription can blow really powerful cavities (see jet power in Fig. 6). The observational data point in Fig. 9 with the highest Bondi power is for M87, for which Russell et al. (2013) use an SMBH mass of $6.6 \times 10^9 M_{\odot}$, much larger than our fiducial choice of $10^9 M_{\odot}$; recalling that $\dot{M}_B \propto M_{BH}^2$ brings this point in agreement with the simulation results. The highest cavity power scales with the halo mass and our halo mass is on the higher side. Also, unlike simulations, observations can miss out fainter parts of cavities, leading to an underestimate of the cavity power.

## 5 DISCUSSION

The problems with Bondi accretion (hot accretion in general) are well recognized (see Section 1). Standard viscous mediated accretion flows on to SMBHs, even with a large $H/R \sim 0.1$, leads to the angular momentum problem; namely, the viscous accretion time-scale for the accretion flow is much longer than the cooling time of the cluster core. Therefore, enough gas cannot accrete fast enough and provide sufficient feedback power to prevent a strong cooling flow.

A number of published studies have invoked the stochastic nature of the condensing cold gas and its fundamentally different behaviour compared to the SAF to overcome the angular momentum problem (Pizzolato & Soker 2010; Hobbs et al. 2011; Nayakshin et al. 2012; Gaspari et al. 2013). Through 3-D hydrodynamic simulations, we explicitly demonstrate that the angular momentum of cold gas condensing out of the hot phase indeed varies stochastically and has an isotropic distribution (especially for the low specific angular momentum, $l \lesssim 10^{28}$ cm$^2$ s$^{-1}$, cold gas). Cancellation of almost randomly directed angular momentum on short time-scales leads to the formation of a compact accretion flow that accretes on to the SMBH on the much shorter dynamical time-scale (compare equations 4 and 8). Recent ALMA observations of infalling cold, clumpy CO clouds in the elliptical galaxy at the centre of Abell 2597 (Tremblay et al. 2016), within 100 pc of the SMBH, provides an almost direct evidence for SCA. Nayakshin et al. (2012) go on to argue that even in the higher $\dot{M}$ regime, SCA is responsible for the emergence of the SMBH $M_{BH} - \sigma$ relation because accretion via a standard extended disc results in more massive SMBHs than observed since a massive extended thin disc does not couple efficiently to an isotropic feedback-powered outflow. Turbulence due to stochastic angular momentum and infall of cold gas (viscous time is much shorter for a turbulent flow with $H/R \sim 0.1$) also prevents the problem of gravitational fragmentation in stochastic accretion flow (Goodman 2003).

Most of the earlier studies of SCA were carried out in numerical setups with turbulence in the cold gas driven artificially (to
realize the stochastic distribution of cold gas angular momentum). However, in this paper we have shown that the turbulence emerges self-consistently from jet-ICM interactions and this behaviour is sustained over time-scales of several Gyr (i.e. much longer than the cluster core cooling time). We analyse the angular momentum distribution of cold gas in our AGN jet-ICM simulations (which appear consistent with most observations of cool cluster cores; see Section 4) and demonstrate that SCA is indeed realized in cool cluster cores stirred by time-dependent AGN jets.

In cluster simulations with limited resolution, the accretion physics input necessarily has to be inserted via a sub-grid model. One such model, which is much better than what is typically used (including ours in this paper where we estimate $\dot{M}$ at $\sim$1 kpc and multiply it by an efficiency parameter $\epsilon$ to calculate the mechanical power injected by the SMBH), is to use the mass accretion rate of the low angular momentum gas (with say $l \lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$; see equation 8 and Fig. 6) crossing the smallest resolved radius (see Section 3.4 for details). Another improvement is the inclusion of star formation and realistic stellar feedback (stellar feedback is increasingly important for lower mass haloes such as clusters and elliptical galaxies; e.g. see fig. 4 in Sharma et al. 2012b) that can prevent artificial accumulation of massive cold gas in cluster cores (e.g. as done by Li et al. 2015). Also, as mentioned earlier, our choice of $\epsilon$ (feedback efficiency parameter) is on the lower side, resulting in somewhat large mass deposition and star formation in the core. More 3-D simulations with larger $\epsilon$s are required to get a complete picture.

## 5.1 Uncertainties in accretion physics

Much more needs to be done in understanding the details of SCA. In our simulations, which focus on larger ($\gtrsim$1 kpc) scales, we have to make assumptions about the underlying stochastic accretion. In particular, we need to understand how the stochastic cold gas condensing out of the ICM is eventually accreted on to the SMBH. Whether a disc forms after the cancellation of stochastic angular momentum and on what time-scale, and at what radii does the transition from SCA to a viscous thin disc (if it does at all) and to an RIAF occurs. While the transition radius from a thin disc to an inner hot/thick RIAF occurs at a smaller radius for a higher accretion rate, the details are not clear even for the standard (non-stochastic) accretion flows with a fixed angular momentum axis (e.g. see the vastly different estimates of disc to RIAF transition radius in Das & Sharma 2013 [Fig. 9b], Liu et al. 1999 [Fig. 1], Yuan & Narayan 2014 [Fig. 7b]).

Strong AGN jets and X-ray cavities (which are powered by jets) are expected only in the radiatively inefficient accretion regime for which $\dot{M}_{BH} \lesssim \alpha^2 \dot{M}_{Edd}$ (Rees et al. 1982; Churazov et al. 2005; Nemmen et al. 2007). Our SMBH accretion rate at most times is lower than this threshold for jet formation (see Fig. 6). Infrequent, short-duration high $\dot{M}_{BH}$ events, in fact, provide a potential basis for short-lived transition from radio mode to quasar mode and for changing the direction of the jet, as observed in numerous clusters (e.g. Babul et al. 2013). Generally the mass accretion rate on to the SMBH is expected to be even lower because of outflows in RIAFs (larger suppression is expected for bigger RIAFs). Therefore, advances in accretion physics are required to estimate the size of RIAFs and the mass accretion rate on to SMBHs fed by SCA.

To sum up, the majority of the cold gas condensing out of the ICM due to thermal instability has large angular momentum such that the viscous accretion time-scale is longer than the core cool-

ing time for $l \gtrsim 10^{28} \text{ cm}^2 \text{ s}^{-1}$ (equation 3). Most of this gas is expected to be Toomre unstable and to form stars, instead of being accreted on to the SMBH. The lower angular momentum cold gas ($l \lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$) has stochastic angular momentum, which implies that cloud collisions (which lead to rapid angular momentum cancellation) will cause the cold, turbulent gas to be accreted rapidly, without much gravitational fragmentation. At innermost radii typically an RIAF with powerful jets is expected, since $\dot{M}_{BH} \ll 0.01 \dot{M}_{Edd}$.

## 6 SUMMARY

We finally summarize our major findings in the following points:

(i) SCA, which was anticipated to play a key role in accretion on to SMBHs from idealized simulations, is realized naturally in our realistic AGN jet-ICM simulations due to the turbulence induced by jets in a non-uniform ICM over several Gyr time-scale.

(ii) We find that the low angular momentum cold gas ($|l| \lesssim 10^{28} \text{ cm}^2 \text{ s}^{-1}$) condensing out of the ICM has an isotropic distribution of angular momentum, which we expect will result in the cancellation of angular momentum on almost a dynamical time. This is unlike a viscosity-mediated SAF, in which an extended thin disc with a long accretion time is expected to form. Another advantage of a stochastic accretion flow is that the net angular momentum of the low angular momentum cold gas varies on a time-scale shorter than the core cooling time, implying that angular momentum cancellation and feedback driven by cold gas accretion can respond fast enough. Furthermore, SCA is independent of disc parameters like $\alpha$ and $H/R$. Elimination of these parameter (usually with a large range) makes SCA a simpler and more robust model.

(iii) Our work suggests an improved feedback prescription based on the angular momentum distribution of cold gas (Section 3.4), which may be more robust than the usual models in which the mass accretion rate is estimated at $\sim$ kpc scales. In the latter approach, one has to fine-tune the feedback efficiency parameter ($\epsilon$) for different radii at which the mass accretion rate is estimated (compare our NFW runs with the inner radius of 1 and 0.5 kpc).

(iv) Most of the features of cold accretion simulations, which show multiphase cooling and heating cycles driven by accretion of cold gas, match observations such as a large scatter of cold gas mass when compared to the SMBH and halo masses (this is not expected for hot/Bondi accretion; see Section 4.2). A large scatter observed between Bondi power and jet/cavity power is also consistent with our cold feedback simulations (see Section 4.3).

## ACKNOWLEDGEMENTS

This work is partly supported by the DST-India grant no. Sr/S2/HEP-048/2012 and an India-Israel joint research grant (6-10/2014[IC]). DP is supported by a CSIR grant (09/079[2599]/2013-EMR-I). AB acknowledges support from NSERC, through the Discovery Grant Program, the Institut Lagrage de Paris, and Pauli centre for Theoretical Studies ETH UZH. He also thanks the Institut d'Astrophysique de Paris (IAP), at the Institute for Computational Sciences and University of Zurich for hosting him. We acknowledge the support of the Supercomputing Education and Research Centre (SERC) at IISc for facilitating our use of Cray XC40-SahrasT cluster on which some of our runs were carried out. We thank the anonymous referee for critical comments that helped improve our paper. PS thanks Sagar Chakraborty for useful discussions.
REFERENCES

Allen S. W., Dunn R. J. H., Fabian A. C., Taylor G. B., Reynolds C. S., 2006, MNRAS, 372, 21

Babul A., Sharma P., Reynolds C. S., 2013, ApJ, 768, 11

Baganoff F. K., et al., 2003, ApJ, 591, 891

Bildfell C., Hoekstra H., Babul A., Mahdavi A., 2008, MNRAS, 889, 1637

Bondi H., 1952, MNRAS, 112, 195

Cavagnolo K. W., Donahue M., Voit G. M., Sun M., 2009, ApJS, 182, 12

Choudhury P. P., Sharma P., 2016, MNRAS, 457, 2554

Churazov E., Sazonov S., Sunyaev R., Forman W., Jones C., Böhringer F., 2005, MNRAS, 363, L91

Das U., Sharma P., 2013, MNRAS, 435, 2431

Fabian A. C., 2012, ARA&A, 50, 455

Frank J., King A., Raine D., 2002, Accretion Power in Astrophysics, 3rd edn. Cambridge Univ. Press, Cambridge

Gaspari M., Ruszkowski M., Sharma P., 2012, ApJ, 746, 94

Gaspari M., Ruszkowski M., Oh S. P., 2013, MNRAS, 432, 3401

Gaspari M., Brighenti F., Temi P., 2015, A&A, 579, A62

Gaspari M., Temi P., Brighenti F., 2017, MNRAS, 466, 677

Goodman J., 2003, MNRAS, 339, 937

Hawley J. F., Krolik J. H., 2001, ApJ, 548, 348

Hayes J. C. et al., 2006, ApJS, 165, 188

Hobbs A., Nayakshin S., Power C., King A., 2011, MNRAS, 413, 2633

Hogan M. T., McNamara B. R., Pulido F., Nulsen P. E. J., Russell H. R., Vantyghem A. N., Edge A. C., Main R. A., 2017a, ApJ, 837, 51

Hogan M. T. et al., 2017b, preprint (arXiv:1704.00011)

Hopkins P. F., Hernquist L., 2006, ApJ, 166, 1

King A. R., Pringle J. E., Livio M., 2007, MNRAS, 376, 1740

Li Y., Bryan G. L., 2014, ApJ, 789, 54

Li Y., Bryan G., Ruszkowski M., Voit G. M., O'Shea B. W., Donahue M., 2015, ApJ, 811, 73

Liu B. F., Yuan W., Meyer F., Meyer-Hofmeister E., Xie G. Z., 1999, ApJ, L17

Loubser S. I., Babul A., Hoekstra H., Mahdavi A., Donahue M., Bildfell C., Voit G. M., 2016, MNRAS, 456, 1565

Macchetto F., Marconi A., Axon D. J., Capetti A., Sparks W., Crane P., 1997, ApJ, 489, 579

Marrone D. P., Moran J. M., Zhao J., Rao R., 2007, ApJ, 654, L57

McCourt M., Sharma P., Quahaert E., Parrish I. J., 2012, 419, 3319

McDonald M., Veilleux S., Rupee D. S. N., Mushotzky R., 2010, ApJ, 721, 1262

McNamara B. R., Nulsen P. E. J., 2007, ARA&A, 45, 117

McNamara B. R., Nulsen P. E. J., 2012, New J. Phys., 14, 055023

McNamara B. R., Rohanizadegan M., Nulsen P. E. J., 2011, ApJ, 727, 39

Mignone A., Bodo G., Massaglia S., Matsakos T., Tesileanu O., Zanni C., Ferrari A., 2007, ApJS, 170, 228

Mittal R., Whelan J. T., Combes F., 2016, 450, 2564

Monin A. S., Yaglom A. M., 1971, Statistical Fluid Mechanics, Dover Edition Dover Press, New York

Nayakshin S., King A., 2007, MNRAS, preprint (arXiv:0705.1686)

Nayakshin S., Power P., King A., 2012, ApJ, 753, 15

Nemmen R. S., Tchekhovskoy A., 2015, MNRAS, 449, 316

Nemmen R. S., Bower R. G., Babul A., Storchi-Bergmann T., 2007, MNRAS, 377, 1652

Peterson J. R., Kahn S. M., Faerels F. B. S., Kaastra J. S., Tamura T., Bleeker J. A. M., Ferrigno C., Jernigan J. G., 2003, ApJ, 590, 207

Pizzolato F., Soker N., 2005, ApJ, 632, 821

Pizzolato F., Soker N., 2010, MNRAS, 408, 961

Prasad D., Sharma P., Babul A., 2015, ApJ, 811, 108 (PSB15)

Quataert E., Narayan R., 2000, ApJ, 528, 236

Rees M. J., Begelman M. C., Blandford R. D., Phinney E. S., 1982, Nature, 295, 17

Russell H. R., McNamara B. R., Edge A. C., Hogan M. T., Main R. A., Vantyghem A. N., 2013, MNRAS, 432, 530

Russell H. R., Fabian A. C., McNamara B. R., Broderick A. E., 2015, MNRAS, 451, 588

Shakura N. I., Sunyaev R. A., 1973, A&A, 24, 337

Sharma P., McCourt M., Quataert E., Parrish I. J., 2012a, MNRAS, 420, 3174

Sharma P., McCourt M., Parrish I. J., Quataert E., 2012b, MNRAS, 427, 1219

Shu F. H., 1992, The Physics of Astrophysics, Vol. II, Gas Dynamics. Univ. Science Books, Sausalito

Soker N., 2016, New Astron. Rev., 75, 1

Soker N., Sternberg A., Pizzolato F., 2009, in Heinz S., Wilcots E., eds, AIP Conf. Ser. Vol. 1201, The Moderate Cooling Flow Model and Feedback in Galaxy Formation. Am. Inst. Phys., Melville, NY, p. 321

Tremblay G. R. et al., 2016, Nature, 534, 218

Voit G. M., Donahue M., Bryan G. L., McDonald M., 2015, Nature, 519, 203

Wong K.-W., Irwin J. A., Yukita M., Million E. T., Mathews W. G., Bregman J. N., 2011, ApJ, 736, L23

Wong K.-W., Irwin J. A., Shcherbakov R. V., Yukita M., Million E. T., Bregman J. N., 2014, ApJ, 780, 9

Yang H.-Y. K., Reynolds C. S., 2016, ApJ, 829, 90

Yuan F., Narayan R., 2014, ARA&A, 52, 529

This paper has been typeset from a TEX/LATEX file prepared by the author.