Galaxy Motions, Turbulence and Conduction in Clusters of Galaxies

Galaxy Motions, Turbulence and Conduction in Clusters of Galaxies


Unopposed radiative cooling in clusters of galaxies results in excessive mass deposition rates onto the central brightest cluster galaxy. However, the cool cores of galaxy clusters are continuously heated by thermal conduction and turbulent heat diffusion due to minor mergers or the galaxies orbiting the cluster center. These process can either reduce the energy requirements for AGN heating of cool cores, or they can prevent overcooling altogether. We perform 3D MHD simulations including field-aligned thermal conduction and self-gravitating particles to model this in detail. Turbulence is not confined to the wakes of galaxies but is instead volume-filling, due to the excitation of large-scale g-modes. We systematically probe the parameter space of galaxy masses and numbers to assess when the cooling catastrophe is prevented. For a wide range of observationally motivated galaxy parameters, we find that the magnetic field is randomized by stirring motions, restoring the conductive heat flow to the core. The cooling catastrophe either does not occur or it is sufficiently delayed to allow the cluster to experience a major merger that could reset the conditions in the intracluster medium. Whilst dissipation of turbulent motions (and hence dynamical friction heating) is negligible as a heat source, turbulent heat diffusion is extremely important; it predominates in the cluster center. However, thermal conduction becomes important at larger radii, and simulations without thermal conduction suffer a cooling catastrophe. Conduction is important both as a heat source and to reduce stabilizing buoyancy forces, enabling more efficient diffusion. Turbulence enables conduction, and conduction enables turbulence. In these simulations, the gas vorticity—which is a good indicator of trapped g-modes–increases with time. The vorticity growth is approximately mirrored by the growth of the magnetic field, which is amplified by turbulence.

conduction – cooling flows – galaxies: clusters: general – galaxies: active – instabilities – X-rays: galaxies: clusters

1 Introduction

The intracluster medium (ICM) in many galaxy clusters has central cooling times shorter then the Hubble time. Radiative cooling should lead to large accumulation of cold material in their centers; however, there is no observational evidence for such gas. This can be understood if some source of heating balances cooling in the ICM. The heating mechanisms invoked to explain this overcooling problem involve AGN “radio mode” heating (e.g., Binney & Tabor (1995); Churazov et al. (2002); Fabian et al. (2003); Ruszkowski et al. (2004a, b); Scannapieco & Brüggen (2008)), preheating by AGN (McCarthy et al., 2008), cosmic rays from AGN (Guo & Oh, 2008; Sharma et al., 2009a), supernovae, turbulent mixing (Kim & Narayan, 2003a; Voigt & Fabian, 2004; Dennis & Chandran, 2005), thermal conduction (Zakamska & Narayan, 2003; Kim & Narayan, 2003b), a combination of thermal conduction and AGN (Ruszkowski & Begelman, 2002) and dynamical friction (El-Zant et al., 2004; Kim et al., 2005; Kim, 2007); see Conroy & Ostriker (2008) and McNamara & Nulsen (2007) and references therein for reviews of the above mechanisms).
Conduction alone is unlikely to offer the complete solution to the overcooling problem for the full range of cluster masses, as its strong temperature dependence implies that it is less effective in lower mass clusters. Furthermore, thermal conduction is well known to be an unstable heating mechanism, either failing to avert a cooling catastrophe, or leading to an isothermal temperature profile (Bregman & David, 1988; Guo & Oh, 2008; Conroy & Ostriker, 2008). Nevertheless, thermal conduction may entirely suppress cooling in non cool-core (NCC) clusters and reduce the constraints on the required energy injection by AGN in cool-core (CC) clusters (Guo et al., 2008); indeed, it is different to stabilize massive clusters with AGN feedback alone (Conroy & Ostriker, 2008), and a second heat source (such as conduction) is generally required. Besides offsetting radiative losses and stemming a cooling catastrophe, conduction can have important implications for establishing the observed bimodality in cluster core entropy (Guo et al., 2008), and the star formation threshold in brightest cluster galaxies (Voit et al., 2008). Indeed, a sudden increase in conduction (due to say, turbulence from an AGN outburst or a merger) could mediate a cool core to non cool-core transition (Guo & Oh, 2009).

Thermal conduction can be strongly suppressed by magnetic fields that are known to be present in the ICM (e.g., Enßlin et al. (2003), Vogt et al. 2003). However, interest in thermal conduction as a potential heating mechanism was revived by Narayan & Medvedev (2001) who suggested that even in the presence of tangled B-fields, the level of conduction can be an appreciable fraction of the Spitzer-Braginskii value. Computing the exact magnitude and distribution of the effective conductivity of the ICM is further complicated by buoyancy instabilities which re-orient the magnetic field. When temperature increases in the direction of gravity, as in the cluster outskirts, the magnetothermal instability (MTI; Balbus (2000); Parrish & Stone (2005)) tends to make the B-fields radial, thereby increasing the effective conduction. On the other hand, in cool cores where temperature decreases in the direction of gravity, the heat flux buoyancy instability (HBI; Quataert (2008); Parrish et al. (2009); Bogdanović et al. (2009)) tends to reorient the fields in the direction perpendicular to that of gravity, effectively shutting down thermal conduction.

However, these instabilities do not operate in a static atmosphere. Chandra and XMM observations show that the cluster gas is rarely in perfect hydrostatic equilibrium. Sloshing motions due to minor mergers, AGN, or galaxy motions can continuously and significantly perturb the gas, as has been repeatedly seen in many disparate numerical simulations (Evrard, 1990; Norman & Bryan, 1999; Nagai et al., 2003; Vazza et al., 2009, 2010; Ascasibar & Markevitch, 2006; ZuHone et al., 2010). Current observational evidence for turbulence ranges from the analysis of pressure maps (Schuecker et al., 2004), its effect on resonant-line scattering (Churazov et al., 2004), and Faraday rotation maps (Vogt & Enßlin, 2005; Enßlin & Vogt, 2006), as well as constraints on turbulent line widths (Sanders et al., 2009). Low levels of turbulence in the ICM can randomize the field configuration set up by the HBI and restore the heat flow to the core (Ruszkowski & Oh, 2010; Parrish et al., 2010). Both of these works modeled the ICM turbulence via a simple driving mechanism to determine the level of turbulence required to effectively restore thermal conduction. This approach did not allow us to link the level of turbulence to the physical properties of the cluster (such as mechanical luminosity of the central AGN or the properties of cluster galaxies). Also, the driving mechanism led, by construction, to volume-filling turbulence which was very effective in randomizing the magnetic field. While the low amplitude of the required subsonic turbulence is eminently feasible (), the realism of volume-filling turbulence is less clear. For instance, both analytic calculations (Subramanian et al., 2006) and numerical simulations (Iapichino & Niemeyer, 2008) predict that turbulence due to galaxy wakes should not be volume-filling (), as turbulence is largely confined to ’streaks’ behind orbiting galaxies.

In this paper, we extend our previous work and perform three-dimensional MHD simulations of the effect of turbulence driven by galaxy motions on the properties of the anisotropic thermal conduction. We show how the trapping of gravity modes excited by the orbiting galaxies can lead to volume-filling turbulence of the right magnitude to restore conductive heat flow. We demonstrate how these subsonic motions generate vorticity and lead to the growth of magnetic field via kinematic dynamo action. We also show that turbulent heat diffusion is an important part of the energy budget.

The paper is organized as follows. In §2 we review basic theoretical expectations for the interaction between turbulence and magnetic fields. In §3 we describe the numerical methods and the setup of the inital conditions. In §4, we describe our results, including the level and volume-filling nature of turbulence, evolution of the gas temperature, generation of vorticity and magnetic fields, and nature of heating mechanisms. Conclusions are summarized in §5.

2 Theoretical Expectations: Turbulence and Trapping of g-modes

