Cold fronts: probes of plasma astrophysics in galaxy clusters

Abstract

The most massive baryonic component of galaxy clusters is the ‘intracluster medium’ (ICM), a diffuse, hot, weakly magnetized plasma that is most easily observed in the X-ray band. Despite being observed for decades, the macroscopic transport properties of the ICM are still not well constrained. A path to determine macroscopic ICM properties opened up with the discovery of ‘cold fronts’. These were observed as sharp discontinuities in surface brightness and temperature in the ICM, with the property that the denser side of the discontinuity is the colder one. The high spatial resolution of the Chandra X-ray Observatory revealed two puzzles about cold fronts. First, they should be subject to Kelvin–Helmholtz instabilities, yet in many cases they appear relatively smooth and undisturbed. Second, the width of the interface between the two gas phases is typically narrower than the mean free path of the particles in the plasma, indicating negligible thermal conduction. It was thus realized that these special characteristics of cold fronts may be used to probe the properties of the cluster plasma. In this review, we will discuss the recent simulations of cold fronts in galaxy clusters, focusing on those which have attempted to use these features to constrain ICM physics. In particular, we will examine the effects of magnetic fields, viscosity, and thermal conductivity on the stability properties and long-term evolution of cold fronts. We conclude with a discussion on what important questions remain unanswered, and the future role of simulations and the next generation of X-ray observatories.

1 Introduction

Galaxy clusters are the largest gravitationally bound objects in the current universe. They reach total masses of ${\sim}10^{13}{-}10^{15}~M_{\odot }$ , only approximately ${\sim}$ 3–5 % of which is in the member galaxies. The other two major components are the ‘intracluster medium’ (ICM), a hot, ionized plasma that comprises ${\sim}$ 15–17 % of a cluster’s mass and the ‘dark matter’ (DM), which forms the majority of the mass ( ${\sim}$ 80 %) and the precise nature of which is at present unknown but is likely comprised of weakly interacting non-baryonic particles.

Table 1. Relevant length scales in the ICM. ^a

a Keeping the astrophysical length unit kpc as the most practical reference, we express smaller lengths in nanoparsec (npc) and femtoparsec (fpc). We assume ICM conditions of $n_{e}\sim n_{p}=0.01~\text{cm}^{-3}$ , $kT=5~\text{keV}$ and $B=1~{\rm\mu}\text{G}$ , except the last two scales, which are global.

The ICM is one of the largest-scale examples of a space plasma and was discovered in early X-ray observations of clusters (e.g. Cavaliere, Gursky & Tucker 1971; Meekins et al. 1971; Forman et al. 1972; Kellogg et al. 1972; Serlemitsos et al. 1977). It consists of fully ionized hydrogen and helium with traces of highly ionized heavy elements at roughly solar to one-third solar abundances. With electron densities of $n_{e}\sim 10^{-4}{-}10^{-1}~\text{cm}^{-3}$ and temperatures of $T\sim 10^{6}{-}10^{8}~\text{K}$ , the ICM radiates X-rays predominantly through bremsstrahlung, 2-photon emission, photo and dielectronic recombinations and line emission due to collisional excitation processes. Except for the brightest emission lines in the cores of clusters, the gas is optically thin to such X-rays. Inverse Compton scattering of the cosmic microwave background (CMB) off of the thermal electrons in the ICM leads to observable distortions in the CMB spectrum known as the Sunyaev–Zeldovich effect (Sunyaev & Zeldovich 1972). Observations of radio sources in clusters and Faraday rotation measurements reveal that the ICM is weakly magnetized (see Carilli & Taylor 2002; Feretti et al. 2012, for reviews), with magnetic field strengths of $B\sim 0.1{-}10~{\rm\mu}\text{G}$ , corresponding to a magnetic pressure which is roughly a few percent of the thermal pressure.

Table 1 shows the relevant length scales of the ICM for typical values of the density, temperature and magnetic field strength. First, we see that the ICM is weakly collisional, since the typical dynamical scales of interest range from ${\sim}$ 100–1000 kpc and the electron and ion mean free paths are on the order of less than a tenth of a kpc (as in the case of the core of the Virgo cluster) up to tens of kpc (as in cluster outskirts). Under these assumptions, the ICM is typically modelled as a collisional fluid using the equations of hydrodynamics (HD) or magnetohydrodynamics (MHD), although structure on scales below several kpc are observable and thus of interest, and where formally the assumption of collisionality is violated. Second, the relevant scales for plasma physics (e.g. the inertial lengths and Larmor radii of the electrons and ions) are many orders of magnitude smaller than these scales and are not directly observable. It should also be noted that although in one sense the effect of the magnetic field is ‘strong’ due to the fact that the electron and ion Larmor radii are orders of magnitude smaller than their mean free paths, in another sense, the magnetization is ‘weak’ since over most of the cluster volume the magnetic pressure is much weaker than the thermal pressure. As we discuss below, accounting for these effects in physical models of the ICM has important and potentially observable consequences.

The launch of the Chandra X-ray Observatory in 1999 provided a unique window into the physics of the ICM, due to its ${\sim}0.5^{\prime \prime }$ spatial resolution. However, even for deep observations of the diffuse ICM emission, realistically the resolution is reduced to a few arcseconds due to the need for sufficient photon statistics. Nevertheless, for nearby clusters at distances of 16–300 Mpc from us, an angular resolution of 3 arcsec still corresponds to a projected spatial resolution on the sky of approximately 0.25–4 kpc, a scale comparable to the particle mean free path. Among the first cluster discoveries by Chandra were sharp edges in surface brightness, first seen in the clusters A2142 (Markevitch et al. 2000) and A3667 (Vikhlinin, Markevitch & Murray 2001a ). In the latter case, the edge had been previously seen by ROSAT, but not with the same detail (Markevitch, Sarazin & Vikhlinin 1999). Such edges may be expected in colliding clusters, which should drive shocks into the ICM. However, spectroscopic temperature measurements revealed these edges were not shocks, since the denser side of the discontinuity is the cooler one. Because of this difference, these edges were dubbed ‘cold fronts’.

An extended review of the properties and theory of cold fronts in clusters was first given by Markevitch & Vikhlinin (2007, hereafter MV07). Since that review, many more cold fronts have been discovered in galaxy clusters and similar objects, e.g. galaxy groups and elliptical galaxies. All of these systems are similar in the sense that in the context of this review all of them can be thought of as DM potentials (of masses ${\sim}10^{13}{-}10^{15}~M_{\odot }$ ) filled with a hot plasma atmospheres (of temperatures $T\sim 10^{6}{-}10^{8}~\text{K}$ ). Despite this large range in mass, for the remainder of the review we will refer to cold fronts in ‘clusters’ as a shorthand.

1.1 Classes of cold fronts

Figure 1. X-ray surface brightness (a) and spectroscopic temperature (b) maps for the elliptical galaxy NGC 1404, which is falling into the Fornax cluster. The galaxy moves towards the northwest (NW) to the Fornax centre. NGC 1404 is a classic example of a remnant-core cold front, which stretches around the galaxy’s atmosphere from northeast (NE) through NW to southwest (SW). Reproduced from Machacek et al. (2005).

Cosmological structure formation, and hence the growth of galaxy clusters, is still ongoing; thus, clusters are dynamic objects. The motion of gas within or between clusters and the presence of entropy gradients, both of which are conditions that are prevalent in clusters, provide the conditions for the formation of cold fronts. There are two classes of cold fronts, ‘remnant-core’ cold fronts and ‘sloshing’ cold fronts, which were first distinguished by Tittley & Henriksen (2005). In the most general terms, remnant-core cold fronts form when a galaxy, group or cluster moves through a hotter ambient plasma and the resulting head wind strips off the outer layer of its atmosphere. This situation arises when a galaxy or group falls into a larger cluster and even when two clusters are passing through each other. This process has also been described as ‘ram-pressure stripping’ in the context of galaxies moving through the ICM of their host clusters (Gunn & Gott 1972, see also Roediger et al. 2015a and references therein). The formation of the cold front is also aided by the fact that the area-proportional drag force on the gas will decelerate the less dense gas of the core more, bringing the densest gas forward (see figure 23 of MV07 for an illustration). These processes lead to an upstream contact discontinuity between the cooler, denser atmosphere of the infalling object and the hotter, less dense ambient gas. This leading edge is the remnant-core cold front (figure 1). It stretches from the stagnation point to the sides of the infalling atmosphere. Prominent examples are 1E0657-56 (Markevitch et al. 2002, the ‘Bullet Cluster’), the infalling elliptical galaxy NGC 1404 in the Fornax cluster (Dosaj et al. 2002; Scharf, Zurek & Bureau 2005; Machacek et al. 2005, shown in figure 1) and A3667 (Vikhlinin et al. 2001a ; Vikhlinin, Markevitch & Murray 2001b ; Datta et al. 2014). In the literature, this type of cold front has also been referred to as a ‘merger’ cold front.

The second class of cold front, sloshing cold fronts, arises when some process offsets the bulk of the central ICM in a cluster from its hydrostatic equilibrium in the cluster potential. The ICM then slowly oscillates – or sloshes – around its hydrostatic equilibrium configuration (Markevitch et al. 2000; Tittley & Henriksen 2005, AM06). The variation of sloshing frequency with cluster radius leads to the formation of arc-shaped contact discontinuities staggered on opposite sides of the cluster core. They are particularly visible if the cluster is of the ‘cool-core’ variety, where the cluster has been relatively undisturbed for some time and the temperature in the centre has slowly decreased via radiative cooling, causing the central density to go up to maintain hydrostatic equilibrium. In such clusters, the entropy profile is a steep function of radius, increasing the likelihood of producing contact discontinuities. Again, since the denser side of the contact discontinuity is the cooler one, these features are cold fronts. The required global perturbation of the ICM is achieved easily in a minor merger, i.e. the passage of a subcluster through the main cluster causes the ICM offset and subsequent sloshing (AM06), although Churazov et al. (2003) demonstrated that similar effects could be caused by the passage of a sufficiently strong shock. An off-centre passage of a subcluster introduces some additional angular momentum into the ICM, resulting in sloshing cold fronts wrapped around the cluster core in a one-armed spiral like pattern (figure 2). The exact combination of linear sloshing, angular momentum and viewing angle determines whether sloshing cold fronts are arranged in a clear spiral or rather in staggered arcs on opposite sides of the cluster core. Sloshing lasts for several Gyrs, such that the sloshing cluster can appear quite relaxed except for the sloshing signatures. Prominent examples of relaxed clusters with sloshing cold fronts are A2029 (Clarke, Blanton & Sarazin 2004, Paterno-Mahler et al. 2013, see figure 2), Virgo (Simionescu et al. 2010; Werner et al. 2016), A2319 (Govoni et al. 2004; O’Hara, Mohr & Guerrero 2004), A2142 (Markevitch et al. 2000; Rossetti et al. 2013) and A1795 (Markevitch, Vikhlinin & Mazzotta 2001; Ehlert et al. 2015). By some estimates, sloshing cold fronts may exist in the cores of ${\sim}1/2{-}2/3$ of relaxed galaxy clusters (Markevitch, Vikhlinin & Forman 2003; Ghizzardi, Rossetti & Molendi 2010).

