A VERY DEEP CHANDRA OBSERVATION OF THE GALAXY GROUP NGC 5813: AGN SHOCKS, FEEDBACK, AND OUTBURST HISTORY

The Chandra observations that were used in the analysis we present here are summarized in Table 1. The aimpoint was on the back-side-illuminated ACIS-S3 CCD for each observation. All data were reprocessed from the level 1 event files using CIAO and CALDB 4.4.3. CTI and time-dependent gain corrections were applied. lc_clean was used to remove background flares.⁹ The mean event rate was calculated from a source-free region using time bins within 3σ of the overall mean, and bins outside a factor of 1.2 of this mean were discarded. There were no periods of strong background flares. The cleaned exposure times are given in Table 1, for a total time of 635 ks.

Diffuse emission from N5813 fills the image field of view (FOV) for each observation. We therefore used the CALDB ¹⁰ blank sky background files appropriate for each observation, normalized to match the 10–12 keV count rate in our observations to account for variations in the particle background. To generate exposure maps, we used a MEKAL model with $kT=0.7$ keV, Galactic absorption, and abundance of 30% solar at a redshift $z=0.006578$.

To study the thermal structure of the ICM, we generated a smoothed temperature map using the method described in Randall et al. (2008) and also used in R11. For each pixel in the temperature map, spectra were extracted from a surrounding circular region that contained 1500 source counts. The temperature map pixel value was then set equal to the best-fitting temperature from fitting to these spectra. Since the extraction regions are generally larger than the pixels, nearby pixels are correlated with one another, such that the temperature map is effectively smoothed on the scale of the local extraction region size. Some advantages of this method are that there is no a priori assumption about the thermal structure of the gas, and experience has shown that such maps can resolve real structure on scales that are smaller than the size of the extraction regions.

For comparison, we also constructed temperature maps using the contour-binning method developed by Sanders (2006). Here, extraction regions are defined based on surface brightness contours. Regions were defined to achieve a signal-to-noise ratio of 38 (corresponding to ∼1400 counts in bright regions) in each region. This method assumes that the temperature distribution follows the surface brightness distribution (which is frequently, but not always, the case). Some advantages of this method are that each extraction region is independent of the others (as long as the region widths are large compared with the local point-spread function [PSF]), and it is relatively cheap computationally (as compared with the smoothed temperature map method). A comparison of the smoothed and contour binned temperature maps gives an indication of the extraction region size as a function of position in the smoothed map.

Smoothed temperature maps over a wider area and a higher-resolution map of the core are shown beside the contour-binned maps in Figures 5 and 6, respectively. Both the smoothed and contour-binned temperature maps reveal relatively hot emission associated with the bright rims surrounding the inner cavities and just inside the prominent 10 kpc surface brightness edges. These features were identified as shocks in R11. The maps suggest additional temperature rises associated with the outer 30 kpc surface brightness edges, particularly to the NW. The temperature profiles across these edges are presented in Section 4.2.2. The temperature maps also show a plume of cooler gas, roughly extending along the line defined by the cavities, from SW to NE. In R11, we showed that this plume is consistent with arising from cool gas that has been lifted by the buoyantly rising X-ray cavities.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We derived pseudo-pressure and pseudo-entropy maps from the smoothed temperature maps. The fitted apec normalization for each spectral map pixel is proportional to the volume integral of the square of the electron density along the line of sight. The square root of the normalization therefore gives an average "projected" density along the line of sight. We define the pseudo-pressure and pseudo-entropy as $kT{{A}^{1/2}}$ and $kT{{A}^{-1/3}}$, respectively, where kT is the fitted projected temperature and A is the normalization scaled by the area of the extraction region. The wide field and core maps are shown in Figures 7 and 8. The pseudo-pressure maps reveal prominent pressure jumps across the 10 kpc shock fronts and increased pressure in the bright central rims (i.e., the 1 kpc shock). There is no obvious equivalent structure in the pseudo-entropy maps. This is not unexpected for weak shocks. For the strongest shock we detect, with $M=1.78$, the Rankine–Hugoniot shock jump conditions for a $\gamma =5/3$ gas give a pressure jump at the shock front of a factor of 3.7, but a jump in the entropy index of only a factor of 1.12. Thus, entropy jumps are expected to be relatively difficult to detect, particularly in projection.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The core pseudo-entropy map (Figure 8) shows some filamentary structure in the center, on the scale of the bright shocked rims. However, unlike the features in the pseudo-pressure map, these features do not directly trace the bright rims. Furthermore, the entropy is lower in these filaments as compared with the surrounding regions, in contrast to the small increase in entropy that is expected at the shock fronts. These structures are likely formed as the central low-entropy gas is pushed and pulled away from the center as the inner cavities expand and rise buoyantly.