It is useful to begin by reviewing some basic theoretical expectations for the behavior of turbulence excited in galactic wakes in clusters. In principle, orbiting galaxies can excite galactic wakes by two means: hydrodynamically (as the ICM collides with the ISM of the galaxy), and gravitationally (similar to dynamical friction for collisionless particles). In practice, we shall conservatively assume that ram pressure stripping is efficient in removing the ISM of galaxies and thus that galaxies only exert gravitational influence.

Volume-filling Turbulence The volume filling factor of galaxy wakes is small (for a simple analytic estimate, see Subramanian et al. (2006)). This might seem to imply that the impact of turbulence excited by galactic wakes is confined to a small fraction of the cluster. However, orbiting galaxies can also resonantly excite g-modes, which from a formal WKBJ analysis have the dispersion relation:


. The Brunt-Väisälä  frequency for buoyant oscillations is:


where is the fluid entropy , and applies if thermal conduction is negligible, while applies if the thermal conduction time is sufficiently short that a displaced blob’s temperature is determined by conductive rather than adiabatic cooling (Sharma et al., 2009b). Note that depend on the entropy and temperature gradient respectively; typically, is about a factor of 2 smaller than .

The above dispersion relation immediately implies that to obtain modes with real , the driving frequency ; otherwise the modes have imaginary radial wavenumber and are evanescent. Physically, one can always achieve a low frequency response by making the mode progressively more tangential, but it is impossible to drive the system at frequencies higher than the maximum response frequency of , corresponding to completely vertical oscillations. This thus implies that waves driven at frequencies can be resonantly excited, and must propagate inward toward the cluster center (as can be seen from their group velocity; Balbus & Soker (1990)), where they will be trapped, reflected and focused inside the resonance radius where . A linear analysis by Balbus & Soker (1990) showed that most of the power in -modes is in the longest wavelengths1. Note that both (which depends on the orbital frequencies of galaxies) and are sensitive to the gravitational potential, which is instrumental in determining if g-modes will be excited.

Isotropic Turbulence Turbulence in the fluid has to compete with buoyancy forces arising from stable stratification. One can show that the ratio of tangential and radial velocities is given by (e.g., see discussion in §2 of Ruszkowski & Oh (2010)):


where is the eddy turnover frequency at a given scale, and Fr is the Froude number, which compares inertial and gravitational forces ( is the Richardson number). If , then turbulence is fundamentally 2D, and for instance it is difficult to rearrange magnetic fields in the radial direction. However, the level of turbulence required to overcome stable stratification is weak; for typical cluster conditions the critical turbulent velocity is (Sharma et al., 2009b):


where is the gravitational acceleration in units of , is a characteristic scale height in units of kpc, and is the critical Richardson number; is typical for hydrodynamic flow.

At first blush, the requirement for might seem to be at odds with the requirement that for modes to be excited. However, note that for homogeneous Kolmogorov turbulence, ; it is therefore conceivable that low frequency g-modes can be excited on large scales, while high-frequency small-scale modes can overcome stabilizing buoyancy forces. Since our background state is not homogeneous, we have to resort to 3D simulations to verify if this expectation is indeed satisfied. This is a major goal of this paper.

Vorticity and B-field growth G-modes excite vorticity. An easy way to see this is to examine the vorticity evolution form of the momentum equation for -waves (i.e., assuming ) (Lufkin et al., 1995):


where is the vorticity, and to note that is in general non-radial, so that is non-zero (indeed, we see in Fig. 2 that since rises toward the center, that g-modes become progressively more tangentially biased there). This implies that vorticity is a good tracer of g-modes, a fact which we shall exploit. It also means that g-modes could conceivably drive an efficient dynamo. There is a well-known analogy between the vorticity equation:


where is the viscosity, and the relation for the magnetic field in the flux-freezing limit,


where is the electricity resistivity. This, together with the fact that the divergence of and both vanish, leads to the expectation that their growth might be related2. There have been a number of studies pointing out that turbulent motions could give rise to magnetic fields in clusters (e.g., Ruzmaikin et al. (1989); Subramanian et al. (2006); Ryu et al. (2008); Cho et al. (2009)). This subject is rich and beyond the scope of this paper; we shall merely compare the growth of vorticity and magnetic fields in our simulations, to see how well they track one another. A reasonable expectation is that the magnetic fields achieve equipartition with turbulence (e.g., Schekochihin & Cowley (2007), and references therein):


where is the rms turbulent velocity on large scales. The above estimate is consistent with observed G fields (Carilli & Taylor, 2002), though there are considerable uncertainties3. The fact that trapping of -modes can give rise to volume-filling turbulence would then be instrumental in allowing volume-filling magnetic fields.

Magnetic tension Magnetic tension can inhibit the HBI (Quataert, 2008). For perturbation scales comparable to the radius (i.e., ) we obtain a critical value:


for suppression, where the fiducial values are measured from our simulated cluster at a radius of kpc. Due to the similarity between the field values in equation (8) and (9), it has been suggested that magnetic fields amplified by turbulence can prevent the onset of the HBI (e.g., see discussion in Kunz et al. (2010). Note that their version of equation (9) yields somewhat lower B-fields than ours, for identical parameters. In any case, equation (9) is only an approximation as the derivation assumes the WKB approximation, while the non-linear saturation of the HBI occurs on global scales). In this paper, we will deliberately ignore this possibility. Observationally, the strength of the magnetic field in the ICM is G and has a large scatter of about an order of magnitude within the ICM and between clusters (Carilli & Taylor, 2002); moreover, there are considerable observational uncertainties in these values, as mentioned above. Numerical simulations show that the HBI still develops for G fields (I. Parrish, priv. comm.), although it can be delayed for increased field strengths. Given the large uncertainty in whether observed field strengths are capable of stabilizing the HBI, past studies of HBI (e.g., Parrish et al. (2009); Bogdanović et al. (2009); Ruszkowski & Oh (2010); Parrish et al. (2010)) focused on the regime where the magnetic tension is unimportant. We also adopt the same approach here, and study if volume-filling turbulence alone can stabilize the HBI. More specifically, we consider plasma and note that, as long as the field is not dynamically important, its exact value does not play a role. In this case, the magnetic field strength scales out of the problem and only serves as a medium to redirect the heat flow via anisotropic thermal conduction. We can therefore study the effects of turbulence alone without the possibly confounding effects of magnetic tension.

Turbulent heating and heat diffusion Turbulence impacts the thermodynamics of the fluid through its effect on thermal conduction, both randomizing and amplifying the magnetic field. Both of these suppress the HBI, and allow thermal conduction at the Spitzer value. However, turbulence can also directly affect the thermal state of the plasma through dissipation of turbulent motions (direct heating), or allowing heat transport via turbulent diffusion (Kim & Narayan (2003b); Dennis & Chandran (2005), and references therein). The heating rate from dissipation of kinetic and magnetic energy is:


where is a dimensionless constant of order unity and , the dominant velocity length-scale, is unknown but almost certainly a function of radius; a reasonable ansatz might be (Dennis & Chandran, 2005), where is some adjustable constant of order unity, and is some minimal lengthscale. On the other hand, the heating rate from turbulent heat diffusion is:


where is the specific entropy, and the turbulent diffusivity is:


where the second factor of takes into account the damping of radial heat transport by buoyancy forces (Dennis & Chandran, 2005). The fact that is of order the hydrodynamic value even for a magnetized plasma was found in MHD simulations by Cho et al. (2003). Nonetheless, equation (11) should be understood to be only approximate, since it assumes that fluid elements are transported adiabatically, which need not be the case when anisotropic conduction is operating. In reality, both the thermal conduction diffusion coefficient and the turbulent heat diffusion coefficient can be comparable, and either could dominate in a specific situation.

Thermal conduction may indirectly assist with turbulent heat diffusion, as it reduces the impact of buoyancy forces (and thus reduces ). Indeed, simulations by Sharma et al. (2009) show that metal mixing in a stratified plasma is much more efficient once conduction is at play, allowing much broader metallicity profiles, for this very reason. Naively, if we think of gas entropy as a scalar to be advected by turbulent motions, similar conclusions should hold, although of course the interaction between heat transport and dynamics requires detailed simulations. We shall investigate the relative role of all these heating processes in our simulations.