Figure 2. Sloshing in A2029. (a) X-ray surface brightness image of the central part of the cluster, showing sloshing cold fronts arranged in a spiral pattern. (b) Surface brightness residual image, highlighting the associated spiral shaped surface brightness enhancement and demonstrating the extent of the sloshing motions out to 400 kpc. Reproduced from Paterno-Mahler et al. (2013).

1.2 Cold fronts and ICM plasma physics

From their initial discovery, it was noted that cold fronts have special properties. First, all of the formation scenarios for cold fronts predict that there should be a velocity shear across the sharp interface. Simple perturbation analysis assuming pure inviscid hydrodynamics shows that such shear flow interfaces should develop Kelvin–Helmholtz instabilities (KHI) on a relatively short time scale (Chandrasekhar 1961), which is a few Myr for cold fronts in the ICM. Thus, a major question in relating cold fronts to ICM properties is whether or not cold fronts are stable or unstable to KHIs.

KHIs at a given shear flow interface can be suppressed in a variety of ways, by, e.g. a sufficiently strong magnetic field aligned with the shear flow interface (Chandrasekhar 1961, Vikhlinin & Markevitch 2002), or by a sufficient ICM viscosity (Roediger et al. 2013b ). In the case of magnetic fields, KHIs are suppressed if the magnetic energy in the component of the field oriented parallel to the interface is roughly in excess of the kinetic energy in the shear flow (Landau & Lifshitz 1960; Vikhlinin & Markevitch 2002). Even if the field is too weak to completely suppress KHI, it may decrease its growth rate. Viscosity suppresses the growth of KHI below a certain length scale by damping out perturbations and smoothing out the velocity gradient at the interface layer by momentum diffusion (Kaiser et al. 2005; Junk et al. 2010; Roediger et al. 2013b ), depending on the Reynolds number of the plasma. It should also be noted in this context that a small initial width of the cold front interface may be sufficient to suppress KHIs for a restricted range of angles along the surface (Churazov & Inogamov 2004), although this is not expected to be true for all observed cold fronts.

At remnant-core cold fronts, the shear flow has a specific structure along the cold front. The shear velocity is zero at the stagnation point, i.e. the most upstream point of the remnant core, and increases along the upstream edge towards the sides of the remnant core. Shear flows also arise along large stretches of sloshing cold fronts. For this reason, we should expect to see evidence of KHIs at the sides of remnant cores and along many sloshing cold fronts, or a washing out of the interface by the KHI. Nevertheless, many (though not all) cold fronts appear relatively undisturbed. This may indicate that the KHIs are being suppressed, or that they are not as easy to recognize.

The second interesting property of cold fronts is the fact of the sharp interface itself. In many cases, the width of the interface is unresolved even by Chandra, and is often smaller than the electron and ion mean free paths. Diffusion of particles and heat should be very efficient across such an interface. In fact, the temperature and density of the gas on either side of the front imply that the thermal conduction time scale across the interface should be very fast, in which case the temperature and density gradients would be smeared out. Nevertheless, sharp cold front interfaces are present in many clusters, indicating that both heat flux and particle diffusion across the interface are both strongly suppressed.

These special properties of cold fronts provide an opportunity to use the combination of numerical simulations with X-ray observations to investigate the physical properties of the cluster plasma. This involves going beyond a simple HD description of the plasma to explicitly including the effects of magnetic fields and transport processes such as viscosity and thermal conduction. Simulations can provide predictions for observable signatures of specific ICM properties, e.g. the appearance or absence of KHIs in X-ray observations or specific shapes for profiles of density and temperature across the fronts. By modelling ICM magnetic fields and transport processes explicitly, simulations may determine what plasma physics of the ICM is responsible for maintaining the smooth and sharp appearance of the cold fronts after their formation and provide a more accurate description of the physics of the ICM. As more simulation works are incorporating these processes to understand their effect on ICM turbulence, their association with AGN feedback and their effects on the thermal state of the cluster as a whole, it is vital to use detailed comparisons between simulations and observations to understand the nature of these processes in the ICM.

In this vein, a number of simulation works have been carried out over the past several years which incorporate various physics beyond pure HD in an attempt to explain the apparent long-term stability properties of cold fronts. These works have provided a number of key insights into ICM physics, but some questions still remain. The purpose of this paper is to review the results of these investigations. These simulations have approached the study of the ICM from the framework of a magnetized fluid, and so the focus of this review is on the constraints that can be placed on parameterized transport coefficients in this framework rather than detailed treatment of the underlying kinetics (see, however, § 4.1.1 for some comments on the latter).

In § 2, we review a number of key observations of nearby cold fronts. In § 3, we review the simulations of cold fronts in clusters which have incorporated physics beyond HD. In § 4 we describe the remaining open questions and the role of future observations and simulations. Lastly, in § 5 we provide a summary.

2 Observations of cold fronts

Since the launch of Chandra, many cold fronts have been discovered in galaxies, galaxy groups and galaxy clusters. In this section we provide an overview of recent X-ray observations of nearby cold fronts that are most relevant to the question of constraining the properties of the cluster plasma. We will first discuss examples of remnant-core cold fronts, followed by sloshing cold fronts.

2.1 NGC 1404/Fornax

An excellent example of a remnant-core cold front is found in the elliptical galaxy NGC 1404 (figure 1), which is falling toward the galaxy NGC 1399 at the centre of the Fornax cluster. Its gas tail, seen already in ROSAT data (Jones et al. 1997; Paolillo et al. 2002), shows that the galaxy is undergoing ram-pressure stripping. Early Chandra observations revealed a sharp leading edge (Dosaj et al. 2002; Machacek et al. 2005; Scharf et al. 2005), and recent deep Chandra observations will allow the width of the remnant-core cold front to be measured in several wedges from the stagnation point towards the sides.

Using stripped galaxies or groups as probes for ICM plasma properties can be complicated by a variety of additional processes superimposed on the infall and stripping, such as central AGN activity as in e.g. M89 (Machacek et al. 2006), a strongly elliptical potential and interactions with other galaxies as in M86 (Randall et al. 2008) or an early stage of infall as for the stripped groups observed in A2142 (Eckert et al. 2014) and A85 (Ichinohe et al. 2015). None of these complications apply to NGC 1404, making it a particularly clean probe of ICM plasma properties.

2.2 A3667

Another nearby example of a remnant-core cold front is in the merging cluster A3667 (figure 3). Line-of-sight velocity measurements of the cluster galaxies show no evidence for multiple components in redshift space, indicating that A3667 is a merger occurring in the plane of the sky (Johnston-Hollitt, Hunstead & Corbett 2008; Owers, Couch & Nulsen 2009). A surface brightness edge in the ICM was identified already in ROSAT and ASCA observations (Knopp, Henry & Briel 1996; Markevitch et al. 1998, 1999), but only Chandra data had sufficient resolution to allow the edge to be identified as a cold front (Vikhlinin et al. 2001a ).

Figure 3. Chandra observation of the cold front in A3667 (from a longer exposure than Vikhlinin et al. 2001a ,   ${\sim}500$  ks). (a) X-ray surface brightness image. The cold front is the sharp edge feature to the SE. (b) X-ray surface brightness profile of the cold front. The blue line shows the best-fit projected spherical density discontinuity, with an infinitely small width. The red line shows the same model smoothed with a Gaussian of width ${\it\sigma}=13$  kpc, corresponding to the collisional mean free path ${\it\lambda}_{in{-}out}$ , which is inconsistent with the data (Reproduced from MV07).

Due to the high spatial resolution of the Chandra data, Vikhlinin et al. (2001a ) were able to show that the surface brightness profiles across the edge were well fit by a projected density discontinuity. They measured an upper limit on the width of the front interface of ${\sim}$ 4 kpc, which is comparable to or smaller than the particle mean free paths on both sides of the cold front, and of those particles moving across the interface (Vikhlinin et al. 2001a ). Particle diffusion across the interface proceeds mostly from the inside of the front (the denser side) to the outside, since the particle flux through the unit area is proportional to $nT^{1/2}$ . Therefore, the most relevant mean free path to consider is ${\it\lambda}_{in\rightarrow out}$ (Spitzer 1962; Vikhlinin et al. 2001a ):

(2.1) $$\begin{eqnarray}\displaystyle {\it\lambda}_{in\rightarrow out}={\it\lambda}_{out}\frac{T_{in}}{T_{out}}\frac{G(1)}{G(\sqrt{T_{in}/T_{out}})}, & & \displaystyle\end{eqnarray}$$

where $G(x)=[{\it\Phi}(x)-x{\it\Phi}^{\prime }(x)]/2x^{2}$ and ${\it\Phi}(x)$ is the error function. For the cold front in A3667, ${\it\lambda}_{in\rightarrow out}\approx 10{-}13$  kpc (MV07). If Coulomb diffusion is not suppressed across the interface, the front width should be at least several times ${\it\lambda}_{in\rightarrow out}$ , yet this is not the case (figure 3), implying that there is a mechanism for such suppression. More recently, Datta et al. (2014) estimated a larger width for the cold front of ${\rm\Delta}R=24\pm 12$  kpc, but their estimate appears to neglect the fact that the projected X-ray emissivity a zero width density jump superficially appears as a surface brightness jump with a finite width. Nevertheless, even with this larger estimate for the width, they were still able to conclude that particle diffusion and thermal conductivity must be strongly suppressed across the interface.

Figure 4. Analysis of cold fronts discovered in A2142 by Markevitch et al. (2000). Panel (a) X-ray image with red regions overlaid showing the regions used for derivation of the temperature profiles shown in (b). Panels (b–e) show the profiles of temperature, surface brightness, density and pressure across the southern (a,b) and northwestern (c–e) edges. Vertical dashed lines show the positions of density jumps. Reproduced from MV07.

The front surface itself appears very smooth, lacking evidence of perturbations within an angular sector of ${\it\phi}\approx 30^{\circ }$ from the direction of the front motion. Outside of this sector, the front becomes less sharp and eventually blends in with the surrounding emission. However, such an interface should be susceptible to the effects of KHI even closer to the axis defined by the direction of motion than this. Vikhlinin & Markevitch (2002) determined that perturbations with wavelengths of ${\approx}20$ –30 kpc would be unstable across the front surface even if the stabilizing effect of the cluster’s gravity is taken into account. However, Churazov & Inogamov (2004) argued that the upstream part of a remnant-core cold front could be stable to KHI if the front has a small initial intrinsic width. Mazzotta, Fusco-Femiano & Vikhlinin (2002) pointed out that the sides of the remnant core in A3667 show a ragged appearance which may indeed be due to KHIs expected at this location.

Assuming that the front is stabilized against KHI by a magnetic field layer, Vikhlinin et al. (2001b ), Vikhlinin & Markevitch (2002) estimated the minimum magnetic field strength that would be required to suppress KHI at small angles from the direction of the front’s motion. They determined that the tangential magnetic field strength at the front surface must be within $6~{\rm\mu}\text{G}<B<14~{\rm\mu}\text{G}$ , which although likely stronger than in most regions of the cluster is still only a small fraction of the thermal pressure, with $p_{m}/p_{th}\sim 0.1{-}0.2$ . Together with the observed sharpness of the front, this is powerful indirect evidence for a magnetic field draping layer at the cold front in A3667.

2.3 A2142