The central regions of the core temperature and pseudo-entropy maps are shown with the Hα contours from R11 overlaid in Figure 9. As shown in R11, the Hα emission follows the plume of cool X-ray-emitting gas up to the location of the intermediate cavities, indicating that these cavities have uplifted cold-phase gas from the core as they buoyantly rise. The Hα emission also clearly traces the SW inner cavity, suggesting that the central cold-phase gas is dynamically coupled to the cavity and is being pushed out as the cavity expands. The SE 1 kpc shock appears to be passing through the Hα filament (at least in projection) without destroying it. Finally, there is a correlation between the Hα emission and the detailed, small-scale, filamentary structure of the central low-entropy gas. This is consistent with the above suggestion that these entropy features likely arise from gas that has been displaced from the core by the cavities, rather than from entropy changes driven by the central shocks.

Zoom In Zoom Out Reset image size

Since each extraction region in the above maps contained only 1500 net counts, the abundance values determined by these fits were not tightly constrained, and the corresponding abundance maps show little to no structure. Abundance measurements in ∼0.7 keV gas are challenging due to the large number of closely spaced emission lines and the resulting degeneracy between line and continuum emission at Chandra's energy resolution. We found that, for N5813, achieving 1σ errors on the order of 5%–10% for the abundance required roughly 30,000 net counts. Fortunately, the new deeper observations provided well over one million photons appropriate for use in spectral fitting, allowing us to map the abundance distribution in the gas. A smoothed abundance map, with 30,000 net counts per extraction region, is shown in Figure 10.

Zoom In Zoom Out Reset image size

The abundance map shows an apparent plume of low-abundance gas that corresponds to the plume of uplifted cool gas seen in Figure 5 (roughly 40 kpc long, extending from SW to NE). It has long been recognized that fitting multi-temperature spectra with a single thermal model can give anomalously low abundance values (the so-called "Fe-bias" effect; Buote 1999). To test for this effect, we fit multi-temperature models to the apparent central minimum and apparent low-metallicity NE plume. We find that the plume is well described by a two-temperature model, with $k{{T}_{{\rm low}}}\approx 0.40$ keV, $k{{T}_{{\rm high}}}\approx 0.75$ keV, and an abundance (set equal for each temperature component) of ${{Z}_{{\rm plume}}}\approx 45$% solar, consistent with the surrounding ICM. In contrast, we find that, while including multiple temperature components for the central region (in a $2.2\times 3.8$ kpc ellipse oriented NE to SW, along the plume) does raise the fitted abundance, it is still significantly lower than the abundance of the surrounding ICM, at ${{Z}_{{\rm center}}}\approx 37$% solar. Thus, we conclude that the apparent low abundance in the extended plume is a fitting artifact arising from the Fe-bias effect, and that the measured abundance in the plume is consistent with that of the surrounding ICM, while the apparent abundance minimum within the central few kpc cannot be easily explained by this effect. We note that central abundance dips that apparently cannot be explained solely by projection effects have been seen in other systems (e.g., Blanton et al. 2003; Panagoulia et al. 2013). In addition to the apparent low-abundance plume, there appear to be increases in the projected abundance at the location of the shock fronts, in contradiction with what is expected from the Fe-bias effect. A more detailed study of these features and the ICM abundance in general in N5813 will be presented in an upcoming paper.

4.2.1. Azimuthally Averaged Profiles

We produced projected radial profiles by fitting spectra extracted from concentric annuli, centered on the centroid of the diffuse emission at large radii. Each annulus was chosen to be at least 1'' wide and contain at least 1000 net counts in the 0.6–3.0 keV band. The inner annuli were limited by the size constraint (and therefore contained more than 1000 net counts), whereas the outer annuli were wider and contained roughly 1000 net counts. The resulting temperature profile is shown in the top panel of Figure 11. The positions of the shock fronts (determined in Section 5) are marked by vertical lines. Even in the projected temperature profile in circular annuli (which effectively smooths out the elliptical shock fronts), temperature enhancements associated with each shock are clearly visible. The central 15 radius region, which contains the central AGN, has been excluded from the fit to the central region. The central temperature rise is due to the presence of the inner shock rims surrounding the inner cavities.

Zoom In Zoom Out Reset image size

As in R11, we determined the 3D structure of the ICM by performing an "onion peeling" deprojection analysis (Fabian et al. 1981) using concentric annuli, each of which contained at least 8000 net counts (in projection) in the 0.6–3.0 keV band. We fit each annulus with an absorbed apec model, with the abundance fixed at 50% solar (consistent with results from spectral fits to the total diffuse emission and with typical values in the projected abundance map shown in Figure 10). We assume spherical shells of uniform density. Given the rich morphology of N5813, with elliptical shock edges, cavities, and other features, it is clear that the assumption of spherical symmetry does not strictly apply. Thus, the deprojected profiles include some level of systematic uncertainty associated with the assumed geometry. The deprojected temperature, density, pressure, and entropy profiles are shown in Figure 11. We see that the deprojected temperature is generally consistent with the projected temperature. We interpret this an indicating that the scale of the azimuthally averaged temperature gradient is generally small compared with the sizes of our extracted annuli, so that projection effects in any given annulus (which are dominated by emission from adjacent annuli) are relatively small. The small discontinuous jumps, or "kinks," in the density, pressure, and entropy profiles around 11 kpc are associated with the sharp surface brightness edges to the NW and SE from the intermediate shock.