3 Methods

3.1 Initial conditions for the gas

The details of the numerical setup are described in Ruszkowski & Oh (2010; hereafter RO10). Here we summarize key differences. The cluster parameters used here are similar to those corresponding to cool-core cluster A2199. In addition to the NFW potential of the cluster halo, we also include the contribution from the central brightest cluster galaxy (BCG), which was not included in RO10. The gravitational potential is described by the sum of the term due to an NFW profile with a softened core


where is the smoothing core radius ( kpc), kpc is the usual NFW scale radius, and the BCG contribution which has a King profile:


where , kpc is the core radius for the BCG and km s is its line-of-sight velocity dispersion. The parameter in equation (13) determines the cluster mass and is of the order of the total cluster mass, . We then solve the equation of hydrostatic equilibrium assuming the entropy distribution as parametrized by Cavagnolo et al. (2009); see equations (16) & (17) of RO10. Note that we do not include the gravitational contribution from other galaxies (§3.2) in our initial conditions, so the system is not initially in full hydrostatic equilibrium. However, after an initial transient, it rapidly relaxes to a new equilibrium configuration.

The addition of the BCG has two effects. Firstly, due to the increased gravitational acceleration, it results in higher gas densities compared to the models we considered in Ruszkowski & Oh (2010). This allows for a more conservative analysis of the effect of cooling. In fact, the central density here is a factor of times higher, which, combined with a slightly lower assumed central temperature, results in a central cooling time which is nearly 5 times shorter. The higher adopted central density in this paper is in line with that observed in A2199 (Johnstone et al., 2002). Given this more stringent setup, some of the stable models in Ruszkowski & Oh (2010) would actually undergo a cooling catastrophe. Secondly, the change in the gravitational potential has consequences for the Brunt-Väisälä  frequency and the trapping of g-modes, as we discuss below.

The initial distribution of density and temperature is shown in Figure 1. The frequency of circular orbits and the Brunt-Väisälä  frequencies for a hydrodynamic and conducting fluid with this initial density and temperature profile are shown in Figure 2. For a mode with a given value of , -modes can be resonantly excited if . This therefore defines an outer trapping radius for such a mode. Note that both and are both strong functions of the gravitational potential. We have directly verified the resonance condition by running simulations both with and without the central cD galaxy; in the latter case, orbiting galaxies fail to excite volume-filling turbulence, which is to be expected since falls inward in this case, and the resonance condition is never satisfied (see also Lufkin et al. (1995); Kim (2007)). Note that fine-turning of the resonance condition is not necessary: the resonance is not very sharp (Balbus & Soker, 1990), and in practice galaxies with non-circular orbits excite modes with a variety of harmonics, some of which can potentially fall below .

The magnetic field setup was identical to that in Ruszkowski & Oh (2010): we generate statistically isotropic random-phase complex fields in Fourier space, with 3D Fourier amplitudes given by:


as appropriate for Kolmogorov turbulence, where and kpc is a smoothing wavelength. We then apply a divergence cleaning operator in space, and then inverse Fourier transform the field back to real space.

Figure 1: Initial electron number density (solid line) and temperature (dashed line) in the ICM of our simulated cluster.
Figure 2: The frequency of circular orbits (solid line), the Brunt-Väisälä  frequency (dashed line) for a hydrodynamic fluid and (dotted line) for a magnetized conducting fluid (all in [Hz]). The frequencies correspond to the initial density and temperature profile shown in Figure 1.

3.2 Initial conditions for the galaxies

The simulations must be initialized with a galaxy population, which has the appropriate spatial distribution, masses and velocities. Rather than relying upon cosmological simulations, we use an empirically grounded approach, which also has the advantage of speed and flexibility. How are the galaxies spatially distributed? From a sample of K-band selected galaxies within 93 clusters and groups, Lin et al. (2004) find that the galaxy number density profile in clusters is well described by the NFW profile (Navarro et al., 1997) with a concentration parameter , with no evidence for cluster mass dependence of the concentration. The theoretical justification for galaxies tracing the NFW profile is somewhat equivocal. If one attempts to use cosmological simulations to set up initial conditions, the radial distribution of subhalos in simulations is well-known to be less concentrated than the dark-matter, or ’anti-biased’ (Nagai & Kravtsov (2005), and references therein). This is due to tidal stripping of sub-halos in the central regions. On the other hand, simulations that include galaxy formation allow subhalos to be selected by stellar mass. This generally shows closer agreement with observed profiles, as the stellar mass (which is tightly bound) remains conserved while the dark matter is stripped from outer regions (Nagai & Kravtsov, 2005). Other investigators find that the fraction of such stellar-dominated halos is small, but caution that numerical resolution effects may preclude robust conclusions at this point (Dolag et al., 2009). Overall, we therefore simply employ the observational result that galaxies trace the NFW profile.

As for the galaxy masses, instead of using a Schechter function, we simply assume (as did, for instance, Subramanian et al. (2006); Kim (2007)) that all galaxies have the same mass. This is for two reasons. Firstly, this allows us to rapidly explore the effect of varying galaxy masses (due, for instance, to different efficiencies of tidal stripping). The assumption of a characteristic mass is reasonable: since dynamical friction scales as , most turbulent motions are induced by galaxies of mass , where most of the mass resides, rather than the more abundant lower mass galaxies. Indeed, we shall find that the induced gas motions are mostly sensitive to the mass of galaxies, and less sensitive to their number (§4.1). Previous hydrodynamic simulations found unchanged results with galaxies drawn from a Schechter distribution, if the characteristic break mass (Kim, 2007). Secondly, it allows us to directly calibrate against lensing estimates for subhalo mass fraction. Unlike K-band surveys, lensing is directly sensitive to total mass, but is generally only sensitive to subhalos with . Natarajan et al. (2009) find from the massive lensing cluster Cl 0024+16 that of the cluster mass can be attributed to substructure with , with typical masses (with a weak radial trend such that galaxies in the outer regions are more massive; see their Fig 6). Their results, including the mass function as a function of radius, is broadly consistent with the results of the Millenium simulation run, except that the typical masses of galaxies is lower in simulations by a factor of . This is subject to the uncertainties of extra binding due to a compact stellar halo mentioned above; note that masses of is also consistent with other observations from lensing (Shin et al., 2008) and galaxy wakes (Sakelliou et al., 2005). Below we explore a grid of models with varying galaxy mass, but never allowing the total substructure mass fraction to rise above . For simplicity in the code, the galaxies are modeled as point masses. Since we are primarily concerned with the excitation of g-modes on scales much larger than kpc galactic scales, we do not expect this simplification to significantly impact our results.

Given these assumptions, the most rigorous way to initialize galaxy velocities is to directly construct the distribution function from the density profile, using Eddington’s formula (Binney & Tremaine, 2008; Kazantzidis et al., 2004). However, velocity anisotropy is only easily incorporated in such models if it has certain parametric forms, as for instance in Osipkov-Merritt models. Instead, we construct a self-consistent velocity model via the local Maxwellian approximation: approximating the velocity dispersion tensor by a multi-variate Gaussian at each point, with dispersions given by the solution of the Jeans equation (Hernquist, 1993). This has the virtue of simplicity and flexibility. Note that such models may not be in strict equilibrium, and can demonstrate evolution (Kazantzidis et al., 2004). However, Springel et al. (2005) find the actual amount of relaxation to be small; furthermore, Faltenbacher et al. (2005), who directly simulate the motion of galaxies in clusters, find their velocity distribution is indeed closely Maxwellian, with good agreement between simulation results and equilibrium Jeans equation solutions. We therefore solve the Jeans equation assuming no rotational support or bulk streaming ():


where is the velocity anisotropy parameter4:


is the galaxy number density and is the combined cluster + cD galaxy gravitational potential. Note that we have not built self-consistent models and ignore the contribution of galaxies to the gravitational potential; for a large sub-halo mass fraction, the system is not in full equilibrium. In practice, this is a small effect, and the galaxy distribution does not evolve significantly over the course of our simulation.
What are appropriate assumptions for ? It may be estimated from observations via Jeans equation modelling, given knowledge of galaxies positions, the cluster potential, and line of sight velocities. A detailed study of 10 clusters using a spectroscopic sample of galaxies from SDSS and 2dF found galaxy orbits to be isotropic within the errors for most clusters (Hwang & Lee, 2008). An earlier paper, using ENACS data, found that the brightest ellipticals do not yield an equilibrium solution, while other ellipticals, SOs and early spirals have isotropic orbits, and late spirals prefer radial to isotropic orbits (Biviano & Katgert, 2004). Overall, we assume isotropy , and regard this as our default model. In passing, we note that one could easily incorporate the effects of the velocity anisotropy by, for example, considering fits to to measurements in simulations (Hoeft et al. (2004)). Given that the evidence for orbital anisotropy in observations is marginal to date, we defer the study of the effect of such orbital distributions to future work. We also note that preferentially radial orbits would enhance the restoration of conduction even further and strengthen our conclusions.

We solve equation (3.2) as an initial value problem, where . Having solved for and , we randomly sample from the multi-variate Gaussian at each position to create a realization of the velocity field. The above procedure allows us to initialize the simulations with a realization of galaxies with both masses and six-dimensional phase space coordinates (position and velocity).

Figure 3: The evolution of the velocity dispersion (left panel) and median velocity within the central 100 kpc (all in [km s]). Blue (red) lines are for 100 (200) galaxies respectively, for equally spaced masses ranging from M to M. The galaxy mass increases gradually from the lightest to the darkest color. The black dashed line is for the pure HBI case, where there is no stirring. The run is halted at early times if a cooling catastrophe occurs.

3.3 Simulation

The simulations were performed using the FLASH code (version 3.2). FLASH is a modular, parallel adaptive mesh refinement magnetohydrodynamic code. Magnetic field evolution was solved by means of a directionally unsplit staggered mesh algorithm (USM; Lee et al. (2009)). The USM module is based on a finite-volume, high-order Godunov scheme combined with constrained transport method (CT). This approach guarantees divergence-free magnetic field distribution. We implemented the anisotropic conduction unit following the approach of Sharma & Hammett (2007). More specifically, we applied monotonized central (MC) limiter to the conductive fluxes. This method ensures that anisotropic conduction does not lead to negative temperatures in the presence of steep temperature gradients. The three-dimensional computational domain was approximately 1 Mpc on each side, enclosing a large fraction of the cluster. The central regions of the cluster had an enhanced refinement level. The maximum spatial resolution for 6 levels of refinement was kpc. The simulations were performed on a 384-processor cluster located at the Michigan Academic Computing Center at the University of Michigan in Ann Arbor and on the Columbia supercomputer at NASA Ames.

4 Results

We performed a total of 16 runs including radiative cooling, anisotropic thermal conduction and self-gravitating particles to emulate the gas “stirring” by galaxies. We considered a uniform grid of parameters: 50, 100, 150, and 200 galaxies characterized by masses of (0.3, 0.6, 0.9, 1.2) M. With our cluster mass of , these parameters corresponds to a mass fraction in galaxies ranging from to a maximum of . For instance, for galaxies with mass , our grid corresponds to . We also performed two control runs: one without the galaxies (and hence without stirring) to isolate the effect of heat buoyancy instability, and one without conduction to isolate the effect of dynamical friction heating by galaxies.

4.1 Gas velocities and volume filling of turbulence

Figure 3 (left panel) shows the evolution of the velocity dispersion measured within 100 kpc from the cluster center. Thin blue (red) lines are for 100 (200) galaxies respectively and for equally spaced masses ranging from to . The mass increases gradually from the lightest to the darkest color. The black dashed line is for the pure HBI case. The HBI case and lighter-colored curves are evolved for shorter times. These runs suffer from overcooling and the central temperatures reaches the low temperature threshold at which point the simulation is stopped. The right panel in Figure 3 shows the median velocity within 100 kpc; the color coding corresponds to that in the left panel. It is clear from these figures that there is a clear trend for the velocity dispersion or the median velocity to increase with the typical galaxy mass. A similar, albeit weaker, trend is seen for the galaxy number. This is consistent with the findings of Kim (2007) who found in pure hydrodynamic simulations, that the gas velocity dispersion scales as , where and are the number and mass of galaxies, respectively. Note that a scaling is consistent with dynamical friction in the linear regime, since for dynamical friction.

As and increase, the cooling catastrophe is delayed, and is completely staved off at the upper envelope of these parameters. In this respect our MHD simulations differ markedly from those of Kim (2007), who found that a cooling catastrophe was inevitable in purely hydrodynamic simulations, for all portions of parameter space. We explore these differences further in §4.5. Note that in our case the velocity dispersion seems to increase more slowly than . Besides the inclusion of MHD in our simulations, differing results could be due to a variety of factors, including the different assumed distribution of galaxies. Note that the stated number of galaxies are distributed over the entire cluster; the number of galaxies in the inner regions which actually result in trapped g-modes is actually considerably smaller, and subject to Poisson fluctuations. Also, the introduction of more galaxies and/or increasing their mass does not cause the velocity dispersion to increase without limit; instead, the growth in velocity dispersion appears to saturate. Kim (2007) also observed this in his hydrodynamic simulations, and attributed it to loss of resonant excitation once density fluctuations become large and the background is nonlinear. We see the same saturation on the same yr timescale, in simulations with driven volume-filling turbulence, where resonant excitation of modes is not an issue (see §3.4 of RO10). The asymptotic velocities of , while generally insufficient for turbulent heating to be important, is enough to restore thermal conduction and enable turbulent heat diffusion.
The comparison between the velocity dispersion and median velocity reveals that that both of these quantities are comparable, as might be expected for volume filling turbulence.

4.2 Bias in magnetic field orientation

The evolution of the anisotropy in the orientation of magnetic fields is shown in Fig. 4. The definition of this parameter is similar to that for galaxy velocity defined in Eq. (16). The only difference is that the velocity dispersions are replaced by magnetic field dispersions. Thus, corresponds to isotropic magnetic fields, whilst corresponds to progressively more tangential (radial) fields respectively. Thin blue lines are for 100 galaxies and for equally spaced masses ranging from M to M. The galaxy mass increases gradually from the lightest to the darkest color. Thick red curves are the corresponding lines for 200 galaxies. The black dashed line is for the pure HBI case. As expected, when stirring is weak, the HBI prevails, leading to a systematic tangential bias in the orientation of magnetic fields. This insulates the core against thermal conduction, leading to a cooling catastrophe. The fields become more tangential with time and the cluster eventually suffers from overcooling. On the other hand, for increasingly vigorous stirring (i.e., increasing the individual masses or number of galaxies), the field becomes increasingly isotropic, and a cooling catastrophe is averted.

These results are consistent with the driven turbulence simulations in RO105, and can be broadly understood in terms of the simple Froude/Richardson number criterion outlined in §2. The main difference in the more realistic scenario we present here is that the discrete nature of the stirrers and resonant excitation process introduces greater stochasticity and time-dependent fluctuations in the velocity field and magnetic anisotropy (e.g., compare the smooth curves in Fig 3 & 4 of RO10 with their noisier equivalents Fig 3 & 4 of this paper). But the main physical conclusions are unchanged. It is also interesting to note that while the velocity dispersion is only weakly dependent on the number of galaxies (depending more sensitively on galaxy masses), the magnetic anisotropy shows somewhat greater sensitivity. In particular, the magnetic field anisotropy cannot simply be predicted from the instantaneous velocity dispersion, as in a naive application of a Froude/Richardson criterion. We saw similar behavior in RO10, where runs with similar asymptotic velocity dispersions had similar velocity anisotropies, but markedly different magnetic anisotropies. The advected magnetic field is sensitive to the integrated past displacement history of a fluid element, and not merely the instantaneous velocity field.