Cold fronts were first discovered by Markevitch et al. (2000) in A2142. They identified two arc-shaped fronts on opposite sides of the cluster core (figure 4). They first proposed that the edges arose from dense subcluster cores which were ram-pressure stripped during a cluster merger – in other words, remnant-core cold fronts. As more cold fronts wrapped around the cores of other clusters were found, it became clear that the cold fronts in A2142 are of the sloshing type (MV07). Rossetti et al. 2013 reported another cold front ${\sim}$ 1 Mpc to the SE, suggesting that sloshing extends to large scales beyond the core region. The hypothesis of sloshing triggered by a minor merger is supported by the observations of Owers, Nulsen & Couch (2011), who identified substructures in the galaxy distribution and their kinematics in A2142 which are candidates for the perturbing subclusters which initiate the sloshing motions.

Ettori & Fabian (2000) noted that the existence of the sharp temperature discontinuities in the A2142 cold fronts indicated that heat flux across the interface must be suppressed by a factor of ${\sim}$ 250–2500, since otherwise they would have been smeared out by thermal conduction on a very short time scale. This was the first indication that cold fronts may be used to constrain transport processes in the ICM.

2.4 Virgo

Figure 5. XMM-Newton EPIC/MOS mosaic image of the central $r\sim 150$  kpc region of the Virgo cluster (Simionescu et al. 2010). The white square indicates the Chandra ACIS-I pointing of Werner et al. (2016) shown in figure 6. Reproduced from Werner et al. (2016).

Figure 6. (a) Chandra image of the Virgo cold front in the 0.8–4.0 keV energy band from Werner et al. (2016). The image was smoothed with a Gaussian function with a $1.5^{\prime \prime }$ window. (b) The residual image, which reveals three X-ray bright quasilinear features separated from each other by ${\sim}15$  kpc in projection. The dark band outside the cold front is an image processing artefact due to the sharp surface brightness discontinuity.

The very nearby Virgo cluster of galaxies (with its central galaxy M87) possesses a pair of sloshing cold fronts (figure 5, Simionescu et al. 2007, 2010). At a distance of only 16.1 Mpc (Tonry et al. 2001), the arcsecond resolution of Chandra corresponds to a projected spatial resolution of 78 pc, making the Virgo cold fronts an ideal target to study cold fronts at the highest resolution. In a recent 500 ks Chandra observation focused on the northern Virgo cold front, Werner et al. (2016) measured the width of the interface at a number of locations along the front surface. In the northwestern part of the front, which appears to be the sharpest, they found an upper limit on the width of the interface to be ${\sim}$ 2.5 kpc at 99 % confidence level, which is approximately 1.5 times the local Coulomb mean free path ${\it\lambda}_{in\rightarrow out}$ for particles diffusing across the interface. In the northern part of the front, the upper limit on the width is found to be somewhat wider, ${\sim}$ 4.3 kpc. They argued that at this location the front may be smeared by the effects of KHI, as had been suggested by the Virgo cluster simulations of Roediger et al. (2013a , hereafter R13; see § 3.4). In several locations, narrow-wedge surface brightness profiles across the Virgo cold front display a multistep structure, a signature for the presence of KHIs, predicted by R13 from mock observations of their simulations. In either case, the Virgo cold front width is much narrower than would be expected if particle diffusion across the front were unsuppressed, which would smooth the front by several ${\it\lambda}_{in\rightarrow out}$ . They also argued from these estimates that thermal conduction is also highly suppressed across the cold front surface.

These observations also revealed a somewhat unexpected feature associated with the cold front. The residual image in figure 6(b) shows the presence of three X-ray bright quasilinear features, extending away from the front surface at an angle which roughly aligns them with the direction of the curvature of the front to the southwest. The features are separated from each other by ${\sim}$ 15 kpc in projection and are not associated with any instrumental artefacts (e.g. chip gaps or readout streaks). Therefore, these features correspond to density enhancements underneath the cold front seen in projection. Werner et al. (2016) argued for an origin of these features that was related to magnetic fields, a possibility we discuss in § 3.3.

2.5 NGC 7618/UGC 12491

The galaxy groups around the elliptical galaxies NGC 7618 and UGC 12491 show classic spiral shaped sloshing cold fronts, arising from their recent mutual encounter (Kraft et al. 2006). However, the cold fronts in these groups are not as smooth as in, e.g. A2142, but display distortions in the form of noses and wings (figure 7, Roediger et al. 2012a ), that resemble the appearance of KHIs as seen in simulations. Cold fronts with similar features exist also in A496 (Dupke, White & Bregman 2007, R12).

Figure 7. The two galaxy groups around NGC 7618 (a) and UGC 12491 (b) passed each other a few 100 Myr ago and caused sloshing in each other’s ICM atmospheres. The spiral shaped cold fronts wrapped around each group centre show distortions in the form of wings and noses, resembling KHIs. Reproduced from Roediger et al. (2012a ).

3 Simulations of cold fronts

Since the discovery of cold fronts, a large number of simulations have been performed to model their formation and evolution in merging galaxy clusters. In this review, we focus on those recent simulations which have attempted to place constraints on the detailed physical properties of the underlying plasma. A table listing the major simulations of cold fronts referenced in this paper, along with the numerical models employed and other properties, can be found in appendix A.

3.1 Basic physics

We will first outline the basic physics underlying these simulations. As previously mentioned, due the fact that the electron and ion mean free paths are smaller than the length scales under consideration, the ICM is typically modelled as a thermal, magnetized plasma in the fluid approximation. However, contrary to the usual MHD ordering of scales, the Larmor radii of the particles are much smaller than their mean free paths, where the Larmor radii are on the order of npc compared to the kpc scale of the mean free path. Therefore, a given particle gyrates many times around the local magnetic field line before encountering another particle. Due to the small Larmor radii and the relative weakness of collisions, we may deduce that the perpendicular and parallel pressures of the particles and the transport coefficients of diffusive processes are anisotropic with respect to the local direction of the magnetic field line. In this situation, the anisotropy of the perpendicular and parallel pressures is produced by the conservation of the first and second adiabatic invariants of the particles.

Under these assumptions, the fluid description of the ICM can be written down in terms of the Braginskii-MHD equations (Braginskii 1965, presented here in conservative form, and in Gaussian units). The equations for conservation of mass, momentum, energy and magnetic field are:

(3.1) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial {\it\rho}}{\partial t}+\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }({\it\rho}\boldsymbol{v})=0, & \displaystyle\end{eqnarray}$$

(3.2) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial ({\it\rho}\boldsymbol{v})}{\partial t}+\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }({\it\rho}\boldsymbol{v}\boldsymbol{v}+\unicode[STIX]{x1D64B})={\it\rho}\boldsymbol{g}, & \displaystyle\end{eqnarray}$$

(3.3) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial ({\it\rho}e)}{\partial t}+\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }[({\it\rho}e\boldsymbol{I}+\unicode[STIX]{x1D64B})\boldsymbol{\cdot }\boldsymbol{v}+\boldsymbol{Q}]={\it\rho}\boldsymbol{g}\boldsymbol{\cdot }\boldsymbol{v}, & \displaystyle\end{eqnarray}$$

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial \boldsymbol{B}}{\partial t}+\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }(\boldsymbol{v}\boldsymbol{B}-\boldsymbol{B}\boldsymbol{v})=0, & \displaystyle\end{eqnarray}$$

where the total (thermal and magnetic) pressure tensor in (3.2) and (3.3) is

(3.5) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D64B}=p_{\bot }\boldsymbol{I}-(p_{\bot }-p_{\Vert })\hat{\boldsymbol{b}}\hat{\boldsymbol{b}}+\frac{B^{2}}{8{\rm\pi}}\boldsymbol{I}-\frac{\boldsymbol{B}\boldsymbol{B}}{4{\rm\pi}}. & & \displaystyle\end{eqnarray}$$

As usual, $\hat{\boldsymbol{b}}=\boldsymbol{B}/B$ is the unit vector in the direction of the local magnetic field and the total thermal pressure satisfies

(3.6) $$\begin{eqnarray}\displaystyle p={\textstyle \frac{2}{3}}p_{\bot }+{\textstyle \frac{1}{3}}p_{\Vert }. & & \displaystyle\end{eqnarray}$$

The terms $B^{2}\boldsymbol{I}/8{\rm\pi}$ and $-\boldsymbol{B}\boldsymbol{B}/4{\rm\pi}$ in (3.5) are the familiar ‘magnetic pressure’ and ‘magnetic tension’, respectively. In (3.3), the total energy per unit volume is

(3.7) $$\begin{eqnarray}\displaystyle {\it\rho}e=\frac{{\it\rho}v^{2}}{2}+{\it\rho}{\it\epsilon}+\frac{B^{2}}{8{\rm\pi}}. & & \displaystyle\end{eqnarray}$$

In these equations, ${\it\rho}$ is the gas density, $\boldsymbol{v}$ is the fluid velocity, ${\it\epsilon}$ is the gas thermal energy per unit volume, $\boldsymbol{B}$ is the magnetic field vector, $\boldsymbol{g}$ is the gravitational acceleration, $\boldsymbol{Q}$ is the heat flux vector and $\boldsymbol{I}$ is the unit dyad. Generally, an adiabatic equation of state of an ideal gas $p=({\it\gamma}-1){\it\rho}{\it\epsilon}$ with ${\it\gamma}=5/3$ is assumed, along with equal electron and ion temperatures, and primordial abundances of H and He with trace amounts of metals, for an average molecular weight of $\bar{A}\approx 0.6$ .

Most of the astrophysical literature on galaxy clusters does not present the MHD equations with a pressure anisotropy, but instead describes this effect in terms of a viscous flux. We will now show how these two descriptions are equivalent. Differences in the two components of the thermal pressure arise from the conservation of the first and second adiabatic invariants for each particle on time scales much greater than the inverse of the ion gyrofrequency, ${\it\Omega}_{g}^{-1}$ (Chew, Goldberger & Low 1956). When the ion–ion collision frequency ${\it\nu}_{ii}$ is much larger than the rates of change of all fields, an equation for the pressure anisotropy can be obtained by balancing its production by adiabatic invariance with its relaxation via collisions (cf. Schekochihin et al. 2005):

(3.8) $$\begin{eqnarray}\displaystyle p_{\bot }-p_{\Vert }=0.960\frac{p_{i}}{{\it\nu}_{ii}}\frac{\text{d}}{\text{d}t}\ln \frac{B^{3}}{{\it\rho}^{2}}, & & \displaystyle\end{eqnarray}$$

where $p_{i}$ is the thermal pressure of the ions (due to their much smaller mass, the effect of the electrons on the viscosity may be neglected).

The factor in front of the time derivative on the right-hand side of (3.8) is simply the dynamic viscosity coefficient for the ions, or the ‘Spitzer’ viscosity, which is given by (Spitzer 1962; Braginskii 1965; Sarazin 1988):

(3.9) $$\begin{eqnarray}\displaystyle {\it\mu}_{Sp} & = & \displaystyle 0.960\frac{p_{i}}{{\it\nu}_{ii}}\nonumber\\ \displaystyle & {\approx} & \displaystyle 2.2\times 10^{-15}\frac{T^{5/2}}{\ln {\it\Lambda}_{i}}~\text{g}~\text{cm}^{-1}~\text{s}^{-1},\end{eqnarray}$$