4.2.2. Shock Front Profiles

To confirm the temperature rises across the surface brightness edges seen in Figures 5 and 6, we derived temperature profiles in annular bin sectors across each of the edges. The profiles were centered on the center of curvature for each feature (such that there was a bin boundary that traced the surface brightness edge), and then the distances were corrected to measure the distance from the center adopted for the azimuthal profiles presented in Section 4.2.1.

The temperature profiles across the 10 and 30 kpc shock fronts are shown in Figures 12 and 13, respectively, with the fitted shock front positions overlaid (see Section 5). Clear temperature rises are seen at each front, confirming that these features are shocks. In the case of the 10 kpc shocks, the deeper observations allow us to measure the temperature profiles with greater angular resolution (on the scale of a few hundred pc) as compared with R11. In addition to the overall shift in normalization due to the updated AtomDB (see Section 4), the peak change in temperature across the fronts is larger as compared with the more coarsely binned profiles in R11. In contrast, the shock Mach numbers that we derive here (Section 5) are consistent with those derived in R11, since the surface brightness profile requires fewer counts per radial bin, allowing for smaller bins. This is consistent with our suggestion in R11 that the temperature jumps we found there were smaller than what was expected based on the Mach numbers derived from the surface brightness jumps due to projection effects.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In some cases, there appears to be shock-heated gas just outside the positions of the shock front (e.g., the two bins to the right of the edge in Figure 12, left). There are at least two factors that can contribute to this effect. First, there is a small uncertainty in the fitted locations of the shock fronts. Second, we show in Section 5 that the shock fronts have finite widths (on the order of 0.5 kpc for the 10 kpc shock). When fitted with a discontinuous model, the preferred position of the edge will naturally be at the center of the smoothed front.

Shocks are expected to heat the ICM as they propagate. While there is a relatively large temperature jump at the shock front, this increase is transient as the gas subsequently expands and cools. The relevant quantity for lasting heating of the gas is the change in entropy. The equivalent amount of heat energy ${\Delta }Q$ imparted to the gas by a shock due to a change in entropy ${\Delta }S$ is given by

$$\begin{eqnarray}&&{\Delta }Q\simeq T{\Delta }S=E{\Delta }{\rm ln} \frac{p}{{{\rho }^{\gamma }}},\end{eqnarray} \tag{ 1 }$$

where $E={{C}_{V}}T$ is the total thermal energy of the gas, γ is the adiabatic index (taken to be 5/3), and p and ρ are the pressure and density, respectively. Thus, each shock contributes a fraction ${\Delta }{\rm ln} (p/{{\rho }^{\gamma }})$ of the total thermal energy in the gas. In the case of weak shocks, this fraction is small (between 0.4% and 12% for the shocks we detect in N5813). To compare with radiative cooling, we are interested in the cumulative effect of shock heating per local cooling time in the gas (i.e., the number of shocks per local cooling time and the total amount of heat energy deposited in the gas). We approximate the local cooling time just outside each shock edge as the time it would take the gas to radiate away all of its thermal energy based on its deprojected emissivity. Based on the shock ages given in Table 2, and keeping in mind that these ages are expected to be underestimates for young shocks and overestimates for older shocks, we assume an outburst repetition rate of one every 20 Myr.

The results from the above comparison of the shock heating and radiative cooling rates are summarized in Table 3. Given the various statistical and systematic uncertainties associated with, e.g., the shock Mach numbers, the shock ages and outburst repetition rate, the deprojected cooling luminosity, etc., these numbers should be considered to be rough estimates. Nevertheless, there is good agreement between the number of shocks expected per cooling time (column 3) and the number required to completely offset radiative cooling (column 4). We conclude that shock heating alone is sufficient to balance radiative cooling in the gas out to at least ∼30 kpc. This extends our earlier result in R11, where we found a balance between heating and cooling out to ∼10 kpc. At 30 kpc, the cooling time (∼2 Gyr) is approaching the expected age of the system, indicating that there is no need for additional heating beyond this radius to completely balance cooling and explain the lack of a strong cooling flow in this system.

Table 3.  Shock Heating and Radiative Cooling

ID	${{t}_{{\rm cool}}}$ a	Shocks/${{t}_{{\rm cool}}}$ b	${\Delta }Q/E$ c
	(108 yr)
1 kpc Shock	1.4	7	9
10 kpc Shock	9.2	46	21
30 kpc Shock	22.2	111	143

Notes.

^aLocal cooling time of the gas, just outside of the shock edge. ^bNumber of shocks per cooling time, assuming an outburst repetition rate of one every 20 Myr. ^cReciprocal of the fraction of the total thermal energy E added by each shock, roughly equal to the number of shocks per cooling time required to offset radiative cooling.

Download table as:  ASCIITypeset image