Figure 4: The evolution of the anisotropy in the orientation of magnetic fields. Vanishing corresponds to isotropic fields. The more negative becomes, the more tangential the fields are. Thin blue lines are for 100 galaxies and for equally spaced masses ranging from M to M. The galaxy mass increases gradually from the lightest to the darkest color. Thick red curves are the corresponding lines for 200 galaxies. The black dashed line is for the pure HBI case. Runs are terminated when a cooling catastrophe sets in.

4.3 Evolution of Gas Temperature and Entropy

In Figure 5 we show the evolution of temperature profiles. This figure is for the models where heating is more efficient. Specifically, it corresponds to the following pairs of parameters (150,1.2), (200,0.9), (200,1.2), where the first number in the parenthesis is the number of galaxies and the second is the galaxy mass in M. Progressively older profiles correspond to systematically brighter colors. The final time corresponds to 5 Gyr and the curves are plotted every 0.1 Gyr. As is clearly seen in this figure, these models do not lead to the cooling catastrophe. Several features are of interest. The temperature profile does not asymptote toward an isothermal profile, as is generically the case when thermal conduction alone offsets cooling (Bregman & David, 1988; Guo & Oh, 2008; Conroy & Ostriker, 2008). Despite the fact that we have not introduced an additional source of central heating such as an AGN, the cluster is able to remain in a thermally stable CC (i.e., with a central temperature which is lower than at the cooling radius) state via heat transport from the outer heater reservoir alone. Without fine-tuning, this is impossible to achieve with thermal conduction alone (when the cluster either becomes isothermal or undergoes a cooling catastrophe). Finally, the temperature profile is not always monotonic, but occasionally increases inward—a situation which is thermodynamically impossible if thermal conduction alone is operating. Note that these fluctuations are transient; such reversals are not present in the later stages of the evolution (progressively lighter blue curves correspond to later times). As we shall see in §4.5, all of these features hint that an additional heat transport process is at play: turbulent heat diffusion.

Figure 5: The evolution of temperature profiles for the strong heating models. The panels correspond to the the following pairs of parameters: (150,1.2), (200,0.9), (200,1.2), where the first number in the parenthesis is the number of galaxies and the second is the galaxy mass in M, from left to right respectively. The final time corresponds to 5 Gyr and the curves are plotted every 0.3 Gyr.
Figure 6: The evolution of temperature profiles for the weak heating models. From left to right are the results for the following sets of parameters: (100, 0.3), (50, 0.6), (50, 0.3), where the first number in the parenthesis is the number of galaxies and the second is the galaxy mass in M. The curves are shown every 0.1 Gyr.

Figure 6 shows the temperature evolution for the parameters where the heating is least efficient. From left to right shown are: (100, 0.3), (50, 0.6), (50, 0.3). Here, the profiles are shown more frequently then in Figure 3 to better capture the evolution of the system just before the imminent cooling catastrophe. For the same reason, we also extend the radial scale to smaller distances from the cluster center to show how the system becomes thermally unstable. It is evident that in all three cases, the cluster quickly evolves toward a cooling catastrophe. In the final stages of the process, the cooling is so fast at the very center that the gas accretion accelerates so much that adiabatic compression in the shells surrounding the center can heat the gas up (e.g., see last profile in the right panel).

Finally, in Fig 7 we show the entropy profiles (where entropy is defined as ) for the strong heating models. The central entropy grows somewhat, consistent with the rise in temperature, and as might be expected if heating by conduction and/or turbulent heat diffusion were taking place. However, these profiles show that turbulent convection/stirring is still a relatively gentle process; we do not see the flat isentropic central profile which might be expected if turbulent convection were extremely efficient. Instead, the fluid always remains stably stratified by entropy, which steadily increases outward at all times.

Figure 7: The evolution of entropy profiles for the strong heating models. From left to right are the results for the following sets of parameters: (150, 1.2), (200, 0.9), (200, 1.2), where the first number in the parenthesis is the number of galaxies and the second is the galaxy mass in M. The curves are shown every 0.1 Gyr.
Figure 8: Time until cooling catastrophe as a function of the number of galaxies and galaxy mass (in units of ). Contours are in Gyr. Models that correspond to 6 Gyr (the maximum duration of the simulations) are thermally stable. See text for details.

As discussed above the cluster will develop a cooling catastrophe if the number of galaxies and/or their masses are too small. The time it takes for the cluster to reach this point (essentially the effective cooling time) is plotted in Figure 8 as a function of galaxy number and mass. The contours are plotted every Gyr. The models that exhibit the effective cooling time of 6 Gyr (the maximum simulation run time) are thermally stable. We point out that in practice the models that possess cooling times to 4 Gyr could be considered stable as they are likely to experience cluster mergers that may reset the conditions in the ICM and further slow down or essentially delay the cooling process. In any case, as can be seen in Figure 5, a substantial fraction of the models shows appreciably long effective cooling times. As a technical note, we add that the reason for the lack of monotonicity in some of the contour lines as a function of galaxy number is that a single random seed was used to generate the conditions for a given number of galaxies and varying galaxy masses.

Figure 9: The evolution of normalized vorticity (left panel) and normalized magnetic pressure. See text for definition of normalization. The curves correspond to the same dataset as that shown in Figure 3 and the meaning of lines is the same as in that figure.

4.4 Generation of vorticity and magnetic fields

As we have previously seen, -modes must be excited for stirring by galaxies to excite volume-filling turbulence. These -modes also induce vorticity (equation 2). Vorticity is therefore an excellent tracer of the growth of -modes. We compute the evolution of vorticity in the central 100 kpc to assess if modes are indeed generated and trapped. Figure 9 (left panel) shows the evolution of the square of the scaled vorticity for the same set of parameters as in Figure 3 that shows velocity dispersion and median velocity. The scaled vorticity is defined as , where kpc and km s are the reference lengthscale and velocity, respectively. Thin blue lines are for 100 galaxies and red ones for 200 galaxies. Galaxy gasses range from M to M and are uniformly sampled (lighter colors are for lighter galaxies). The black dashed line corresponds to the pure HBI case. A clear trend for the vorticity to increase with time is seen in this figure, suggesting that modes are present and at least partially trapped, leading to the volume-filling turbulence seen.

As discussed in §2, a growth in vorticity might also lead to growth in the magnetic field; the possibility that magnetic fields could be turbulently amplified in clusters has been repeatedly raised (e.g., Ruzmaikin et al. (1989); Subramanian et al. (2006); Ryu et al. (2008); Cho et al. (2009)). Whilst a detailed study is beyond the scope of this paper, we check whether these theoretical expectations are satisfied in our simulations. Fig. 9 shows that the magnetic energy density indeed grows in tandem with vorticity, with more vigorous stirring corresponding to greater field amplification. However, the characteristic growth time appears to be somewhat longer. Note that the simulations were initialized with extremely small magnetic fields: the initial plasma beta , where is typically measured in the ICM. These small initial fields were for computational convenience (since the MHD approximation is satisfied with a trivially small magnetic field), and to ensure that magnetic fields never become dynamically important6. Hence, despite growing by a factor of  50, the magnetic energy density has not yet reached its saturated state, and is not yet in equipartition with turbulence. Nonetheless, the turbulent amplification of the B-fields, which mirrors the growth of vorticity, is a robust result.

4.5 Relative contribution to gas heating

Figure 10: Heating to cooling ratios for the case of steady volume-filling turbulence driven by a source function (see text and RO10 for more details). Top row: Conductive (solid line) and convective heating to cooling ratios as a function of time for the ICM within 100 kpc (left panel) and 50 kpc (right panel) from the cluster center, respectively. Bottom row: Conductive (left panel) and convective (right panel) heating to cooling ratios as a function of radius. Progressively lighter blue color denotes later times. The curves are plotted every Myr.
Figure 11: Same as for Fig 10 for the model of turbulence stirred by galaxy motions. Top row: Conductive (solid line) and convective heating (dashed line) to cooling ratios as a function of time for the ICM within 100 kpc (left panel) and 50 kpc (right panel) from the cluster center, respectively. Bottom row: Conductive (left panel) and convective (right panel) heating to cooling ratios as a function of radius. The curves are plotted every Myr; progressively lighter blue color denotes later times.