where the temperature $T$ is in Kelvin, ${\it\nu}_{ii}$ is the ion–ion collision frequency and $\ln {\it\Lambda}_{i}$ is the ion Coulomb logarithm, which is a weak function of ${\it\rho}$ and $T$ ; for conditions in the ICM, it may be approximated as $\ln {\it\Lambda}_{i}\approx 40$ . Using (3.1) and (3.4) to replace the time derivatives of density and magnetic field strength with velocity gradients, equation (3.8) may be written as

(3.10) $$\begin{eqnarray}\displaystyle p_{\bot }-p_{\Vert }={\it\mu}_{Sp}[3\hat{\boldsymbol{b}}\boldsymbol{\cdot }(\hat{\boldsymbol{b}}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}\boldsymbol{v})-\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\boldsymbol{v}]. & & \displaystyle\end{eqnarray}$$

which, upon plugging back into (3.2), shows that the effect of the pressure anisotropy manifests itself as an anisotropic viscous flux.

Next, we consider the effect of thermal conduction. Due to the effect of the small Larmor radii of the electrons, heat only flows along the field lines and the only component of the temperature gradient which affects the heat flux is that which is parallel to the field direction:

(3.11) $$\begin{eqnarray}\displaystyle \boldsymbol{Q}=-{\it\kappa}_{Sp}\hat{\boldsymbol{b}}\hat{\boldsymbol{b}}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}T, & & \displaystyle\end{eqnarray}$$

where ${\it\kappa}_{Sp}$ is the familiar ‘Spitzer’ thermal conductivity coefficient (Spitzer 1962; Braginskii 1965; Sarazin 1988):

(3.12) $$\begin{eqnarray}\displaystyle {\it\kappa}_{Sp} & = & \displaystyle \frac{3.2k_{B}p_{e}}{m_{e}{\it\nu}_{ee}}\nonumber\\ \displaystyle & {\approx} & \displaystyle 1.84\times 10^{-5}\frac{T^{5/2}}{\ln {\it\Lambda}_{e}}~\text{erg}~\text{cm}^{-1}\text{s}^{-1}~\text{K}^{-1},\end{eqnarray}$$

where $p_{e}$ is the thermal pressure of the electrons, ${\it\nu}_{ee}$ is the electron–electron collision frequency and $\ln {\it\Lambda}_{e}$ is the electron Coulomb logarithm, which may also may be approximated as $\ln {\it\Lambda}_{e}\approx 40$ . Due to their much larger mass, the heat flux due to ion collisions is negligible, so it is neglected.

So far, most simulations of the ICM have not worked with the full set of these equations, but rather in certain physically motivated limits (though Dong & Stone 2009, Kunz et al. 2012, Parrish et al. 2012, Suzuki et al. 2013 and ZuHone et al. 2015 are notable exceptions). In particular, there have been two approximations of these equations in common use. In the limit that the pressure anisotropy $p_{\bot }-p_{\Vert }$ is small (in other words, collisions are relatively strong and the viscosity is low), we recover standard MHD (or simply HD if it is assumed that the magnetic field can also be neglected). In this review we will describe a number of simulations that work in this MHD limit, which emphasize the effect of the magnetic field on cold fronts. Since the Prandtl number of the Spitzer description of the ICM is very small, the thermal diffusivity dominates over the momentum diffusivity:

(3.13) $$\begin{eqnarray}\displaystyle Pr\equiv \frac{{\it\nu}}{{\it\chi}}=0.5\frac{\ln {\it\Lambda}_{e}}{\ln {\it\Lambda}_{i}}\left(\frac{2m_{e}}{m_{i}}\right)^{1/2}\simeq 0.02. & & \displaystyle\end{eqnarray}$$

Therefore, another limit of these equations considers MHD with thermal conductivity added. We will also consider a number of simulations that have operated in this limit.

A second limit of these equations that is often modelled is that in which the magnetic field is very weak and turbulent throughout the volume of the cluster. In this case, although the effects of magnetic pressure and tension are neglected, the anisotropy of the heat and momentum fluxes manifest as an effective suppression of the transport coefficients when averaged over a small volume, assuming the magnetic field is isotropically tangled. We may then rewrite (3.2) and (3.3) as:

(3.14) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial ({\it\rho}\boldsymbol{v})}{\partial t}+\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }({\it\rho}\boldsymbol{v}\boldsymbol{v}+p\boldsymbol{I}+\unicode[STIX]{x1D6D1})={\it\rho}\boldsymbol{g}, & \displaystyle\end{eqnarray}$$

(3.15) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial ({\it\rho}e)}{\partial t}+\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\{[({\it\rho}e+p)\boldsymbol{I}+\unicode[STIX]{x1D6D1}]\boldsymbol{\cdot }\boldsymbol{v}+\boldsymbol{Q}\}={\it\rho}\boldsymbol{g}\boldsymbol{\cdot }\boldsymbol{v}, & \displaystyle\end{eqnarray}$$

and the conductive and viscous flux terms may be written as

(3.16) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{Q}=-f_{c}{\it\kappa}_{Sp}\boldsymbol{{\rm\nabla}}T, & \displaystyle\end{eqnarray}$$

(3.17) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6D1}=-f_{v}{\it\mu}_{Sp}\boldsymbol{{\rm\nabla}}\boldsymbol{v}, & \displaystyle\end{eqnarray}$$

where the ‘suppression factors’ $f_{c}$ and $f_{v}$ account for the reduced conductivity and viscosity due to the averaging over the random direction of the magnetic field, or other microscale processes that suppress the effective conductivity and/or viscosity below the Spitzer value. This HD limit of the equations is the approximation that has been assumed in a number of simulation works, including some important ones considered in this review. Consideration of this limit is relevant, since the finite Larmor radius effects on the transport coefficients may often only reveal themselves in a spatially averaged way due to projection effects and the finite resolution of X-ray observations.

It is also important to note that the Braginskii MHD equations are ill-posed in situations where the absolute value of the pressure anisotropy exceeds certain bounds (on the order of the magnetic pressure), and in this regime finite Larmor radius effects become important. We will describe these issues further and the consequences of these effects for simulations of the ICM in § 4.1.1.

3.2 Results from simulations of the formation and global evolution of cold fronts

The basic physics of the formation and global evolution of cold fronts can be understood in the pure hydrodynamics framework for both the remnant core and the sloshing variety. The results of these studies, in particular their points of divergence from observed cold fronts, indicated the need for additional physics beyond HD that is discussed in subsequent sections.

Stripping of gas from one gravitational potential falling into another, usually larger, gas-filled potential has been simulated in various contexts, ranging from stripping of gas clouds inside galaxy halos (e.g. Murray et al. 1993; Pittard, Hartquist & Falle 2010), to elliptical galaxies in clusters (e.g. Lea & De Young 1976; Takeda, Nulsen & Fabian 1984; Portnoy, Pistinner & Shaviv 1993; Stevens, Acreman & Ponman 1999; Acreman et al. 2003; Takizawa 2005, etc.) to mergers between clusters (Roettiger, Loken & Burns 1997; Ricker & Sarazin 2001; Ritchie & Thomas 2002; Poole et al. 2006, 2007, 2008). Such simulations routinely produce a remnant-core cold front at the upstream side of the stripped atmosphere. Strictly speaking, several of these idealized simulations do not prove that gas stripping forms an upstream cold front because they already start with a cold front, a contact discontinuity between the cooler to-be-stripped atmosphere and the hotter ambient medium. However, simulations of high mass ratio cluster mergers such as Poole et al. (2006) and ZuHone (2011) demonstrate that remnant merger core cold fronts form readily in such situations. Heinz et al. (2003) pointed out an additional effect that may help to form and maintain the upstream cold front. They suggest that the flow of ambient gas around the stripped atmosphere causes a slow forward motion of the gas inside the stripped potential, bringing cooler gas from the centre of the infalling atmosphere towards the upstream edge. However, this does not seem to be a universal process and may depend on the depth and steepness infalling potential as well as on the time scale available for the stripping during the cluster infall (Roediger et al. 2015a , hereafter R15A). All purely hydrodynamical simulations of the formation of such remnant-core cold fronts that have sufficient resolution predict the growth of KHIs along the sides of the stripped atmosphere.

Ascasibar & Markevitch (2006, hereafter AM06) carried out simulations of sloshing cold fronts with the Gadget SPH code, and first described in great detail the formation and evolution of these features. Roediger et al. (2011) and Roediger et al. (2012b , hereafter R12) demonstrated that the positions of the sloshing cold fronts and the contrasts in density and temperature across them can be used to determine the merger history of a given cluster. In the simulations of AM06, the cold fronts which formed were stable against KHI in spite of the presence of shear flows along them, since the SPH method employed at the time was unable to resolve the KHIs. Simulations with grid codes, however, routinely predict KHIs at sloshing cold fronts. ZuHone, Markevitch & Johnson (2010) performed sloshing simulations using identical initial conditions as AM06, except that in their case the simulations were carried out using the FLASH grid-based code. Although the same cold fronts were formed in both sets of simulations, the FLASH simulations produced KHI at the interfaces (see figure 9). KHIs also developed naturally in the sloshing simulations for the Virgo cluster and A496 (Roediger et al. 2011, R12, R13). Their appearance in synthetic X-ray observations is discussed in § 3.4.

Many cold fronts are believed to be free of KHIs. The following sections discuss ICM plasma properties that could suppress them.

Figure 8. Formation of a merger-remnant cold front in a purely hydrodynamic simulation from Heinz et al. (2003). Slices are of gas entropy, with velocity vectors overlaid. Note the rapid development of KHI.

Figure 9. Sloshing cold fronts (shown with slices of gas temperature) from an SPH simulation by AM06 (a) and a grid-based simulation of the same cluster set-up by ZuHone et al. (2010) (b). The grid-based simulation resolves KHI that are unresolved by the SPH simulation.

3.3 Magnetic fields

Throughout most of the volume of a galaxy cluster, the magnetic field strength is dynamically negligible. The influence of the magnetic field is typically parameterized by the plasma parameter ${\it\beta}=p_{th}/p_{B}$ , the ratio of the thermal pressure $p_{th}$ and the magnetic pressure $p_{B}=B^{2}/8{\rm\pi}$ . Faraday rotation measurements (RM) constrain the value of ${\it\beta}$ to be as low as ${\sim}$ 50 or as high as several hundred, with the magnetic field strength decreasing in radius from the cluster centre, and turbulent on scales ranging from 100 kpc down to 100 pc, as determined from RM maps (see Carilli & Taylor 2002; Eilek & Owen 2002; Govoni et al. 2006; Laing et al. 2008; Bonafede et al. 2010; Guidetti et al. 2008, 2010; Feretti et al. 2012).

However, the existence of the discontinuous changes in density and temperature observed in cold fronts leads to the expectation that the magnetic field in these regions may not be typical of the cluster as a whole. This was first shown for the case of remnant-core cold fronts using analytic arguments by Lyutikov (2006). As a cold front moves through the magnetized medium, the field is amplified and aligned with the cold front surface as the flow stretches and wraps or ‘drapes’ the field lines around the interface. This effect occurs whether the cold front is moving subsonically or supersonically, so long as the velocity is super-Alfvénic, which is easily satisfied for the weak magnetic fields in the ICM. The draping layer that forms is very thin compared to the radius of curvature of the cold front, and has a magnetic pressure nearly equal to the surrounding pressure.