We can estimate the average cavity power ${{P}_{{\rm cav}}}$ by dividing the total enthalpy of the cavities (estimated as $4pV$ for each cavity, for a total of $5.7\times {{10}^{57}}$ erg s⁻¹) by the age of the oldest cavity pair ($\sim {{10}^{8}}$ yr; see R11, Table 3). This gives ${{P}_{{\rm cav}}}\approx 1.8\times {{10}^{42}}$ erg s⁻¹. Spectral fitting gives a total X-ray luminosity of ${{L}_{X}}\approx 5.5\times {{10}^{41}}$ erg s⁻¹ within a radius enclosing all of the cavities (∼26 kpc). Thus, in N5813 we find that a few times the pV work required to inflate the cavities is in principle large enough to offset radiative cooling, as has generally been found in other systems (e.g., Bîrzan et al. 2004; Rafferty et al. 2006; Nulsen et al. 2007; Hlavacek-Larrondo et al. 2012).

Since we have shown that, in N5813, shock heating alone is sufficient to balance cooling, one might ask whether the total heating from cavities and shocks leads to an "overheating" of the ICM, breaking the balance between heating and cooling that is required for AGN feedback to operate over long timescales. Although cavities are expected to deposit some of their enthalpy in the ICM as they rise buoyantly (McNamara & Nulsen 2007), this heating is not directly observed, and the details on how and where cavities heat the ICM are poorly understood. Furthermore, cavities can only directly heat the ICM locally, as they rise radially, whereas roughly isotropic heating must take place close to the central AGN to regulate feedback. Mixing may help to distribute locally heated gas (e.g., Zhuravleva et al. 2014), although the gas motions must be small enough to not completely destroy the observed X-ray cavities. Shocks, in contrast, heat the gas roughly isotropically and most strongly at smaller radii, and the thermal effect of shocks on the ICM is, in some cases, directly observed. We suggest that AGN outburst shocks may generally play a significant role in AGN feedback, particularly at small radii where the Mach numbers are higher. The cavities are then free to rise buoyantly and release their energy to heat the ICM at larger radii. The fact that cavities are observed at larger radii demonstrates that they do not necessarily release their enthalpy to the ICM at small radii, consistent with this picture. Weak shocks are generally difficult to detect, since the shock fronts are thin and easily masked by projection effects. They are more easily detected in N5813 due to its proximity, low gas temperature, and regular morphology. The fact that one of the handful of examples of confirmed outburst shocks known in clusters, in Abell 2052 (Blanton et al. 2011), was confirmed by a measured temperature jump only after very deep Chandra observations and careful deprojection is consistent with this scenario. Similarly, no temperature jump is found associated with the presumed outburst shock in the Perseus cluster using extremely deep Chandra observations, even though it is the X-ray-brightest cluster (Graham et al. 2008).

Measurements of the cavities and shocks allow us to estimate the total mechanical outburst energy for each of the three outbursts in N5813. Since the mechanical luminosity dominates the total AGN energy output in N5813 (R11), as is generally the case for kinetic mode AGN feedback (Fabian 2012), this gives us information on the total outburst energy history of the AGN. We take the cavity internal energy to be $3pV$, where the pressure p is taken from the azimuthally averaged pressure profile. The assumed cavity dimensions and derived internal energies are given in Table 4.

Table 4.  Properties of the X-Ray Cavities

ID	aa	bb	rc	${{E}_{\operatorname{int}}}$ d
	(kpc)	(kpc)	(kpc)	(1056 erg)
Inner, SW	0.95	0.95	1.3	0.5
Inner, NE	1.03	0.93	1.4	0.6
Middle, SW	3.9	3.9	7.7	13.1
Middle-1e, NE	2.9	2.2	4.9	5.5
Middle-2e, NE	2.8	2.4	9.3	3.8
Outer, SW	5.2	3.0	22.2	3.9
Outer, NE	8.0	4.4	18.0	15.6

Notes.

^aSemi-major axis. ^bSemi-minor axis. ^cDistance from central AGN. ^dThe internal energy of the cavity, estimated as $3PV$. ^ePart of a "split" cavity.

Download table as:  ASCIITypeset image

We use the cavity internal energy rather than the enthalpy (which is the sum of the internal energy and the work required to inflate the cavity, and roughly equal to $4pV$ for cavities filled with relativistic particles) since the work done during cavity expansion goes into driving shocks, the energy of which we account for separately. As in R11, we estimate the outburst energy in shocks as

$$\begin{eqnarray}&&{{E}_{s}}\approx {{p}_{1}}{{V}_{s}}\left( {{p}_{2}}/{{p}_{1}}-1 \right),\end{eqnarray} \tag{ 2 }$$

where p₁ and p₂ are the pre- and post-shock pressures, respectively, and ${{V}_{s}}$ is the total volume enclosed by each roughly ellipsoidal shock surface. Although Equation (2) is expected to slightly underestimate the shock energy as the Mach numbers were larger at earlier times (giving larger pressure jumps), in R11 we showed that this estimate agrees well with shock energies derived from simple 1D hydrodynamic simulations of a central point explosion by matching the simulated surface brightness profiles with the observations. The total mechanical energy output of each outburst is simply the sum of the total cavity and shock energies. We note that the ratio of the energy in cavities to that in shocks is between 0.15 and 0.3 for each outburst, and that in the case of the central outburst, where the cavities occupy a significant fraction of the volume enclosed by the 1 kpc shock surface, this ratio is somewhat lower than the expected value of ∼1. This is likely a result of underestimating the central gas pressure, and hence the central cavity energies, due to the complicated morphology in this region, which leads to systematic uncertainties associated with projection effects and a lack of knowledge on the volume filling factor of the X-ray-emitting gas. However, we note that such an underestimate would only moderately affect our total outburst energy, well within the factor of a few at which we expect our simple estimates to be accurate, and hence would not affect our main conclusions.