In §4.3 we noted a number of interesting features in the temperature profiles of our stable clusters. They remained stable CC clusters, neither becoming isothermal nor developing cooling catastrophes, as clusters stabilized solely by thermal conduction generally do. Furthermore, the central temperature showed time-dependent oscillations, sometimes becoming hotter than gas further out. A temperature inversion would not happen if only thermal conduction was at play. This behooves us to take a closer look at what actually stabilizes the customary thermal runaway. We have already discussed the effect that turbulence can have on thermal conduction, by tangling field lines and countering the HBI. However, turbulence itself can be a source of heating, either via viscous dissipation of turbulence, or turbulence diffusion of high entropy gas into low entropy regions (e.g., Dennis & Chandran (2005), and references therein). Let us examine these in turn.

As long as there is sufficient separation of scales that an inertial range can develop (such that the energy per unit mass per unit time is independent of scale), the heating rate from dissipation of turbulent motions is independent of the nature of viscosity. In particular, it is unimportant if our numerical viscosity is different from the actual physical viscosity in the ICM. The heating rate per unit volume due to dissipation of such motions is (Dennis & Chandran, 2005):


where is the dominant lengthscale, is the energy density in turbulence, and is the eddy-turnover time on the dominant lengthscale. To estimate , we can note that vorticity has units of , and that our scaled vorticity in Fig 9 is . This implies


Consistently, note that the vorticity in Fig 9 indeed takes Gyr to rise to its asymptotic value. This implies that the heating time for turbulent dissipation of motions is:


where is the thermal energy density, and we have defined the turbulent Mach number (note that our quoted velocities are in 3D). While there are factors of order unity uncertainty, it is clear that the mild subsonic motions we explore are a negligible source of heating via viscous dissipation (and consistent with other estimates; Dennis & Chandran (2005)). This also implies that dynamical friction heating due to galaxy motions (El-Zant et al., 2004; Kim et al., 2005; Kim, 2007; Conroy & Ostriker, 2008; Birnboim & Dekel, 2010) is not the source of heating which averts the cooling catastrophe in these simulations.

On the other hand, turbulent heat diffusion is not negligible. One can estimate its contribution from a simple mixing length prescription as in equation 11; this shows that it can be at least comparable to and may exceed the thermal conduction contribution. However, the coefficient of turbulent diffusivity, , is only approximate and subject to order unity corrections. Since we have full knowledge of the density, velocity, temperature and magnetic fields in our simulations, we can attempt to directly compute the heating contributions from thermal conduction and turbulent heat convection. In particular, at a radius we can calculate the inward heat flux due to conduction:


where is a unit vector pointing in the direction of the magnetic field and is the Spitzer-Braginskii conduction coefficient given by erg scmK, as well as the convective heat flux (Parrish et al., 2008):


where is the spatial average of quantity in the shell and is the local deviation of that quantity from its average; generally the second term is dominant. We can then compare these to the total rate of energy loss within radius due to radiative cooling. We can also compute the volumetric heating rate due to these two processes, via , , although these are of course much noisier quantities.

It is useful perhaps to begin by considering a case where the properties of turbulence are well known: the ’strong’ driven turbulence case of RO10, which has volume filling turbulence by construction, and rms velocities of . Conductive (solid line) and convective (dashed line) heating to cooling ratios as a function of time for the ICM within 100 (50) kpc from the cluster center are shown in the upper left (right) panel of Fig 10. It is clear that conduction only contributes of the heat necessary to overcome cooling and convective heat flow is an important part of the energy budget; indeed, in the central regions turbulent advection of heat is the dominant heating process (note that the cluster is not in complete equilibrium, so the sum of the two ratios is not necessarily unity). The convective heat flow shows time-dependent fluctuations, as might be expected. In the bottom panels, we show the volumetric heating to cooling ratios as a function of radius, for conduction (left panel), and turbulent convection (right panel). The curves are plotted every  100 Myr; progressively lighter colors denote later times. The divergence of heat fluxes is a much noisier quantity, as reflected in the plots. Nonetheless, it is clear from the plots that convective heating dominates near the center, whilst conductive heating dominates further out. This reminiscent of stable hybrid AGN+conduction heating models (Ruszkowski & Begelman, 2002) where the AGN heats the cluster center and conduction is important further out.

In Fig. 11, we show the same plots, but for the case where turbulence is due to stirring by galaxies. All results presented in this figure are for the case of 200 galaxies, each with M; this is stable against a cooling catastrophe. As before, conduction is only a fraction of the energy budget. However, in this case convective heating shows dramatic oscillations as a function of time; the amplitude of the oscillations near the center is much larger than in the driven turbulence case . The reason for this is that the dominant lengthscales of motion are comparable or larger than the depicted radii, as might be expected if -modes are excited (since most of the energy in -modes are in the largest lengthscales, comparable to the trapping radius). For instance, from §4.4, a typical lengthscale on which vorticity is excited is . Since fluctuations in the velocity field span larger scales than the ones under interest, our calculation of will show strong time dependence (however, the calculations of the conductive heat flux are of course still valid. Note that has to be unity on average, since the cluster is stabilized against a cooling catastrophe). As noted earlier, Poisson fluctuations in the number of galaxies in the core will also drive time-dependent fluctuations in the velocity field. The gas is sloshing in the potential well; we observe this directly too in the simulations, as the gas pressure maximum wanders in time from the center of the potential well. Nevertheless, despite the breakdown of equation 22 in a rigorous sense, it is clear from the amplitude of fluctuations in the bottom panels of Fig 11 that (as in the driven turbulence case) conductive heating increases outwards in radius, while convective heating is more important near the center. In particular, the dominance of convective heating near the center, and its positive and negative fluctuations, allow us to understand the fluctuating temperature profiles seen in Fig 5. Since conduction is only a part of the energy budget, there is no reason for the stabilized temperature profile to approach isothermality. Furthermore, the reason why the central temperature gradient can occasionally become inverted (with the center hotter than its surroundings) is clear: if a high-entropy fluid element is compressed at the center, this will result in higher central temperatures. While thermal conduction seeks to make the fluid isothermal (since heat flows down the temperature gradient), turbulent heat diffusion seeks to make the fluid isentropic (since heat flows down the entropy gradient). In this sense, the subsonic turbulence induced by galaxies results in only mild convection, since as seen in Fig 7, the gas remains convectively stable with entropy increasing monotonically outward. Whilst we have not directly calculated the diffusion of metals directly, this also suggests that metal mixing to larger radii will be somewhat enhanced (so that metals will have a broader distribution than the galaxies), but not greatly so. Indeed, a mixing-length theory calculation of metal dispersal via turbulent diffusion by Rebusco et al. (2005), who assumes levels of turbulence very similar to those we have simulated, shows excellent agreement with observations.

Could turbulent heat diffusion alone stabilize a thermal runaway? We tested this hypothesis by running purely hydrodynamic simulations of our most vigorous stirring case (200 galaxies of mass ), without thermal conduction. The results are shown in Fig. 12. The cluster rapidly undergoes a cooling catastrophe (consistent with the results of Kim (2007)), demonstrating that thermal conduction is an essential element in stabilizing the cluster. Note that this model already represents the most extreme choice in parameter space of the amount of stirring possible by galaxies. Also in line with these results, purely hydrodynamic simulations of ’sloshing’ZuHone et al. (2010) show that the cooling catastrophe can be delayed by not disrupted. Besides providing the dominant source of heating in the outer regions, thermal conduction also reduces stabilizing buoyancy forces (as discussed in §2) and thus enables more rapid, efficient mixing and turbulent heat diffusion. Passive scalars such as metals are more efficiently advected in the presence of conduction (Sharma et al., 2009a), and the same is likely true of the advection of entropy. Thus, intriguingly, whilst neither process alone can stabilize the cluster, the interplay between turbulence and conduction does permit stability: turbulence enables conduction, and conduction enables turbulence.