Using an analytical stability analysis and simulations, Dursi (2007) showed that a thin magnetized layer located at a cold front-like interface could stabilize interfaces even against long-wavelength perturbations of KHIs and Rayleigh–Taylor instabilities (RTIs). In particular, they showed that given a magnetic field layer with thickness $l$ , to stabilize wavelengths ${\it\lambda}\approx 10l$ requires that the Alfvén speed in the layer is comparable to the relevant destabilizing velocity – namely, the shear velocity in the case of KHI. This was followed up with simulations by Dursi & Pfrommer (2008), who performed a detailed study of the propagation of cold clumps of gas through a magnetized medium. They established that the magnetic field strength in the draping layer that forms is set by a competition between the ‘plowing up’ of field lines around the cold front and the slipping of field lines around the sides of the front. An animation showing this effect can be seen in figure 10. Another early set of MHD simulations of the formation of remnant-core cold fronts that were explicitly used to investigate the effects of the draping layers were carried out by Asai, Fukuda & Matsumoto (2004, 2005, 2007, hereafter together referred to as ‘Asai-MHD’). These simulations used a similar set-up to the simulations of Dursi & Pfrommer (2008), but also included gravity and had various initial configurations for the magnetic field lines ranging from uniform to turbulent. The Asai et al. (2004) simulations were carried out in two dimensions (2-D), and the Asai et al. (2005) and Asai et al. (2007) simulations were carried out in 3-D.

Figure 10. Animation of the draping of magnetic field lines around a cold cloud in the simulations of Dursi & Pfrommer (2008). The cutting plane is coloured by magnetic energy, as are the magnetic field lines.

Figure 11. Slices through the temperature (a), temperature gradient (b) and magnetic field strength (c) for the formation of a cold front in a magnetically turbulent medium from Asai et al. (2007). The solid curves in (a) show the contours of the magnetic field strength. Arrows in (a) and (c) show the velocity vectors and arrows in (b) show the gradients of temperature. A length unit of ‘1’ corresponds to 250 kpc.

In all of these simulations, the magnetic draping layers that form suppress KHI at the interface (see figure 11 for examples of the temperature and magnetic field distributions from one of the Asai et al. (2007) simulations), as predicted by Lyutikov (2006). Importantly, these works showed that such layers also formed regardless of the initial field configuration, whether uniform or turbulent. They also find that the convergence of the flow behind the subcluster amplifies the field in this location as well. The Asai-MHD simulations also included anisotropic thermal conduction, which will be discussed in § 3.5.

The MHD cluster merger simulations of Takizawa (2008) also produced magnetic draping layers at the cold front surfaces, but these simulations had low resolution and were not used to investigate the properties of cold front stability. It should also be noted that magnetic draping effects may also be responsible for the stability of AGN-blown bubbles; which has been demonstrated by a number of simulation works (Robinson et al. 2004; Jones & De Young 2005; Ruszkowski et al. 2007; Dong & Stone 2009; O’Neill, De Young & Jones 2009). In particular, Ruszkowski et al. (2007) pointed out that the ability of the magnetic fields to stabilize a buoyantly rising bubble depends on the coherence length of the fields. If the latter is smaller than the bubble radius, no useful draping layer can form, and the bubble dissolves by KHIs and RTIs. The simulations of Asai et al. (2007) indicate that for remnant-core cold fronts a draping layer forms quickly also for a turbulent magnetic field with a coherence length smaller than the remnant atmosphere.

Figure 12. Slices of gas temperature (a) and magnetic field strength (b) in an MHD simulation of gas sloshing. Clicking on the figure shows an animation of the development of the cold fronts and the magnetized layers over ${\sim}$ 2 Gyr. Along the front surfaces, the magnetic field strength is increased and the magnetic field within the cold fronts remains high as the cold fronts expand.

Figure 13. Slices through the plasma ${\it\beta}$ for the simulations from ZML11 with varying initial ${\it\beta}$ . Each panel is 500 kpc on a side. Major tick marks indicate 100 kpc distances.

Figure 14. Slices through the gas temperature for the simulations from ZML11 with varying initial ${\it\beta}$ . Each panel is 500 kpc on a side. Major tick marks indicate 100 kpc distances.

Figure 15. Comparison of simulations of gas sloshing with and without magnetic fields. Panels show slices of gas temperature. In the magnetized simulation, ${\it\beta}\sim 100$ . Clicking on the figure shows an animation of the sloshing cold fronts developing over ${\sim}$ 2 Gyr. KHI develop more readily in the unmagnetized simulation.

Though similar magnetic field layers to the remnant-core cold fronts are also expected to appear in sloshing cold fronts, they are also somewhat different. In the remnant-core case, the movement of the cold, dense gas through the surrounding ICM builds up a layer of amplified magnetic field on the outside of the interface (Lyutikov 2006). In the case of sloshing cold fronts, in simulations we often find a velocity shear across the front surface due to the velocity difference between the fast moving cold gas beneath the front, moving tangential to the interface, and the slowly moving hot gas above the front (see also Keshet 2012). This shear amplifies the magnetic field via line freezing, $\text{d}(\boldsymbol{B}/{\it\rho})/\text{d}t\simeq [(\boldsymbol{B}/{\it\rho})\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}]\boldsymbol{v}$ , on the inside of the interface (Keshet et al. 2010). Figure 12 shows an animation of sloshing cold fronts using slices of gas density and magnetic field strength from an MHD simulation.

The first simulations of sloshing cold fronts to include magnetic fields were those of ZuHone, Markevitch & Lee (2011, hereafter ZML11). They performed a parameter space exploration over the initial magnetic field configuration. The initial magnetic field was either turbulent or purely tangential, and the average field strength and its coherence length scale were varied. The former was set by enforcing the ${\it\beta}$ parameter to be the same everywhere. This produces a declining magnetic field radial profile which is consistent with existing observations (Bonafede et al. 2010).

In these simulations, the shear flows initiated by the gas sloshing produce highly magnetized (low- ${\it\beta}$ ) layers just underneath the cold front surfaces. Figure 13 shows slices of the plasma ${\it\beta}$ for the simulations with varying magnetic field strength. These layers are produced in all of the simulations, but the strength of the layer depends on the strength of the initial field. For example, figure 13(d) shows a magnetic layer in the cold front to the southwest which is nearly in equipartition with the thermal pressure, but the corresponding layers in the simulations with initial weaker field (increasing ${\it\beta}$ ) are not nearly as strong. Roughly, the magnetic pressure in the layers is increased by an order of magnitude over the initial field. These strongly magnetized layers are transient, only occurring directly underneath the front surface as the fronts expand. Within the increasing volume defined by the expanding fronts, the magnetic field that develops in the core behind the front is stronger than the initial field by a factor of a few. This field is typically very turbulent, except exactly below the front surface.

Figure 14 shows slices through the gas temperature for the same simulations at the same epoch. For simulations with stronger initial field, the magnetic field in the layers becomes stronger, and more capable of suppressing KHI. Figure 15 shows a side-by-side animation of gas sloshing in an unmagnetized versus a magnetized medium (ZuHone et al., in prep), showing again that the magnetic field suppresses the growth of KHI considerably compared to a case where there is no magnetization. For essentially complete suppression of KHI, the initial magnetic pressure must be at least roughly 1 % of the thermal pressure, or ${\it\beta}\lesssim 100$ . In contrast, ZML11 did not find that the suppression of KHI has a strong dependence on the coherence length of the magnetic field fluctuations or the magnetic field geometry (whether turbulent or tangential), provided that the field strength was the same (figure 22 of ZML11).

However, they did find that these conclusions were somewhat dependent on resolution – for a simulation with initial ${\it\beta}=100$ and a finest resolution of 2 kpc, KHI is highly suppressed, but for an otherwise identical simulation with a finest resolution of 1 kpc, small-wavelength KHI were able to grow (figures 38 and 39 of ZML11). This is likely due to a confluence of two factors: the decreased numerical viscosity associated with higher spatial resolution, and the decrease in the wavelength at which KHI can be suppressed due to the thinning of the magnetized layer which occurs with decreased resolution (see Dursi 2007). Given these facts, and the uncertainties associated with the magnetic field strength in clusters, it is uncertain that magnetic fields provide the observed suppression of KHI in all clusters where such suppression is observed.

Magnetic fields may also explain the peculiar linear features seen in the observations of the Virgo cold front by Werner et al. (2016) (figure 6). Figure 16 shows slices of density, temperature and magnetic field strength, and residuals of surface brightness for a simulation designed to reproduce the Virgo cluster cold fronts (based on the simulations from ZuHone et al. (2015)) using the initial conditions from R12. The initial magnetic field strength for the simulation was set by ${\it\beta}\sim 100$ . The images show a region of the cluster approximately located at the same place where the observations of Werner et al. (2016) were taken. They show narrow, quasilinear, dense, cold features, which correlate with narrow channels of weak magnetic field in between wide bands of strong magnetic field which are less dense and hotter. The dense features produce surface brightness features similar to those seen in figure 6.

Since the magnetic field in the simulation is turbulent before the sloshing begins, the local magnetic field strength will vary from place to place. When the sloshing motions begin, the field is stretched and amplified in the direction of the sloshing, which erases much of the small-scale structure in the field. However, we may expect that perpendicular to the direction of motion, some memory of the variations in magnetic field strength will be preserved. Since the simulations demonstrate that the final amplified field strength does depend to some extent on the initial field strength, these variations will produce regions of slightly weaker magnetic field sandwiched between regions of stronger field. The total pressure across these bands (magnetic and thermal) will be continuous.

Figure 16. Sloshing in a Virgo-like simulated cluster from Werner et al. (2016). (a) Slice through the centre plane of electron density. (b) Slice of gas temperature. (c) Slice of magnetic field strength. (d) Projected surface brightness residuals.

The narrow regions of slightly weaker field are correlated with regions of higher density and lower temperature compared to the surrounding regions with higher magnetic field strength. The anticorrelation between density and temperature is somewhat peculiar, as the simulations exhibiting these features are adiabatic (without sources of heating or cooling), and for an individual fluid element the entropy should be constant. This implies that a decrease in the gas thermal pressure due to an increase the magnetic pressure in the wide bands should correspond to a decrease in both the gas density and temperature in these bands, yet the temperature is higher. This may be explained by low-entropy, low-temperature high-density gas being advected from the cluster centre to larger radii by the sloshing, but it is not completely clear why such gas should be associated with regions of low magnetic field strength. In real clusters, the density enhancement and temperature decrement will be enhanced by gas cooling. A more careful investigation of this phenomenon from simulations is required.

3.4 Cold fronts and viscosity