Using the above estimates, we find total outburst energies of $8.5\times {{10}^{56}}$, $9.9\times {{10}^{57}}$, and $8.9\times {{10}^{57}}$ erg for the outbursts associated with the 1, 10, and 30 kpc shocks, respectively. As in R11, we find the outburst energy of the 1 kpc outburst to be more than a factor of 10 less than that of the 10 kpc outburst. Based on this, we concluded that either the 1 kpc outburst is ongoing, with the AGN actively inflating the cavities, or that the mean jet power varies significantly over timescales on the order of 10 Myr. Here, we find remarkable agreement between the total outburst energies of the 10 and 30 kpc outbursts. This suggests that the mean outburst (or jet) power, averaged over the outburst interval of a few times 10 Myr, is fairly stable, and that the total energy of the 1 kpc outburst is lower because this outburst is ongoing.

In R11, we estimated the outburst power as the total energy divided by the age of the outburst. Here, we instead estimate this power as the total energy divided by the duration of each outburst, which we take to be the shock age for the youngest outburst (Table 2) and the outburst interval (∼20 Myr, Section 6.1) for the remaining outbursts. This gives mean outburst powers of $1.6\times {{10}^{43}}$, $1.6\times {{10}^{43}}$, and $1.4\times {{10}^{43}}$ erg s⁻¹ for the 1, 10, and 30 kpc outbursts, respectively. Thus, we find that the mean outburst power is indeed roughly constant, certainly within a factor of a few, which is the level at which we expect our rough estimates to be accurate. We note that in R11 we found a lower power for the 1 kpc outburst. This difference arises mainly from using results from our numerical model in R11, which predicts a lower shock energy and larger shock age as compared with the estimates based on observations alone (see Table 2 in R11). The source of this difference is likely to be a systematic error that arises when applying the spherically symmetric model to the 1 kpc shock, which clearly has a more complicated morphology (see Figures 1 and 6).

We stress that our result on the consistency of the outburst power only applies to the mean power averaged over the outburst interval of a few times 10 Myr. AGN luminosities are known to vary by up to several orders of magnitude over very short (observable) timescales (e.g., Harris et al. 2009). However, in principle the outburst signatures in N5813 could be created with a constant jet power. In this case, the cavities expand rapidly and drive shocks just after they are formed. The expansion rate drops as the fractional energy input rate decreases, and the shocks separate from the cavities and propagate outward. As the cavities rise buoyantly, they eventually disconnect from the central jets (once the buoyant speed exceeds the expansion rate), which then begin inflating a new pair of cavities and the process repeats (similar models are discussed, e.g., in Fabian et al. 2003). Therefore, we use the term "outburst" to refer to the creation of a cavity pair and its associated shock, rather than a rapid increase in jet power, since these features are consistent with a kinetic jet power that either also cycles in outbursts or is roughly constant.

We conclude that the most likely scenario is that N5813 is in a "steady state" kinetic feedback mode, with the AGN outbursts roughly equally spaced in time, and each outburst depositing a similar amount of total energy into the ICM, such that the mean AGN outburst power is roughly constant. The lower total power of the 1 kpc outburst is an indication that this outburst is ongoing and has yet to deposit the bulk of its energy in the ICM. The detection of 1.4 GHz radio emission filling the inner cavities (from young, high-energy, non-thermal particles injected by the AGN), along with shock-heated cavity rims that have not yet separated from the cavities (as they have clearly done for the intermediate and outer cavities), is consistent with this interpretation (R11). This is also consistent with results from Allen et al. (2006), who find evidence for accretion flows around central AGNs that are stable over a few million years. The fact that N5813 shows little evidence of a recent merger or other complicated ICM "weather" probably contributes to it being able to maintain such steady-state feedback over long timescales.

Fits to the surface brightness profiles across the shock fronts reveal that the SE and NW 10 kpc fronts have non-zero widths that cannot be explained by Chandra's PSF (Section 5). Measurements of the widths of surface brightness edges in cluster X-ray observations have been used to place constraints on transport processes in the ICM by comparing these widths with the local collisional mean free paths of the particles. Edges that are significantly narrower than the particle mean free path indicate that Coulomb diffusion is suppressed across the edge. In the case of cold fronts in clusters, the front widths are found to be significantly smaller than the particle mean free paths (e.g., Vikhlinin et al. 2001a; Russell et al. 2012), implying that Coulomb diffusion is suppressed across the fronts. However, in cold fronts the suppression of transport processes (and Kelvin–Helmholtz instabilities) is likely due to magnetic draping, where magnetic field lines are stretched along the cold front edge as the relatively cool, dense gas moves through the ICM (Vikhlinin et al. 2001b; Asai et al. 2004, 2005, 2007; Lyutikov 2006). In this case we consider shock fronts that propagate through the ICM, and therefore magnetic draping is not expected to occur.