Figure 12: Evolution of temperature profiles for a purely hydrodynamic run (no thermal conduction) with 200 galaxies of mass ; our strongest heating model. Curves are shown every 0.1 Gyr. The cluster rapidly undergoes a cooling catastrophe, demonstrating that thermal conduction is an essential element in stabilizing the cluster; purely hydrodynamic turbulent heat diffusion alone is insufficient.

5 Conclusions

Using three-dimensional MHD simulations, we have studied the effect of anisotropic thermal conduction and stirring motions due to galaxies orbiting in the cluster potential on the effective cooling rate in cluster cool cores. Such galaxies excite mild subsonic turbulence with . We find that a combination of thermal conduction and turbulent heat transport can stabilize the cluster, for realistic parameter choices consistent with gravitational lensing observations of substructure in clusters. Unlike much previous work, there is no subgrid physics in our simulations: we do not invoke sub-grid prescriptions for the topology of the magnetic field (which affects the effective thermal conductivity), the magnitude and volume-filling factor of turbulence, which is calculated directly from the gravity/hydro solver (unlike previous work (Ruszkowski & Oh, 2010; Parrish et al., 2010) in which volume-filling turbulence is inserted by hand), or turbulent heat diffusion (which is directly simulated). We have also simulated a cluster with significantly higher central density than in RO10, and still found it to be thermally stable. Other salient points include:

  • In order for galaxies to excite volume-filling turbulence, rather than have turbulence confined to galactic wakes, they must excite -modes, which requires that , where is the Brunt-Väisälä  frequency appropriate when thermal conduction timescales are rapid (equation 2). On the other hand, overwhelming the stabilizing buoyancy forces to randomize the magnetic field requires that . These two requirements can be simultaneously satisfied since for Kolomogorov turbulence; hence, the low-frequency, large scale modes can be trapped, while the high-frequency, small scale modes overcome the HBI.

  • We observed strong growth in vorticity, which is a good tracer of the growth of -modes. We also observed turbulent amplification of B-fields in tandem with vorticity.

  • Thermal conduction provided about of the heating budget, with the rest coming from turbulent heat diffusion. Viscous dissipation of turbulent motions (and hence dynamical friction heating) is negligible. Turbulent heat diffusion tends to be more important in the center of the cluster, while conduction plays a greater role further out. The predominance of turbulent heat diffusion in the center—which is powered by motions on large scales—implies that it exhibits oscillations about the equilibrium temperature profile, and can occasionally exhibit small temperature inversions as high entropy fluid elements are compressed near the center. However, conduction plays a crucial part of the story; our most extreme stirring case still suffered a cooling catastrophe if thermal conduction was omitted. Besides supplying heat further out in the cluster, conduction also reduces stabilizing buoyancy forces and enables more efficient turbulent heat diffusion. It appears that turbulence enables conduction to operate, as well as vice-versa. The details of the interplay between turbulence and conduction, as well as the diffusion of metals in our stirring simulations, are interesting topics for future work.

In this paper, we have focussed on a time-steady source of turbulence—stirring by galaxy motions—but we stress that other intermittent sources of turbulence such as mergers or AGN outbursts, can also contribute. Indeed, a sudden rise in heat transport processes such as conduction and turbulent heat diffusion due to an increase in turbulence could effect a CC to NCC transition (Guo & Oh, 2009; Ruszkowski & Oh, 2010; Parrish et al., 2010). Other processes which could reorient field lines in galaxy cluster include rising bubbles, which could amplify and straighten magnetic fields in their wake (Ruszkowski et al., 2007; Guo et al., 2008; Bogdanović et al., 2009). In the future, observations of Faraday rotation by SKA (Bogdanovic et al., 2010) or magnetic draping around galaxies orbiting the cluster center (Pfrommer & Dursi, 2010), could probe the topology of magnetic field lines and test these ideas. Finally, these ideas about the interplay between between the thermal conduction, the HBI and turbulence in the inner regions of the cluster also apply with equal force to the interplay between conduction, the MTI and turbulence in the outer regions of the cluster, which we present elsewhere (Ruszkowski et al 2010).