Small-scale KHIs will be suppressed if the shear velocity across the interface does not change discontinuously, but is smoothed over a certain thickness $\pm d$ . In this case, KHIs of length scale ${\lesssim}10d$ are suppressed (Chandrasekhar 1961, § 102). In a viscous fluid, momentum diffusion across a shear flow interface will establish exactly such a smoothed out shear flow layer. Thus, we can naturally expect that with increasing viscosity KHIs are suppressed at ever larger length scales. Assuming an isotropic viscosity, Roediger et al. (2013b ) showed that the suppression of KHIs occurs roughly at Reynolds numbers $Re=LU/{\it\nu}$ smaller than 30–100, where $L$ is the length scale of the KHI, $U$ the shear velocity and ${\it\nu}={\it\mu}/{\it\rho}$ is the kinematic viscosity. The detailed value of the critical Reynolds number depends on the nature of viscosity (e.g. constant versus temperature dependent ${\it\nu}$ ) and density contrast across the shear interface.

At full Spitzer viscosity (3.9), the critical Reynolds number for suppressing KHIs is easily reached for both remnant core and sloshing cold fronts. Furthermore, a sufficient ICM viscosity prevents the mixing of gas stripped from the infalling subcluster or galaxy with the ambient ICM in the wake (Roediger et al. 2015b , hereafter R15B). Recent deep Chandra and XMM-Newton observations resolve the structure of the cold fronts and of galaxy or subcluster wakes. However, the interpretation of the observations on their own is not straightforward because the appearance of individual objects depends on several other parameters such as the infall and/or merger history, viewing geometry, initial gas contents, or internal perturbations due to, e.g. AGN activity.

R15A and R15B simulated the infall of the elliptical galaxy M89 into the Virgo cluster and followed the progressive gas stripping process as the galaxy is crossing its host cluster. The authors showed that progressive gas stripping of a galaxy or subcluster leads to the formation of a ‘remnant tail’ of the atmosphere (figure 17). Gas from the galaxy’s atmosphere is removed predominantly upstream and from its sides, but the downstream part of the atmosphere is shielded from the ICM head wind (see animation in figure 17). Consequently, the stripped atmosphere develops a head–tail structure, where the tail is simply the unstripped downstream atmosphere. Making the analogy to the flow around a blunt body, the wake of the galaxy starts only downstream of the remnant tail. These global flow patterns, i.e. the upstream edge or cold front of the remnant core, its remnant tail, and the downstream wake in the ICM are independent of viscosity.

R15B focused on observable signatures of an isotropic but Spitzer-like temperature-dependent viscosity (figure 18). In inviscid stripping, KHIs at the sides of the remnant atmosphere take the shape of horns or wings that are observable in X-ray images, and mixing of stripped colder galactic gas leads to fainter, warmer wakes. Interestingly, already at a level of $1/10\text{th}$ of the Spitzer viscosity, KHIs are suppressed, horns or wings at the sides of the remnant core are absent and the stripped-off colder and denser galactic gas leads to bright, cool wakes far behind the galaxy, implying that the viscosity is low in systems such as M89 which possess evidence of ‘horned’ cold fronts (figure 18).

Figure 17. Snapshots and animation from simulations of progressive gas stripping of an elliptical galaxy, from R15A and R15B. (a,b) and the associated animation show density slices through the stripped galaxy for inviscid stripping and with a viscosity of 0.1 of the Spitzer level. During the infall into the host cluster, characteristic flow regions develop, which are labelled in (c). A remnant-core cold front exists at the upstream side, stretching around to the sides of the atmosphere. Here the interface between the galactic gas and the ICM becomes KHI unstable at sufficiently low viscosity. The downstream atmosphere is shielded from the ambient flow for a long time and takes the appearance of a remnant tail. Only in the wake can stripped gas mix with the ambient gas, unless the mixing is suppressed by viscosity.

Figure 18. Chandra X-ray image of the stripped elliptical galaxy M89 in the Virgo cluster (a), compared to mock X-ray images of inviscid and a viscous stripping simulations for the same galaxy (b,c). The observation and simulations clearly show the remnant-core cold front, labelled ‘upstream edge’ here. Comparing the simulations and the observations identifies the observed near tail of M89 as the remnant tail, and its wake is only found at the edge of the observed field of view. In the inviscid simulation, KHIs peel off galactic gas from the sides of the remnant core. If the Virgo ICM has a viscosity roughly 0.1 of the Spitzer viscosity, KHIs at the sides of M89 would be suppressed, and it should have a cold bright wake out to 10 times the remnant-core radius downstream of the galaxy. Reproduced from R15B.

Figure 19. The effects of Braginskii viscosity on remnant-core cold fronts, from S13. (a,b) Slices through the gas temperature in the $y$ – $z$ plane from the inviscid (a,c) and Braginskii viscosity (b,d) simulations. In these units, $1=250$  kpc. (c,d) Profiles of the gas density and temperature across the cold front interface, from the inviscid (a,c) and Braginskii viscosity (b,d) simulations. The cold front interface is located at $z\approx -0.5$ ( $-125~\text{kpc}$ ). The growth of KHI in both simulations is similar, demonstrating that in this scenario Braginskii viscosity has only a marginal effect on suppressing KHI.

Figure 20. Snapshot and animation from simulations of inviscid and viscid sloshing in Virgo cluster, from R13. (b,c) Shows temperature slices through the cluster centre in the orbital plane of the perturber. Clicking on the images shows an animation of the final 0.5 Gyr of the simulations. Up to this point both simulations had been run at zero viscosity. A viscosity at 0.1 of the Spitzer level erases pre-existing KHIs along the northern front, which continued to grow and form in the inviscid run. (a) Shows a mock X-ray image of the inviscid sloshing, zoomed in on the northern cold front. In X-ray images, KHIs appear as multiple edges and give the cold front a ragged appearance. The multiple edges also appear in narrow-wedge surface brightness profiles across the cold front (three edges are marked in the profile in d). At a sufficient viscosity, however, the cold front appears as a single edge in the surface brightness profile (e).

Suzuki et al. (2013, hereafter S13) carried out the first investigation of remnant-core cold front evolution under the influence of an anisotropic (Braginskii) viscosity. They used weak magnetic fields so that their effect on KHI suppression would be negligible. The magnetic fields were uniform in direction, oriented either perpendicular to the axis of motion of the cold front or inclined by a $45^{\circ }$ angle. They showed that anisotropic viscosity is not as efficient as an isotropic viscosity at suppressing KHIs because the damping effect of viscosity is reduced by the dependence on the field line direction (see figure 19). In the plane parallel to the initial field line direction, KHI grows very fast, since the field lines are perpendicular to the velocity gradient at the interface. Conversely, in the plane perpendicular to the initial field line direction, the growth rate of KHI was somewhat reduced, due to the fact that the field line had components parallel to the velocity gradient. A similar dependence on the initial field direction on the growth of KHI at the surfaces of AGN-blown bubbles was found by Dong & Stone (2009), who also included Braginskii viscosity in their simulations. This points to the need for future simulations to include more realistic magnetic field configurations.

The first simulation of sloshing cold fronts with an isotropic viscosity was carried out by ZuHone et al. (2010), who studied how sloshing motions in a galaxy cluster provide a source of heat to offset cooling via mixing of higher entropy gas from larger radii with the lower entropy gas from the core. They compared simulations without viscosity and with a constant kinematic viscosity which roughly corresponded to the Spitzer value in the core region. They found that gas mixing is very efficient if the ICM is inviscid, whereas it is negligible if the plasma is viscous. They noted that the viscosity greatly suppressed KHIs along the cold fronts.

Figure 21. Simulations of the Virgo cold fronts with different prescriptions for viscosity, including anisotropic (Braginskii) viscosity. From top left counterclockwise, simulations are: inviscid, Spitzer anisotropic viscosity, Spitzer isotropic viscosity, 1/10th Spitzer isotropic viscosity. (a) Slices through the gas temperature in keV. (b) Simulated 300 ks Chandra observations. Cold fronts with an appearance that has been visibly modified by viscosity are marked with arrows; for a close up of the NW cold front see figure 22.

R13 showed that the presence of KHIs at sloshing cold fronts leads to characteristic observable signatures such as a ragged appearance of the surface brightness edge and a multistep structure of the cold front in surface brightness profiles (figure 20). In contrast, suppressing KHIs with an isotropic, temperature-dependent Spitzer-like viscosity with a suppression factor of $f_{v}\sim 0.1$ , led to smooth sloshing cold fronts and a clean, single step edge in surface brightness profiles. Furthermore, the region below the cold fronts shows much stronger surface brightness fluctuations in the inviscid case than in the viscous one. The deep Virgo observations of Werner et al. (2016) show both surface brightness fluctuations below the cold front and a multistep structure of the cold front itself, which may indicate the presence of KHIs and a strongly suppressed ICM viscosity in the Virgo cluster. Distorted or ‘ragged’ sloshing cold fronts have also been identified in the merging groups NGC 7618/UGC 12491 (Roediger et al. 2012a , figure 7) and in A496 (Dupke et al. 2007, R12).

Figure 22. Close up of the NW cold front in figure 21 for simulations with inviscid and anisotropic (Braginskii) viscosity. (a–c) Show the inviscid simulation, whereas the (d–f) show the anisotropic viscosity simulation. (a–f) The panels show the simulated 300 ks exposure, simulated 300 ks residuals, and slice through the magnetic field strength. In the inviscid simulation, KHI distort the cold front and tangle the magnetic field lines.

The first simulations of sloshing cold fronts to employ anisotropic viscosity were those of ZuHone et al. (2015, hereafter Z15). They compared magnetized simulations with anisotropic viscosity to those with isotropic viscosity with various suppression factors $f_{v}$ . Similar to S13, they used weak magnetic fields (with ${\it\beta}\sim 1000$ ) to isolate the effect of viscosity on KHI. Using synthetic X-ray observations, they found that anisotropic viscosity is capable of suppressing KHI at cold front surfaces to a degree that is qualitatively similar to an isotropic viscosity with a Spitzer fraction of $f_{v}=0.1$ (see figures 21 and 22). They also compared simulations with isotropic viscosity with and without magnetic fields, and noted that even their weak seed field has effects on cold front stability that can yield observational consequences. This demonstrates that it will be essential for future studies to consider the effects of viscosity and magnetic fields together.

3.5 Cold fronts and thermal conductivity

Due to the small electron Larmor radius, thermal conduction in the ICM should be greatly reduced perpendicular to magnetic fields. Given that magnetic fields readily align with cold fronts, we expect thermal conduction across cold fronts to be greatly reduced. Xiang et al. (2007) used pure HD simulations of remnant-core cold front formation with isotropic thermal conduction and various suppression factors for the thermal conductivity to show that to reproduce the width of the cold front in A3667 the heat flux across the interface must be suppressed by a factor of ${\sim}$ 67. They also found that the thickness of the diffusive layer near the stagnation point of the cold front does not depend on the distance along the front from this point, implying that the front does not widen with distance from this point.

The previously described Asai-MHD simulations were the first to investigate the effects of anisotropic thermal conduction and magnetic fields on a remnant-core cold front. These simulations demonstrated that sufficient magnetic draping would occur around the developing cold front to suppress thermal conduction strongly enough that a sharp cold front edge could be maintained for time scales much longer than that required to smear out the cold front interface by uninhibited conduction, regardless of whether or not the initial magnetic field was uniform or turbulent. Figure 23 shows profiles of gas density, temperature and pressure across the cold front for models with anisotropic and isotropic thermal conduction, where in the former case the cold front is shielded from conduction by the magnetic field, but in the unmagnetized latter case conduction has completely smoothed out the front surface.