In the case of merger shocks, Markevitch & Vikhlinin (2007) find a bow shock width of roughly 35 kpc, on the order of the local particle mean free path, for 1E 0657-56 (the Bullet cluster), although this width is only marginally preferred over a zero-width front. Russell et al. (2012) examine shock widths in the merging cluster Abell 2146, which contains a leading bow shock and a trailing reverse shock. They find finite width shocks in each case, with the reverse shock width being significantly smaller, and the bow shock width being marginally smaller, than the particle mean free path. They conclude that transport processes are suppressed across these merger shocks. In the case of AGN outburst shocks, Croston et al. (2009) consider the thickness of the northeastern shock in Centaurus A, but their results are inconclusive.

It is therefore of interest to compare our measured shock widths with the Coulomb mean free path of particles in the ICM. In the region of the shock fronts there are four relevant mean free paths: that in the pre-shock region, ${{\lambda }_{{\rm in}}}$, that in the post-shock region, ${{\lambda }_{{\rm out}}}$, that of particles crossing the front from the post-shock to the pre-shock region, ${{\lambda }_{{\rm in}\to {\rm out}}}$, and that of particles crossing the front from the pre-shock region to the post-shock region, ${{\lambda }_{{\rm out}\to {\rm in}}}$. In our case, the largest and most relevant mean free path is ${{\lambda }_{{\rm in}\to {\rm out}}}$, which is given by

$$\begin{eqnarray}&&{{\lambda }_{{\rm in}\to {\rm out}}}=15{{\left( \frac{{{T}_{{\rm out}}}}{7\;{\rm keV}} \right)}^{2}}{{\left( \frac{{{n}_{e,{\rm out}}}}{{{10}^{-3}}\;{\rm c}{{{\rm m}}^{-3}}} \right)}^{-1}}\frac{x\;G(1)}{G(\sqrt{x})}\;{\rm kpc},\end{eqnarray} \tag{ 3 }$$

where ${{n}_{e,{\rm out}}}$ is the pre-shock electron density, $x={{T}_{{\rm in}}}/{{T}_{{\rm out}}}$, ${{T}_{{\rm in}}}$ is the post-shock gas temperature, ${{T}_{{\rm out}}}$ is the pre-shock gas temperature, $G(y)=[\phi (y)-y{{\phi }^{\prime }}(y)]/2{{y}^{2}}$, and $\phi (y)$ is the error function (Spitzer 1962). When calculating ${{\lambda }_{{\rm in}\to {\rm out}}}$, pre-shock temperatures were taken from the projected temperature profiles just outside the shock fronts. Post-shock temperatures were calculated from the pre-shock temperatures by applying the Rankine–Hugoniot jump conditions, rather than taken from the observed projected post-shock temperatures. We note that this is a conservative assumption from the perspective of finding a mean free path that is smaller than the shock width, as it gives a higher post-shock temperature (that is not diminished by projection effects), which gives a larger value for ${{\lambda }_{{\rm in}\to {\rm out}}}$.

The relevant particle mean free path is given for each shock front in Table 2. The uncertainty in the mean free path is dominated by the uncertainty in the pre-shock temperature, which is on the order of 10%, leading to a mean free path uncertainty of roughly 20%. In all cases, the particle mean free path is significantly smaller than the shock width, on the order of 10 times smaller than the 90% confidence lower limit for the 1 and 10 kpc shocks (where lower limits on the shock widths can be placed). For a weak shock, the shock width is expected to be roughly $w\approx \lambda /(M-1)$, where λ is the effective particle mean free path (Landau & Lifshitz 1987; McNamara & Nulsen 2007), or 1.4–2 times λ for the weak shocks we consider here. Thus, the apparent shock widths are too large to be explained by particle diffusion alone. We note that this result likely precludes using these shock fronts to place accurate constraints on the electron-ion equilibration timescale, as has been done, e.g., by Russell et al. (2012).

As noted in Section 5, a measured finite shock width may arise due to deformations in and broadening of a shock front as it propagates through a "clumpy" and/or turbulent ICM. Nulsen et al. (2013) provide an estimate of the expected shock width due to turbulence as a function of radius as the shock propagates through a uniformly turbulent ICM. The shock width (w, defined as the rms of the displacement of sections of the shock front due to turbulence) is determined by the distance traveled (${{r}_{s}}$), the shock speed ${{v}_{s}}$, the coherence length of the turbulence ($\ell$), and the rms turbulent speed (${{\sigma }_{t}}$). Here, we can invert this relation to find the turbulent speed implied by the observed shock width:

$$\begin{eqnarray}&&{{\sigma }_{t}}\approx \frac{w\;{{v}_{s}}}{\sqrt{{{r}_{s}}\;\ell }}.\end{eqnarray} \tag{ 4 }$$