The software used in this work was in part developed by the DOE-supported ASC/Alliance Center for Astrophysical Thermonuclear Flashes at the University of Chicago. MR thanks Jeremy Hallum for his invaluable help with maintaining the computing cluster at the Michigan Academic Computing Center where most of the computations were performed. We are indebted to Dongwook Lee, the author of the MHD module in FLASH for letting us use a proprietary three-dimensional set of MHD modules. MR thanks Justin Nieusma for technical assistance with performing the simulations. We thank Ian Parrish, Eliot Quataert, Elena Rasia, Prateek Sharma, Min-Su Shin, Maxim Markevitch, John ZuHone, Paul Nulsen, Aneta Siemiginowska, Christine Jones, Larry David, Bill Forman, Renato Dupke, Milos Milosavljević, John ZuHone, Matt Kunz, and Alex Schekochihin for discussions. SPO acknowledges support by NASA grant NNG06GH95G, and NSF grant 0908480. MR and SPO thank Institute of Astronomy, Cambridge, UK and Max Planck Institute for Astrophysics, Garching, Germany for their hospitality.


  1. Although a WKBJ analysis formally breaks down in this regime, a subsequent numerical study (Lufkin et al., 1995) showed that many of the linear theory results are still valid.
  2. Note, however, that this analogy is imperfect, since , which leads to a nonlinear coupling in the equations, whereas no such relation exists between and .
  3. In general, estimates based on rotation measure (RM) lead to stronger magnetic fields, while those based on synchrotron and inverse Compton (IC) analysis give weaker fields. However, RM methods may overestimate fields if single-scale magnetic field correlation length is used (Newman et al., 2002) or when the small-scale fluctuations in density and magnetic field are correlated in a turbulent medium (Beck et al., 2003). Moreover, these estimates depend on whether radio sources used to probe the field strength are embeded in the ICM, with smaller values inferred when background sources rather than embedded ones are used (Carilli & Taylor, 2002).
  4. From cosmological simulations, Benson (2005) finds that radial and tangential velocities can be correlated (at least at the time of merger), a detail we ignore.
  5. Although note that all but one of the simulations in RO10 were adiabatic simulations; by contrast, all the simulations presented here simultaneously include radiative cooling.
  6. Thus, as least in these simulations, the HBI is stabilized by the stirring motions and not by magnetic tension.


  1. Ascasibar Y., Markevitch M., 2006, ApJ, 650, 102
  2. Balbus S. A., 2000, ApJ, 534, 420
  3. Balbus S. A., Soker N., 1990, ApJ, 357, 353
  4. Beck R., Shukurov A., Sokoloff D., Wielebinski R., 2003, A&A, 411, 99
  5. Benson A. J., 2005, MNRAS, 358, 551
  6. Binney J., Tabor G., 1995, MNRAS, 276, 663
  7. Binney J., Tremaine S., 2008, Galactic Dynamics: Second Edition. Princeton University Press
  8. Birnboim Y., Dekel A., 2010, ArXiv e-prints
  9. Biviano A., Katgert P., 2004, A&A, 424, 779
  10. Bogdanovic T., Reynolds C., Massey R., 2010, ArXiv e-prints
  11. Bogdanović T., Reynolds C. S., Balbus S. A., Parrish I. J., 2009, ApJ, 704, 211
  12. Bregman J. N., David L. P., 1988, ApJ, 326, 639
  13. Carilli C. L., Taylor G. B., 2002, ARA&A, 40, 319
  14. Cavagnolo K. W., Donahue M., Voit G. M., Sun M., 2009, ApJS, 182, 12
  15. Cho J., Lazarian A., Honein A., Knaepen B., Kassinos S., Moin P., 2003, ApJ, 589, L77
  16. Cho J., Vishniac E. T., Beresnyak A., Lazarian A., Ryu D., 2009, ApJ, 693, 1449
  17. Churazov E., Forman W., Jones C., Sunyaev R., Böhringer H., 2004, MNRAS, 347, 29
  18. Churazov E., Sunyaev R., Forman W., Böhringer H., 2002, MNRAS, 332, 729
  19. Conroy C., Ostriker J. P., 2008, ApJ, 681, 151
  20. Dennis T. J., Chandran B. D. G., 2005, ApJ, 622, 205
  21. Dolag K., Borgani S., Murante G., Springel V., 2009, MNRAS, 399, 497
  22. El-Zant A. A., Kim W.-T., Kamionkowski M., 2004, MNRAS, 354, 169
  23. Enßlin T. A., Vogt C., 2006, A&A, 453, 447
  24. Evrard A. E., 1990, ApJ, 363, 349
  25. Fabian A. C., Sanders J. S., Allen S. W., Crawford C. S., Iwasawa K., Johnstone R. M., Schmidt R. W., Taylor G. B., 2003, MNRAS, 344, L43
  26. Faltenbacher A., Kravtsov A. V., Nagai D., Gottlöber S., 2005, MNRAS, 358, 139
  27. Goldshmidt O., Rephaeli Y., 1993, ApJ, 411, 518
  28. Guo F., Oh S. P., 2008, MNRAS, 384, 251
  29. Guo F., Oh S. P., 2009, MNRAS, 400, 1992
  30. Guo F., Oh S. P., Ruszkowski M., 2008, ApJ, 688, 859
  31. Hernquist L., 1993, ApJS, 86, 389
  32. Hoeft M., Mücket J. P., Gottlöber S., 2004, ApJ, 602, 162
  33. Hwang H. S., Lee M. G., 2008, ApJ, 676, 218
  34. Iapichino L., Niemeyer J. C., 2008, MNRAS, 388, 1089
  35. Johnstone R. M., Allen S. W., Fabian A. C., Sanders J. S., 2002, MNRAS, 336, 299
  36. Kazantzidis S., Magorrian J., Moore B., 2004, ApJ, 601, 37
  37. Kim W., 2007, ApJ, 667, L5
  38. Kim W., El-Zant A. A., Kamionkowski M., 2005, ApJ, 632, 157
  39. Kim W., Narayan R., 2003a, ApJ, 596, L139
  40. Kim W.-T., Narayan R., 2003b, ApJ, 596, 889
  41. Kunz M. W., Schekochihin A. A., Cowley S. C., Binney J. J., Sanders J. S., 2010, ArXiv e-prints
  42. Lee D., Deane A. E., Federrath C., 2009, in N. V. Pogorelov, E. Audit, P. Colella, & G. P. Zank ed., Astronomical Society of the Pacific Conference Series Vol. 406 of Astronomical Society of the Pacific Conference Series, A New Multidimensional Unsplit MHD Solver in FLASH3. pp 243–+
  43. Lin Y.-T., Mohr J. J., Stanford S. A., 2004, ApJ, 610, 745
  44. Lufkin E. A., Balbus S. A., Hawley J. F., 1995, ApJ, 446, 529
  45. Markevitch M., Gonzalez A. H., David L., Vikhlinin A., Murray S., Forman W., Jones C., Tucker W., 2002, ApJ, 567, L27
  46. McCarthy I. G., Babul A., Bower R. G., Balogh M. L., 2008, MNRAS, 386, 1309
  47. McNamara B. R., Nulsen P. E. J., 2007, ARA&A, 45, 117
  48. Nagai D., Kravtsov A. V., 2005, ApJ, 618, 557
  49. Nagai D., Kravtsov A. V., Kosowsky A., 2003, ApJ, 587, 524
  50. Narayan R., Medvedev M. V., 2001, ApJ, 562, L129
  51. Natarajan P., Kneib J., Smail I., Treu T., Ellis R., Moran S., Limousin M., Czoske O., 2009, ApJ, 693, 970
  52. Navarro J. F., Frenk C. S., White S. D. M., 1997, ApJ, 490, 493 (NFW)
  53. Newman W. I., Newman A. L., Rephaeli Y., 2002, ApJ, 575, 755
  54. Norman M. L., Bryan G. L., 1999, in H.-J. Röser & K. Meisenheimer ed., The Radio Galaxy Messier 87 Vol. 530 of Lecture Notes in Physics, Berlin Springer Verlag, Cluster Turbulence. pp 106–+
  55. Parrish I. J., Quataert E., Sharma P., 2009, ApJ, 703, 96
  56. Parrish I. J., Quataert E., Sharma P., 2010, ApJ, 712, L194
  57. Parrish I. J., Stone J. M., 2005, ApJ, 633, 334
  58. Parrish I. J., Stone J. M., Lemaster N., 2008, ApJ, 688, 905
  59. Petrosian V., 2001, ApJ, 557, 560
  60. Pfrommer C., Dursi J. L., 2010, Nature Physics, 6, 520
  61. Quataert E., 2008, ApJ, 673, 758
  62. Rebusco P., Churazov E., Böhringer H., Forman W., 2005, MNRAS, 359, 1041
  63. Ruszkowski M., Begelman M. C., 2002, ApJ, 581, 223
  64. Ruszkowski M., Brüggen M., Begelman M. C., 2004a, ApJ, 611, 158
  65. Ruszkowski M., Brüggen M., Begelman M. C., 2004b, ApJ, 615, 675
  66. Ruszkowski M., Enßlin T. A., Brüggen M., Heinz S., Pfrommer C., 2007, MNRAS, 378, 662
  67. Ruszkowski M., Oh S. P., 2010, ApJ, 713, 1332
  68. Ruzmaikin A., Sokolov D., Shukurov A., 1989, MNRAS, 241, 1
  69. Ryu D., Kang H., Cho J., Das S., 2008, Science, 320, 909
  70. Sakelliou I., Acreman D. M., Hardcastle M. J., Merrifield M. R., Ponman T. J., Stevens I. R., 2005, MNRAS, 360, 1069
  71. Sanders J. S., Fabian A. C., Smith R. K., Peterson J. R., 2009, ArXiv e-prints, 0911.0763
  72. Scannapieco E., Brüggen M., 2008, ApJ, 686, 927
  73. Schekochihin A. A., Cowley S. C., 2007, Turbulence and Magnetic Fields in Astrophysical Plasmas. Springer
  74. Schuecker P., Finoguenov A., Miniati F., Böhringer H., Briel U. G., 2004, A&A, 426, 387
  75. Sharma P., Chandran B. D. G., Quataert E., Parrish I. J., 2009a, ApJ, 699, 348
  76. Sharma P., Chandran B. D. G., Quataert E., Parrish I. J., 2009b, ArXiv e-prints
  77. Sharma P., Colella P., Martin D. F., 2009, ArXiv e-prints
  78. Sharma P., Hammett G. W., 2007, Journal of Computational Physics, 227, 123
  79. Shin M., Strauss M. A., Oguri M., Inada N., Falco E. E., Broadhurst T., Gunn J. E., 2008, AJ, 136, 44
  80. Springel V., Di Matteo T., Hernquist L., 2005, MNRAS, 361, 776
  81. Subramanian K., Shukurov A., Haugen N. E. L., 2006, MNRAS, 366, 1437
  82. Vazza F., Brunetti G., Gheller C., Brunino R., 2010, New Astronomy, 15, 695
  83. Vazza F., Brunetti G., Kritsuk A., Wagner R., Gheller C., Norman M., 2009, A&A, 504, 33
  84. Vogt C., Enßlin T. A., 2005, A&A, 434, 67
  85. Voigt L. M., Fabian A. C., 2004, MNRAS, 347, 1130
  86. Voit G. M., Cavagnolo K. W., Donahue M., Rafferty D. A., McNamara B. R., Nulsen P. E. J., 2008, ApJ, 681, L5
  87. Zakamska N. L., Narayan R., 2003, ApJ, 582, 162
  88. ZuHone J. A., Markevitch M., Johnson R. E., 2010, ApJ, 717, 908
Comments 0
Request Comment
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 minumum 40 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description