S13 also claimed to find that the magnetic field layer also suppressed conduction across the cold front interface for their nearly identical set-up, although they noted that some smoothing out of the temperature discontinuity occurred at early times before the magnetic draping layer was able to develop. Figure 24(a,b) show slices through the gas temperature of their magnetized simulations with anisotropic thermal conduction (with and without Braginskii viscosity). The interface does not appear quite as sharp as in their non-conductive simulations (see figure 19). Figure 24(c,d) highlight this fact in more detail, showing the profiles of density and temperature across the cold front surface at $z\approx -0.5$ .

Figure 23. Profiles of the gas density, temperature and pressure across the cold front interface from the simulations of Asai et al. (2007). (a) Profiles from the magnetically turbulent simulation with anisotropic thermal conduction. (b) Profiles from the unmagnetized simulation with isotropic thermal conduction. The cold front interface that is preserved by the magnetic field layer at $x\approx -1$ ( $-250~\text{kpc}$ ) in (a) are wiped out by conduction in the (b).

Figure 24. The effects of anisotropic thermal conduction and viscosity on remnant-core cold fronts, from S13. (a,b) Slices through the gas temperature in the $y$ – $z$ plane from the inviscid (a,c) and Braginskii viscosity (b,d) simulations. In these units, $1=250$  kpc. (c,d) Profiles of the gas density and temperature across the cold front interface, from the inviscid (a,c) and Braginskii viscosity (b,d) simulations. The density and temperature jump across the cold front interface at $z\approx -0.5$ ( $-125~\text{kpc}$ ) appears to be smeared out by conduction, regardless of the presence or absence of Braginskii viscosity (compare with the sharp jumps in figure 19 c,d).

ZuHone et al. (2013, hereafter Z13) were the first to carry out simulations of the formation of sloshing cold fronts in the ICM with magnetic fields and anisotropic thermal conduction. They used the same set of cluster initial conditions as ZML11, performing simulations with different thermal conductivities and magnetic field geometries, with the initial magnetic field strength set by enforcing that ${\it\beta}=400$ everywhere. The two most important cases were a simulation with full Spitzer conduction along the field lines, and an otherwise identical simulation with $1/10\text{th}$ Spitzer conduction along the field lines, modelling a situation where thermal conduction is suppressed along the field line, by strong curvature of the field lines, stochastic magnetic mirrors, or microscale plasma effects (Chandran et al. 1999; Malyshkin & Kulsrud 2001; Narayan & Medvedev 2001; Schekochihin et al. 2008, also see § 4.1.1 and references therein).

Figure 25. The effect of anisotropic thermal conduction on sloshing cold fronts (reproduced from Z13). (a–c) Slices of gas temperature through the centre of the domain in the $x$ – $y$ plane from the simulations of Z13 at the epoch $t=3.25$  Gyr. Each panel is 500 kpc on a side. Black lines indicate the locations of profiles taken across the cold fronts. (d–f) Profiles of the gas temperature taken across the lines shown in (a–c).

These simulations revealed an unexpected result: despite the presence of magnetic field layers, thermal conductivity was able to smear out the front surfaces to a certain degree. Figure 25 shows slices of temperature through the midplane of the sloshing cold fronts (a–c), along with profiles taken across the cold fronts in particular places (shown in (d–f), marked by black lines in (a–c)). In the case of full anisotropic Spitzer conduction (c), the temperature jumps have been smeared out and the interior of the cold fronts have considerably higher temperature than the case of no conductivity. In the simulation with 1/10th Spitzer conduction (b), the temperature jumps are essentially identical to those of the case without conduction. Figure 26 shows that thermal conductivity would have a significant effect on the appearance of cold fronts as seen in X-rays – full Spitzer conduction, even when restricted to the field lines, smears out the interface to such a degree that the cold fronts are essentially eliminated, as shown in figure 26. For 1/10th Spitzer conduction, the cold fronts appear roughly the same as in the case without conduction.

Figure 26. Projected X-ray emission along the $z$ -axis of the simulation domain for the simulations from Z13 with a Chandra X-ray image of A2319 included for comparison. White bars indicate 100 kpc distances. Conduction smooths out the surface brightness jumps, making them barely discernible in comparison to the sharp jumps in emission seen in A2319.

Figure 27. The effect of anisotropic thermal conduction on a simulated cluster similar to Virgo, from ZuHone et al. (2015). (a–c) Projected spectroscopic temperature for three different models. Coloured regions show the locations of profiles in (d,e). (d,e) Profiles of surface brightness and projected temperature across the cold front surfaces at the locations given by the regions in (a–c).

These results are in conflict with the initial expectation that a magnetic layer completely protects the cold front from heat conduction. Not surprisingly, the devil is in the details – the full simulations reveal two reasons why the simple picture, drawn from analytic considerations, breaks down. First, although the magnetic field is stretched tangential to the front surface by the shear flow, KHI may develop rapidly along parts of the front, tangling the field lines at the interface once again, allowing some heat to diffuse between the hot and cold phases (cf. figure 8 of Z13). This effect is also seen in the plane-parallel simulations of Lecoanet, Parrish & Quataert (2012), which showed that perturbations in the magnetic field direction across an interface similar to a cold front can cause a similar smearing of the interface. Only in a scenario where the amplified magnetic field were strong enough to prevent KHI completely would this effect be absent. Secondly, the 3-D nature of the problem and the geometry of the cold fronts leave the cold gas underneath the front magnetically connected to hot gas above the front and perpendicular to the sloshing plane, providing yet another avenue for the transfer of thermal energy (cf. figures 9 and 10 of Z13). We will discuss the discrepancies between the different simulations of cold fronts with anisotropic thermal conduction in § 4.1.3.

In the previously described work by Z15, simulations with anisotropic thermal conduction were also included, with and without Braginskii viscosity. Regardless of whether or not the ICM is viscous, the results from these simulations are the same as in Z13: anisotropic thermal conduction smears out cold fronts to such a degree that they would be unobservable, even by Chandra (see figure 27). This implies that these constraints on thermal conduction may apply generally to clusters, as the cluster model resembling the Virgo cluster from ZuHone et al. (2015) has a much colder average temperature than the more massive cluster model from Z13, and hence much weaker thermal conductivity, due to the strong dependence of the latter on the temperature ( $\propto T^{5/2}$ ).

4 Discussion

4.1 Open questions

4.1.1 The role of plasma instabilities

The Braginskii MHD equations outlined in § 3.1 do not have an explicit dependence on the Larmor radius; it is simply assumed to be much smaller than the mean free path. However, a number of recent works have shown that the effects of a finite Larmor radius may be relevant on larger scales in the ICM. In particular, when the pressure anisotropy violates the approximate inequalities

(4.1) $$\begin{eqnarray}\displaystyle -\frac{B^{2}}{4{\rm\pi}}\lesssim p_{\bot }-p_{\Vert }\lesssim \frac{B^{2}}{8{\rm\pi}}, & & \displaystyle\end{eqnarray}$$

a situation which is expected to occur readily in the cluster plasma (Schekochihin et al. 2005; Lyutikov 2007; Kunz et al. 2011), rapidly growing Larmor-scale instabilities (namely, the firehose instability on the left-hand side of the equation and the mirror instability on the right-hand side) are triggered and act to regulate the pressure anisotropy back to within its stability boundaries. The effect of these microscale instabilities on the pressure anisotropy and, conversely, the transport of momentum, is an open question.

Collisionless kinetic simulations of the driven firehose instability (e.g. Matteini et al. 2006; Hellinger & Trávníček 2008; Kunz, Schekochihin & Stone 2014) have shown that the collisionality of the plasma, supplemented by the anomalous scattering of particles off the microscale fluctuations, adjusts to maintain a marginally firehose-stable pressure anisotropy. This effectively reduces the parallel viscosity of the plasma, making it possible for KHI to develop more easily at surfaces like cold fronts. In addition, with the pressure anisotropy microphysically pinned at the firehose stability threshold (the left inequality of (4.1)), the resulting viscous stress would effectively cancel the magnetic tension, again easing development of KHI. In the case of the mirror instability (see Kunz et al. 2014; Riquelme, Quataert & Verscharen 2015), the pressure anisotropy appears to be regulated not just by anomalous particle scattering but also by an increasing population of resonant particles becoming trapped in magnetic mirrors where the pressure is naturally less anisotropic.

Consideration of these plasma effects requires choices to be made when employing the Braginskii MHD equations in simulations of cold fronts. In such simulations, when the pressure anisotropy exceeds the bounds of (4.1), these instabilities will be triggered, but since the Larmor radius is severely underresolved, they will occur on the smallest scale available, that of the cell size. The effects of the instabilities will then be damped out by the numerical diffusivity associated with the grid scale. Alternatively, as suggested by the kinetic simulations, the pressure anisotropy may be limited within the bounds suggested by (4.1), which amounts to an increase of the parallel Reynolds number of the ICM.

In regions near cold fronts, these effects may not be of significant concern, as the magnetic field is amplified to near-equipartition scales, and the boundaries in (4.1) will be more difficult to transgress. However, future simulation studies of cold fronts using the Braginskii MHD equations should investigate the potentially observable effects on cold fronts of limiting the pressure anisotropy or not.

4.1.2 Distinguishing between the effects of magnetic fields and viscosity

All of the cold front simulations we have detailed in this review have shown that magnetic fields and viscosity are both capable of suppressing KHI at cold front surfaces, provided that their effects are strong enough. The question then remains as to which effect is most responsible for this suppression in the cold fronts that lack evidence of KHI. Does one effect dominate over another, or are both effects important?

The major difficulty in answering this question are the lack of strong independent constraints on either effect apart from the lack of KHI at cold fronts itself. From the observational side, we do not currently have an independent check on the strength of the magnetic field in the layer at the cold front surface. In principle, an estimate of the field strength along the line of sight could be made using RM observations in the vicinity of the front. However, the strength of this component may or may not be comparable to the orthogonal component in the plane of the sky where the cold fronts are being observed, especially if the latter component is increased by shear amplification and/or draping effects at the front surface.

The most promising way to constrain the magnetic field strength in the layer would be to observe a deficit in the thermal pressure at a front surface, providing an essentially direct measurement of the magnetic pressure at that location, since the sum of the thermal and magnetic pressures should be continuous, provided the measurement is made on the side of the front where the magnetized layer is located. This would be easiest to check in sloshing cold fronts, which are isobaric across the front surface. Recently, Reiss & Keshet (2014) used deprojected pressure profiles across a number of sloshing cold fronts from the literature and found that many of the fronts in their sample have pressure discontinuities consistent with the existence of enhanced magnetic pressure underneath the front surface. The average ratio of thermal to total pressure in these regions was ${\it\xi}\simeq 0.8$ , which implies ${\sim}$ 20 % of the total pressure is in the form of magnetic fields. The average value is statistically significant, although for many of the cold fronts in their sample the significance is modest (up to a few ${\it\sigma}$ ). They also did not examine the characteristics of each individual front in detail to determine if the presumed magnetization in that front would be sufficient to prevent the development of KHI. However, these initial results are promising, and deserve future study and longer exposures by Chandra.