The rms turbulent speed implied by measurements of each shock front we consider is given in Table 2. Following Nulsen et al. (2013), we take the coherence length to be $\ell \sim 0.1r$ (Rebusco et al. 2005). Remarkably, the implied turbulent speed is consistent with roughly 70 km s⁻¹ in each case. We note that, in principle, the turbulent speeds calculated here are upper limits on the true ICM turbulence, as this calculation does not include the contribution of density inhomogeneities to shock broadening, projection effects, and the effects of an irregular shock front geometry (which is approximated as spherically symmetric in sectors), all of which will act to increase the apparent shock width. However, as the calculation is only approximate, our results should not be taken as hard upper limits.

The implied turbulent velocity of roughly 70 km s⁻¹ is reasonable compared with results from simulations (100–300 km s⁻¹; Lau et al. 2009) and observations (e.g., de Plaa et al. 2012; Fabian et al. 2013; Zhuravleva et al. 2014). While this value is on the low end of the reported range, we note that N5813 is a relatively isolated galaxy group (R11), which has a regular morphology and shows no signs of a recent major merger. Thus, the level of merger-driven turbulence is expected to be low. Additionally, we note that observational measurements will naturally be biased toward detecting systems with larger turbulent velocities. In several cases only upper limits can be placed, with all reported upper limits significantly larger (a few to several 100 km s⁻¹) than the ∼70 km s⁻¹ detection we report here (Sanders et al. 2011; Bulbul et al. 2012b; Fabian et al. 2013). Finally, this value is consistent with the range of turbulent velocities of 43–140 km s⁻¹ recently reported by Zhuravleva et al. (2014) for the (relatively low mass) Virgo cluster.

By considering the effects of resonant scattering on Fe xvii line ratios, de Plaa et al. (2012) constrain the turbulent velocity in N5813 to be $140\lt {{V}_{{\rm turb}}}\lt 540$ km s⁻¹, significantly larger than our value of 70 km s⁻¹. However, we do not view this as a serious discrepancy for several reasons. First, as mentioned above, the limits provided by Equation (4) are expected to be rough estimates. In particular, the coherence length $\ell$ is unknown and the only hard constraint is that $\ell \lt r$. Second, as noted by the authors themselves, the atomic data used to calculate the line ratios in de Plaa et al. (2012) suffer from significant systematic uncertainty. For example, by comparing results derived with the spex ¹¹ code versus AtomDB ¹² , they find inferred turbulent velocities for N5813 that differ by almost a factor of 2. Finally, we note that there are plausible physical origins for the higher level of turbulence found by de Plaa et al. (2012). For example, they place constraints for the central $r\lesssim 5$ kpc only, where the turbulent velocity may be larger than the total ICM average due to the influence of the AGN, whereas the estimates given here are averages over the entire region interior to each shock. The trend of the inferred turbulent speed decreasing with radius (Table 2), although not statistically significant, is consistent with this picture. It is interesting to note that the turbulent velocity inferred from the width of the central r ≈ 1 kpc shock is consistent with the range given by de Plaa et al. (2012), although, as noted above, the errors are large and the effects of the larger relative PSF are unclear. Additionally, we note that the resonant line ratio may also be affected by bulk flows in the ICM (e.g., the radial gas flows behind the shock fronts), possibly leading to an over-estimate of the turbulence.

We conclude that, while the measured shock widths are too large to be explained by Coulomb diffusion, they are consistent with arising from the shocks propagating through a turbulent ICM. The implied ratio of turbulent to thermal pressure support, estimated as the ratio of the turbulent energy to the thermal energy (de Plaa et al. 2012), is on the order of a few percent ($\lesssim 5\%$) for a turbulent speed of 70 km s⁻¹. For each shock we consider, the local thermal diffusion coefficient $\nu \approx \lambda {{c}_{s}}$, where ${{c}_{s}}$ is the sound speed, is significantly smaller (by up to an order of magnitude) than the implied turbulent diffusion coefficient $\eta \approx \ell {{\sigma }_{t}}$. Thus, turbulent diffusion dominates the microscopic transport in the ICM (at least in this case, and likely in other similar-size groups, which are expected to have similar bulk ICM properties).

Although the X-ray cavities must be close to pressure equilibrium with their surroundings, radio observations show that generally the pressure support from relativistic non-thermal electrons, under the assumption of equipartition, is insufficient to balance the thermal pressure in the ICM (Bîrzan et al. 2008). Additional pressure support might be provided by heavy non-thermal particles, or by high-temperature thermal gas, with deep X-ray observations providing constraints on the contribution from the latter (e.g., Schmidt et al. 2002; Blanton et al. 2003). Here we place similar limits for the low-mass galaxy group N5813.

To place lower limits on the temperature of thermal gas in the cavities, we use a method similar to that employed in previous studies of other systems (e.g., Sanders & Fabian 2007; Rafferty et al. 2013). Due to the complex thermal structure of the gas, we do not use comparison regions outside of the cavities as in Sanders & Fabian (2007). Rather, we place limits based on fits to spectra extracted from regions at the cavity centers, where the contrast is highest and the relative contribution from high-temperature gas inside the cavities is expected to be the largest.

Each spectrum is fit with two absorbed apec model components. The first component models the contribution from the total projected emission, while the second models hotter thermal emission from within the cavities. To derive conservative estimates, the temperature, abundance, and normalization of the first component are allowed to vary. We find that within the cavity regions this single thermal component is sufficient to adequately model the total projected emission (additional thermal components do not significantly improve the fits and are not well constrained). For the second component, the temperature is fixed at some value, the abundance is tied to that of the first component, and the normalization is allowed to vary. For an assumed cavity volume, the upper limit on the normalization gives an upper limit on the density. If we further assume that the cavity is in pressure equilibrium with its surroundings, that the pressure is equal to the average azimuthal pressure at the cavity radius (Figure 11), and that the thermal gas dominates the pressure support inside the cavities, an upper limit on the volume filling fraction of the gas in the cavities (for the given temperature) is obtained. A lower limit on the gas temperature is obtained by increasing the fixed temperature of the hot gas until the upper limit on the volume filling fraction is equal to unity.

We concentrate on the cavities with the largest contrast and the most regular morphology. Since they are closely spaced (and therefore in a similar environment), we fit the inner cavities together, assuming a cylindrical geometry with a radius of $r=0.408$ (0.488) kpc and a depth of $d=0.689$ (0.815) kpc for the northeastern (southwestern) central cavity. For the intermediate cavities, only the southwestern cavity is considered (with $r=1.82$ kpc, $d=3.01$ kpc) as it has a larger contrast and more regular shape compared with the northeastern intermediate cavity, which may be split into two cavities and may be connected to the northeastern central cavity (R11). Similarly, we choose the northeastern outer cavity over its southwestern counterpart due to its larger contrast and more regular shape, and we assume an oblate spheroidal geometry with major axis $a=5.89$ kpc and minor axis $b=1.90$ kpc in the plane of the sky. We find limits on the temperature of any volume filling thermal gas in the cavities of $\gt 3.4$, $\gt 16$, and $\gt 4.5$ keV for the inner, intermediate, and outer cavities, respectively (for comparison, the limit placed by Rafferty et al. 2013). Our most stringent constraint of $kT\gt 16$ keV for the intermediate cavities implies a gas density in the cavities that is at least 20 times lower than the azimuthal average at that radius (${{n}_{e}}\lt 9\times {{10}^{-4}}$ cm⁻³).

We note that this method relies on several assumptions. For instance, the hot temperature component is assumed to arise solely from gas contained within the cavities (it may be the case that hot thermal gas supplies only part of the pressure support, while, e.g., heavy non-thermal particles or magnetic fields make up the difference). Additionally, the azimuthally averaged pressure profile, which is used to estimate the mean external pressure at each cavity, is derived under the assumption of spherical symmetry (Section 4.2.1), which is clearly violated at some level for this system (Figure 2). Furthermore, some assumptions about the cavity geometries (or, more specifically, their volumes) must be made to calculate densities and pressures. The largest source of error is likely from the latter assumption. We find that an assumed uncertainty in cavity volume of 20% leads to an uncertainty of ∼10% on the temperature limit (a similar result was found by Rafferty et al. 2013).

The residual image shown in Figure 2 reveals a "channel" of decreased surface brightness to the north, apparently connected to the NE outer cavity. To test the significance of this feature, we extracted the surface brightness profile in azimuthal bins roughly 9° wide and 64'' long across the channel. The profile, shown in Figure 15, shows a clearly significant dip across the channel, indicating that this is a real feature.

Zoom In Zoom Out Reset image size

Since the northern channel appears to connect to the NE outer cavity, one interesting possibility is that this channel is the result of energy (presumably mostly in the form of non-thermal particles) leaking through the bright cavity rim and heating the ICM. Indeed, the temperature maps shown in Figure 5 show relatively high temperatures in the region of the northern channel, although the correlation between the high-temperature region in the temperature map and the channel in the residual map is not exact. Although cavities are expected to heat the ICM (see Section 6.1), such a direct detection of this heating has not yet been observed. We therefore extracted temperatures both in (wider) azimuthal bins across the northern channel and in identical elliptical regions within and to either side of the channel (at roughly the same distance from the center). In neither case did we find a significant temperature enhancement associated with the channel, with temperatures of all regions agreeing within roughly 1σ. We conclude that, although it is possible that the northern channel is a result of energy leaking from the NE cavity and heating the ICM, causing the gas in this region to expand to pressure equilibrium (and subsequently cool adiabatically) and hence the density (and surface brightness) to drop, we find little evidence to support this interpretation. In particular, we cannot rule out the possibility that this feature is simply due to large-scale structure in the ICM, e.g., associated with previous mergers or bulk motion of the group, even in this relatively relaxed system. Another possibility is that this deficit is a leftover cavity from an even older, fourth outburst that lies off of the axis defined by the other three cavity pairs, although the fact that it is at least partially within the 30 kpc shock front is difficult to reconcile with this interpretation. Deep, low-frequency radio observations would be useful both to detect the putative low-energy non-thermal particles in the NE outer cavity and to see if this emission extends out along the northern channel, consistent with it arising from a leaking cavity.