For the viscosity, if we assume its origin is tied to the ion collisions in the ICM, it will have the Spitzer form (3.9), which can be easily determined from the plasma properties. However, the effective viscosity within a small region will be dependent on the magnetic field direction, and as noted in the previous section there are still large theoretical uncertainties regarding the effect of microscale plasma instabilities on the collisionality of the plasma.

To break this apparent degeneracy between magnetic fields and viscosity, simulations should investigate other aspects of the cold fronts and their near vicinity and make predictions for what observable effects may arise from either aside from the suppression of KHI. For example, R13 showed in their unmagnetized simulations of the sloshing cold fronts in Virgo that in the absence of viscosity, the cold gas underneath the front surface should exhibit multiple edges in surface brightness as seen in projection. These edges are absent if the plasma is viscous. It is not yet known whether the same features would be present in a cold front with a magnetized layer. In the context of remnant-core cold fronts, the characteristics of the edges of the cold front and the stripped tail and wake of the subcluster explored in detail by R15A and R15B may be compared to otherwise identical simulations with magnetic fields to determine which effect best matches the observations.

4.1.3 How efficient is thermal conduction in the ICM?

As detailed in § 3.5, MHD simulations of cold fronts with anisotropic thermal conduction have given mixed results as to whether or not magnetic fields can effectively protect the sharp density and temperature jumps from being smeared out by an efficient thermal conductivity. In the simulations reviewed here, magnetic fields seem to largely prevent the sharp temperature and density gradients of remnant-core cold fronts from becoming smeared out due to conduction, but cannot do so in the case of sloshing cold fronts (though we again note as we did in § 3.5 that the cold front jumps do appear to be somewhat smeared out by conduction in the S13 simulations; see figures 19 and 24). Why is this the case?

Several factors may be at play. As noted above, the behaviour of magnetic fields at the cold front surface is rather different in the case of remnant-core cold fronts and sloshing cold fronts: in the former, the field is ‘draped’ on the outside of the front, in the latter, the field is amplified by shear amplification on the inside of the front. The draping layers on the outside of the fronts tend to be more long lasting in the simulations than the layers on the inside of the fronts, which are more transient (see § 3.3). The geometry of the cold front is also an important consideration. In the Asai-MHD and S13 simulations, the magnetic field layer appears to wrap essentially all the way around the cold front, strongly restricting conductive heat flux across the interface. As noted above, in the sloshing simulations of Z13 and Z15, tangling of field lines along the front surface by KHI and the presence of other regions of hot gas which are connected to the cold gas by magnetic fields permit a flux of heat to the cold phase. However, the KHI in these simulations were able to develop due to the initial magnetic field strength being set by ${\it\beta}=400/1000$ (Z13/Z15); in a more strongly initially magnetized simulation ( ${\it\beta}\lesssim 100$ ), the draping layers may have been more effective.

It should also be noted that the remnant-core cold front simulations from the Asai-MHD and S13 works have poorer resolution ( ${\sim}$ 7–10 kpc, in the 3-D simulations) than the sloshing cold front simulations by ZuHone ( ${\sim}$ 1–2 kpc). The lower resolution simulations have higher numerical dissipation, preventing some KHI modes from growing and limiting the development of turbulence, both of which could retangle magnetic field lines if they were not strong enough to prevent the growth of either effect in the first place. For these reasons, the remnant-core scenario with anisotropic thermal conduction should be re-examined at higher resolution.

The question of whether or not the narrow widths of cold front interfaces are even compatible with anisotropic thermal conduction is important for the study of the thermodynamic properties of the ICM as a whole. A considerable amount of attention in recent years has been devoted to the effect of anisotropic thermal conduction on the stability properties of the ICM. In the classical picture, the ICM is convectively stable since the entropy in clusters is stratified with $\text{d}S/\text{d}r>0$ . However, a number of studies (Balbus 2000, 2001; Parrish & Quataert 2008; Parrish, Stone & Lemaster 2008; Quataert 2008; Bogdanović et al. 2009; Parrish, Quataert & Sharma 2009) have demonstrated that with efficient thermal conductivity in the ICM, the entire cluster atmosphere is unstable to convection provided that $\text{d}T/\text{d}r\neq 0$ , which is observed to be true for nearly all clusters, especially those with cool cores (though other works have emphasized that these instabilities are affected strongly by the inclusion of Braginskii viscosity; see Kunz et al. 2011, 2012). Other works have emphasized the possible role of thermal conduction in setting the temperature profiles of clusters (Zakamska & Narayan 2003; Guo, Oh & Ruszkowski 2008; McCourt, Quataert & Parrish 2013), partially offsetting cooling in cores (Bertschinger & Meiksin 1986; Bregman & David 1988; Ruszkowski & Begelman 2002; Kim & Narayan 2003; Voigt & Fabian 2004; Parrish, Quataert & Sharma 2010; Ruszkowski & Oh 2010, 2011) and influencing the condensation of multiphase gas in cluster cores and giant elliptical galaxies (McCourt et al. 2012; Sharma et al. 2012; Wagh, Sharma & McCourt 2014; Voit & Donahue 2015; Voit et al. 2015; Yang & Reynolds 2016). Therefore, determining if cold fronts demonstrate that the thermal conductivity in the ICM is very low would have a significant impact on our understanding of the thermodynamics of the cluster core. Such a suppression of thermal conductivity has been suggested by previous theoretical studies, either due to strong curvature of magnetic field lines at small scales or by microscale plasma instabilities (Chandran et al. 1999; Malyshkin & Kulsrud 2001; Narayan & Medvedev 2001; Schekochihin et al. 2008).

4.2 Future observations

More observations are needed to increase the number of cold front observations with high-quality, long-exposure images and spectra for simulations to compare against. To best constrain the plasma properties, observations of nearby clusters with Chandra are most useful. For instance, other portions of the largest Virgo cold front than the one observed by Werner et al. (2016) should also be observed with similar exposure ( ${\sim}500$  ks). Simulations predict (R13) that the velocity shear may change along the length of the front, indicating that the magnetic field strength in the layer and/or the presence of KHI may be different at other locations along the front.

An important leap forward in the use of cold fronts for studying the physics of the cluster plasma will be provided by the next generation of X-ray observatories, which will have effective areas that far surpass what is currently available. There is no doubt that the increased effective area of Athena (Nandra et al. 2013), with over an order of magnitude improvement in area over that of current instruments, will greatly improve the statistics of measuring the temperature on either side of cold fronts, but this improvement will be blunted by the poorer spatial resolution of Athena ( ${\sim}5^{\prime \prime }$ ) with respect to Chandra’s superb subarcsecond resolution. The mission concept X-ray Surveyor (Weisskopf et al. 2015) would provide both a similar increase in effective area, coupled with Chandra-like spatial resolution.

The next generation of X-ray observatories will also possess high resolution microcalorimeters which will be able to measure the properties of the velocity field of the ICM from the shifting and broadening of emission lines in the X-ray spectrum. The first such mission will be Astro-H (Takahashi et al. 2014), to be launched in early 2016. A recent study carried out by ZuHone et al. (2016) used synthetic observations to show that Astro-H will be able to measure the basic kinematic properties of sloshing cold fronts, namely the shift of the spectral lines produced by the bulk motion of the sloshing gas (provided that our line of sight is at least partially within the sloshing plane) and the line broadening produced by the variation in this bulk motion and sloshing-driven turbulence. Unfortunately, they also showed that due to Astro-H’s poor spatial resolution of ${\sim}$ 1’, it will be unable to place any meaningful constraints on the ICM viscosity from studies of cold fronts. Meaningful constraints on the cluster plasma from kinematics will require a mission with both high spatial and spectral resolution, such as Athena, or (more likely) X-ray Surveyor, both of which will have calorimeters with spectral resolution equivalent to or higher than Astro-H.

5 Summary

Cold fronts are important probes of plasma physics in galaxy clusters. We presented a review of recent simulation works, along with examples of relevant observations, which have provided the strongest constraints on the properties of the fronts themselves as well as the ICM physics implied by these features. Many insights have been gained from these simulation studies, but a number of open questions remain. To summarize:

1. (i) Recent simulations conclusively demonstrate that the surfaces of cold fronts should be locations of strong magnetic field layers oriented parallel to the front surface, whether on the outside of the front (as in remnant-core cold fronts) or the inside of the front (as in sloshing cold fronts). Such layers may suppress KHI, explaining the fact that KHI do not appear to be present in many observations of cold fronts. However, whether or not this is the case depends strongly on the initial magnetic field strength of the surrounding medium.

2. (ii) In the absence of strong magnetic fields, viscosity may play an important role in explaining the lack of evidence of KHI in many cold fronts. Recent simulations have demonstrated that even a small viscosity may be enough to explain the observed features of some cold fronts. On the other hand, a number of cold fronts do show evidence of KHI, in the form of distortions such as the ‘horns’ seen in the cold front of the remnant-core cold front in M89, and the ‘ragged’ appearance of some cold fronts as seen in NGC 7618, UGC 12491 and A496. In these systems, such features imply a low viscosity for the cluster plasma. More detailed comparisons between observations and simulations for a variety of systems are necessary to determine what range of Reynolds number is permissible in the ICM.

3. (iii) Where KHI suppression is observed in cold fronts, it is not clear which effect is likely to be more responsible – magnetic fields or viscosity? The most straightforward way to break this degeneracy is to use long-exposure Chandra observations of cold fronts to attempt to measure the thermal pressure deficit at the interface arising from the presence of a strongly magnetized layer. More quantitative comparisons of the observable signatures of cold fronts predicted from simulations with magnetic fields and viscosity are also needed.

4. (iv) Whether or not magnetic fields suppress thermal conduction across cold fronts appears to depend on the type of cold front. Simulations of remnant-core cold fronts suggest that the highly magnetized layers effectively prevent the flow of heat to the colder side of the front. In contrast, for sloshing cold fronts, simulations predict that the regions above and below the cold fronts are connected with magnetic fields along which sufficient heat conduction can occur, erasing the sharp cold fronts. Although there are geometrical differences between the two types of cold fronts which may be responsible for this difference, the remnant-core simulations of this effect performed so far have low resolution. More MHD simulations of remnant-core cold fronts with anisotropic thermal conduction are needed to determine if these differences between these two types of cold fronts persist at higher resolution or is numeric in nature.

5. (v) The focus of future simulations should be to employ the most detailed and computationally feasible physical models available. This includes the full modelling of the Braginskii MHD equations, with possible Larmor-radius-scale corrections as provided by future kinetic simulations of the cluster plasma. From these, more accurate hydrodynamic approximations may be derived. The role of magnetic reconnection near cold front surfaces should also be investigated.

6. (vi) In future comparisons between simulations and observations of cold fronts, careful attention needs to be paid to predicting the distinct observable signatures of cold fronts under different assumptions for the underlying physics, since only if we can recognize and characterize these signatures accurately can they be used as a probe of the ICM physics. The most accurate comparisons to the observations are facilitated by synthetic X-ray observations of cold fronts, which include the effects of Poisson noise and instrumental responses.

7. (vii) A more accurate determination of the plasma properties of the ICM from cold fronts will be provided by the next generation of X-ray telescopes with high effective area and high spectral resolution, such as Athena and X-ray Surveyor, although the latter will be essential for examining the cold front interfaces with the same spatial resolution as Chandra.