Spatially resolved Spitzer-IRS spectral maps of the superwind in M82 Free

Abstract

We have mapped the superwind/halo region of the nearby starburst galaxy M82 in the mid-infrared with Spitzer − IRS. The spectral regions covered include the H₂ S(1)–S(3), [Ne ii], [Ne iii] emission lines and polycyclic aromatic hydrocarbon (PAH) features. We estimate the total warm H₂ mass and the kinetic energy of the outflowing warm molecular gas to be between M_warm ∼ 5 and 17 × 10⁶ M_⊙ and E_K ∼ 6 and 20 × 10⁵³ erg. Using the ratios of the 6.2, 7.7 and 11.3 μm PAH features in the IRS spectra, we are able to estimate the average size and ionization state of the small grains in the superwind. There are large variations in the PAH flux ratios throughout the outflow. The 11.3/7.7 and the 6.2/7.7 PAH ratios both vary by more than a factor of 5 across the wind region. The northern part of the wind has a significant population of PAH's with smaller 6.2/7.7 ratios than either the starburst disc or the southern wind, indicating that on average, PAH emitters are larger and more ionized. The warm molecular gas to PAH flux ratios (H₂/PAH) are enhanced in the outflow by factors of 10–100 as compared to the starburst disc. This enhancement in the H₂/PAH ratio does not seem to follow the ionization of the atomic gas (as measured with the [Ne iii]/[Ne ii] line flux ratio) in the outflow. This suggests that much of the warm H₂ in the outflow is excited by shocks. The observed H₂ line intensities can be reproduced with low-velocity shocks (v < 40 km s⁻¹) driven into moderately dense molecular gas (10² < n_H < 10⁴ cm⁻³) entrained in the outflow.

1 INTRODUCTION

Galactic outflows, ‘superwinds’, are driven by hot gas and momentum generated by the combined mechanical and radiative energy and momentum output from stellar winds and supernovae (e.g. Heckman, Armus & Miley 1990). They have been invoked as a source for the heating and metal-enrichment of both the intracluster and intergalactic medium (Adelberger et al. 2003; Loewenstein 2004), as a major driver of galaxy evolution (e.g. Croton et al. 2006), and as playing a crucial role in establishing the mass–metallicity relation among galaxies (Tremonti et al. 2004). Superwinds have also been hypothesized to be the mechanism through which a dust enshrouded ultraluminous infrared galaxy (ULIRG) transitions to an optically identified quasar (e.g. Sanders et al. 1988).

The dynamical evolution of a starburst-driven outflow has been extensively discussed (e.g. Chevalier & Clegg 1985; Suchkov et al. 1994; Wang 1995; Tenorio-Tagle & Muñoz-Tuñon 1998; Strickland & Stevens 2000). Superwinds are generated when the kinetic energy in the outflow from massive stars and supernovae is thermalized, generating a region of very hot (T ∼ 10⁸ K) low-density gas in the (interstellar medium) ISM of a starburst galaxy (Chevalier & Clegg 1985). This hot gas will expand, sweeping-up and entraining ambient gas and creating a hot cavity within the overall structure of the ISM. This region of hot gas is highly overpressurized relative to its surrounding ambient ISM and will expand along the steepest pressure gradients. As it expands, it sweeps up and shocks ambient material, creating an outflowing ‘superwind’ (Castor, McCray & Weaver 1975; Weaver et al. 1977). Contributing to this overall energy driven outflow is the radiation pressure on dust grains that are likely coupled to the ionized gas within the flow (Murray, Quataert & Thompson 2005).

There is a significant amount of morphological, physical and kinematic evidence for the existence of superwinds in nearby starburst and infrared luminous galaxies (e.g. Heckman et al. 1990). Large- scale optical emission-line and associated X-ray nebulae are ubiquitous in starbursts (Heckman et al. 1990; Lehnert, Heckman & Weaver 1999; Strickland et al. 2004), and these increase in size and luminosity from dwarf starbursts to ULIRGs, where they can be tens of kpc in size (Grimes et al. 2005). In addition, the winds are observed in molecular lines such as CO (Walter, Weiss & Scoville 2002; Bolatto et al. 2013) and blueshifted absorption features (e.g. Na D) seen projected against the starburst nuclei (Heckman et al. 2000; Martin 2005; Rupke & Veilleux 2005) indicating that the neutral, cold ISM is being swept-up in the outflow. Recently, vigorous outflows have also been detected in ULIRGs and LIRGs in molecular absorption lines with Herschel-PACS which suggest very high mass outflow rates (e.g. Fischer et al. 2010; Sturm et al. 2011; Spoon et al. 2013; Veilleux et al. 2013). ULIRGs which have higher AGN bolometric fractions appear to have higher terminal velocities and shorter gas depletion time-scales than those with lower AGN bolometric fractions, suggesting that at least in ULIRGs, the AGN power is related to the feedback on the cold molecular gas (Veilleux et al. 2013; Cicone et al. 2014).

M82 is the brightest, nearest, and best-studied starburst galaxy with an outflow (Heckman et al. 1990; Shopbell & Bland-Hawthorn 1998; Lehnert et al. 1999). At a distance of about 3.5 Mpc (Jacobs et al. 2009), M82 affords an unparalleled opportunity to study a dusty outflow in great detail. The outflow in M82 has extensive and complex extraplanar polycyclic aromatic hydrocarbon (PAH) and H₂ emission, evident from the Spitzer-Infrared Array Camera (IRAC) and Spitzer-Infrared Spectrograph (IRS) observations presented in Engelbracht et al. (2006). Extensive PAH emission is also seen at 3.3 μm in the AKARI data Yamagishi et al. (2012). PAHs are thought to be the origin of strong, broad emission features arising in the mid-infrared spectra of star-forming galaxies (e.g. Brandl et al. 2006; Smith et al. 2007a). M82 has also been observed in the [C ii] and [O i] emission lines by Contursi et al. (2013) who find a low-velocity outflow in the atomic gas. Complex extraplanar warm H₂ knots and filaments extending more than ∼5 kpc above and below the galactic plane of M82 have been imaged in the near-infrared (NIR) using narrow-band filters by Veilleux, Rupke & Swaters (2009). Using wide-field imaging, they found that H₂ is widespread, suppressed relative to PAHs in the fainter regions, and concluded that H₂ is not a dynamically important component of the outflow. However, Veilleux et al. (2009) could only measure the very hot (∼1000 K) H₂ gas, and the mass of very hot H₂ gas is much smaller than the mass of the H₂ gas probed by this study, as we show in Section 3.3.

In this paper, we study the M82 outflow with Spitzer–IRS. With Spitzer/IRS spectroscopy, we can study in detail all the PAH features, derive the grain properties through analysis of the PAH ratios, and compare these to the warm (200 K) molecular gas and the ionized gas (via the strong Ne lines) as a function of position above and below the plane of M82 for the first time. The rotational lines of H₂, observable in the infrared, unlike the NIR ro-vibrational lines, probe most of the warm molecular gas mass. We then use the ratio of pure rotational H₂ lines to constrain the physical and excitation conditions of the gas. We use PAH 11.3/7.7 and 6.2/7.7 ratio maps to study PAH sizes and ionization across the outflow, [Ne iii]/[Ne ii] to study the variation of the hardness of the ionization field, H₂/PAH to diagnose the origin of the H₂ emission, and the H₂/[Ne ii] and PAH/[Ne ii] line flux ratios to help understand the survivability of the PAH grains in the wind. For the first time we can make robust maps of these ratios throughout the outflow of M82, and compare all of these to the spatially-resolved X-ray emission to understand the properties of the dust and molecular gas and ionized gas in the wind.

The paper is organized as follows: in Section 2 we describe the observations and the data reduction; in Section 3 we present and analyse the spectral maps built from the observations, using the PAH 11.3/7.7 and 6.2/7.7, [Ne iii]/[Ne ii], H₂/PAH, H₂/[Ne ii], and PAH/[Ne ii] line ratios, and also H₂ excitation diagrams to diagnose the properties of the ionized gas, PAHs and the H₂ emission; in Section 4 we diagnose the origin of H₂ emission using X-ray images and shock excitation models; and in Section 5 we present the conclusions.

2 OBSERVATIONS AND DATA REDUCTION

The outflow of M82 was observed with Spitzer–IRS between 2006 November 10–14 and 2008 May 9–22, using all low-resolution modules. In Fig. 1 we show the footprints of the IRS observations on the outflow of M82. We extracted only the maps that contain both modules SL1 (7.8-14 μm) and SL2 (5.5-7.8 μm) or LL1 (20-35 μm) and LL2 (14-20 μm). In total we mapped an area of approximately 5 arcmin² above and below the plane of M82.

We used CUBISM (Smith et al. 2007a) to build SL and LL spectral cubes from the Spitzer–IRS Basic Calibration Data (BCDs). In order to create and subtract a local background from the wind spectra, we observed a region 20 kpc North of the M82 disc with both IRS modules. This region was observed with the same parameters (integration time per pixel and cycles per position) that were used for the wind. After verifying that no line, PAH or dust continuum was visible in the off field, we subtracted an average image created from this region from the wind data at each location and for each IRS module before building the spectral cubes.

We used PAHFIT (Smith et al. 2007b) to create maps of spectral lines and PAH features. PAHFIT is a spectral fitting routine that decomposes IRS low-resolution spectra into broad PAH features, unresolved line emission, and continuum emission from dust grains. PAHFIT allows for the deblending of overlapping features, in particular the PAH emission, silicate absorption and fine-structure atomic, and warm H₂ emission lines. Because it performs a simultaneous fit to the emission features and the underlying continuum, PAHFIT works best if it can fit combined SL and LL spectra across the full 5-40 μm range. In order to provide a stable long wavelength continuum from which to estimate the true PAH emission, including the broad feature wings, we formed a combined SL + LL cube. In order to investigate the spatial variations of the H₂/PAH ratio at the maximum resolution we decided to rotate and regrid the LL cubes to match the orientation and pixel size of the SL cubes (1.85 arcsec × 1.85 arcsec). Thus we provide a stable long-term extrapolation of the LL fit to the SL fit. By calculating the flux in the new LL cubes using a bilinear interpolation (in surface brightness units), we assure that the flux over the native SL pixel scale is conserved.

Finally, we applied a small correction factor to each individual SL spectrum, necessary to match the continuum in the corresponding LL spectrum. This scale factor ranges between 0.65 and 0.99, and was calculated by comparing pixels in the spectral overlap region between SL and LL at approximately 14 μm. This scale factor is applied before running PAHFIT on the combined SL + LL cubes. Such scaling is commonly observed in spatially extended sources as a result of the different entrance slit dimensions of the SL and LL spectrograph modules. In Fig. 2 we present an example of the spectral decomposition performed at every location of the map with PAHFIT. The output of PAHFIT is saved for each feature, and is used to create all the maps (PAH, H₂, [Ne ii] and [Ne iii]) at the native resolution of the SL data as shown in Fig. 3. Although care was taken to avoid the FIR-bright nucleus of M82 during the spectral mapping, the LL1 data exhibit a ‘ripple’ pattern that varies spatially and spectrally, and the LL2 data have a slightly elevated continuum level in some pixels (compared to SL1). We believe both of these to be the effect of light scattered into the nuclear IRS entrance slits. Therefore, we have opted not to use the LL1 data in the fitting and analysis, except for the derivation of an upper limit on the H₂ S(0) line flux, calculated as the 3σ noise level at 28 μm. We have excluded regions of the LL2 data cubes where the implied scale factor in the continuum between SL1 and LL2 is greater than 2.0.

3 RESULTS

3.1 Morphology of the dust, ionized gas and warm H₂ emission

In Fig. 3 we present the total PAH emission map calculated by summing over all the PAH bands between 5 and 14 μm, the [Ne ii]12.8 μm flux map derived from the SL data cubes, the total H₂ map derived from the SL and LL2 cubes, and the [Ne iii]15.6 μm emission map. We find widespread emission with significant structure in all of the line maps. To the north, both PAH and [Ne ii] exhibit two prominent, broad radial filaments. To the south, the PAH map shows a complex structure with a bubble, or loop. This loop is also seen in the [Ne ii] map. The overall morphologies of the warm H₂ and Hα are generally very similar, as seen in the NIR by Veilleux et al. (2009). Maps of the cold molecular gas made in CO (J = 1–0) by Salak et al. (2013) are at lower spatial resolution, but also show extended, extraplanar gas.

We can study the physical properties of the ionized gas and PAHs by mapping the ratios between ionized lines such as [Ne iii]/[Ne ii] and between the different PAH bands. If the ionized gas is photoionized, the Ne line ratio is sensitive to the hardness of the radiation field and the ionization parameter. The Ne ratios observed in the outflow are standard values for the hardness of the radiation field found in starburst galaxies. The ratio of the 6.2, 7.7 and 11.3 μm PAH features can be used to estimate the ionization state and the size of the small dust grains (e.g. Draine & Li 2001). Even the overall absence of PAH emission may be related to the ionization state of the gas since we do not expect gas emitting strongly in [Ne iii] to show much PAH emission, as the small grains are easily destroyed by energetic photons (e.g. Engelbracht et al. 2006; Madden et al. 2006; Wu et al. 2006). The 11.3 μm/7.7 μm ratio is high when PAHs are more neutral and larger, while the 6.2 μm/7.7 μm PAH ratio is high when PAH grains are smaller. However, we caution that variations in the PAH line ratios are also influenced by the structural properties and elemental abundances in small particles giving rise to the PAH features (e.g. Yamagishi et al. 2012).

In Fig. 4, we present maps of the PAH 11.3/7.7 and 6.2/7.7 ratios, the H₂/PAH ratio and the [Ne iii]/[Ne ii] ratio. The H₂/PAH ratio map was produced using a H₂ [S(1)–S(3)] total flux map and PAH (7.7 + 8.6 μm) flux map. Although there is a great deal of structure in the PAH ratio maps, the general trend is that the 11.3/7.7 and the 6.2/7.7 ratios increase as we move away from the starburst disc. The PAH 11.3/7.7 is enhanced by a factor of ∼4 and the PAH 6.2/7.7 is enhanced by a factor of ∼3 at the northern edge of the outflow compared to the regions near the centre of the starburst disc. There is an enhancement of the 11.3/7.7 and 6.2/7.7 ratios by ∼2 in the southeast edges of the outflow.

The H₂/PAH ratio varies between 0.01 and 1, and the [Ne iii]/[Ne ii] ratio varies between 0.2 and 0.9. In the region between the PAH emission ‘fingers’ (as particularly evident in the Northern outflow region in Fig. 3), both maps show an enhancement of the ratio.

In Fig. 4, we also mark the 12 regions from which we have extracted 5-20 μm spectra and performed multitemperature fits to the H₂ lines. We show these spectra in Fig. 5. These regions were chosen to represent the full range of H₂/PAH ratios and to approximately sample the full spatial extent of the outflow. For comparison, we also show an integrated spectrum of the centre of M82 (Fig. 5; Beirão et al. 2008), which we use as a reference. In addition to the 12 regions, we also show a composite spectrum of the entire IRS map in Fig. 5. We measure the fluxes of the PAH features between 5 and 14 μm, the fine-structure lines [Ne ii] and [Ne iii], and the rotational excitation lines of H₂ from S(0) to S(3) and present the results in Tables 1 and 2. Regions 1, 5, 6, 7 and 9 are the most similar to the reference spectrum in slope and the relative strength of PAH features, which is not surprising since they are located close to the centre. PAHs are also strong in regions 3, 8, 10, 11 and 12, but their spectrum also exhibits a flatter continuum, probably indicating the presence of colder dust. Region 3 is the most distant from the central starburst and it has the weakest PAH emission. Region 4 exhibits a flat continuum, but a comparatively strong PAH 11.3 μm feature.

As we have seen in the maps, the spectra show strong variations of the H₂/PAH ratios between the regions. The H₂S(1) line is very prominent in the spectrum of Region 3 compared to the PAH bands, while in the regions close to the disc, such as Regions 7 and 9, this line is much fainter compared to the PAH features. Indeed, using the values in Tables 1 and 2, we calculate in Region 3 a very high H₂/PAH ∼ 0.08, while in Region 9, H₂/PAH ∼ 0.01.

In Fig. 6, we present a plot of the PAH 11.3/7.7 ratio as a function of the 6.2/7.7 PAH ratio. The range of 6.2/7.7 PAH ratios in the outflow is much larger and the values are, on average, smaller than the disc. The outflow exhibits 6.2/7.7 and 11.3/7.7 ratios ranging between ∼0.2 and 0.4. On average, the bulk of the superwind shows larger, more ionized grains than in the starburst disc, although there are regions where the ratios overlap. As an example of the different PAH ratios found in the outflow, we present in the upper-right corner of Fig. 6, two spectra from different regions of the outflow, and an average spectrum of the disc. In general, the outflow has stronger 17 μm PAH and H₂ emission than the SB disc, and the spectrum in the region 2 kpc north of the disc (the blue spectrum shown in Fig. 6) has higher 17 μm feature and H₂ S(1) flux, indicating the prevalence of larger PAH grains and enhanced H₂ emission further away from the disc. In summary, we see a much larger range in the PAH ratios in the outflow than is seen in the starburst disc. The average spectrum of the wind has smaller 6.2/7.7 ratios and larger 11.3/7.7 ratios than the average spectrum of the starburst, suggesting the grains may be larger and more ionized, on average, than those in the starburst disc.

3.2 H₂/PAH variations in the M82 wind

In Fig. 7, we plot the variation of H₂S(1)–S(3)/PAH7.7 + 8.6 μm as a function of [Ne ii] surface brightness. Since the [Ne ii] emission arises predominantly from the starburst disc, this traces the relative H₂ and PAH emission as a function of distance from the starburst. This ratio has been previously used by Ogle et al. (2010) and Guillard et al. (2012) to show the relative importance of the UV versus mechanical heating of the gas. We can see that the H₂/PAH ratio decreases with increasing [Ne ii] flux. The vast majority of the mapped regions lie well above the values seen in star-forming galaxies in the local Universe (the dashed lines; see Roussel et al. 2007). The H₂/PAH ratio values span two orders-of-magnitude, 0.01–1. Most of the M82 wind points have H₂/PAH between 0.01 and 0.10, but ∼4 per cent have extreme values, H₂/PAH >0.5. The fact that the H₂/PAH ratios seen in the wind are well above those seen in the SB disc, the normal galaxies observed by Roussel et al. (2007), and PDR models of Guillard et al. (2012, see section 4) suggests that the H₂ is not solely excited by stellar radiation from the starburst.

We can also investigate possible correlations of the H₂/PAH ratio with diagnostics of the atomic gas and the size and ionization state of the grains. In Fig. 8, we show H₂/PAH as a function of the PAH 11.3/7.7 ratio, the PAH 6.2/7.7 ratio, the H₂ S(3)/S(1) ratio, and the [Ne iii]/[Ne ii] ratio. We see that in the disc the H₂/PAH ratio increases with increasing PAH 11.3/7.7 ratio. While the trend of increasing H₂/PAH in the wind is similar to that seen in the disc, the wind values extend to much larger values of both 11.3/7.7 PAH and H₂/PAH than seen in the starburst. The scatter to high H₂/PAH and a much larger range in the 6.2/7.7 ratio is seen in the upper-right panel of Fig. 8 as well. In the lower left, we observe an anticorrelation between H₂/PAH and S(3)/S(1) in the disc as warmer gas has lower H₂/PAH. The wind points seem consistent with this trend but the uncertainties are large. This general lack of very warm molecular gas and high H₂/PAH ratios in the wind may be a function of inability of the MIR diagnostic features to probe the hottest gas. The [Ne iii]/[Ne ii] also increases with H₂/PAH in the disc, and the wind points to the North seem consistent with the trend, but clearly there are many wind points in a few regions in the northern outflow which have much higher H₂/PAH ratios than seen in the disc at a given [Ne iii]/[Ne ii] value.

To further help us understand what drives the variations of the H₂/PAH flux ratios, we present maps of H₂/[Ne ii] and PAH/[Ne ii] in Fig. 9. In general, the regions with high H₂/PAH correspond to high PAH/[Ne ii], suggesting that it is not simply a decreasing PAH that causes the rise in H₂/PAH off the plane of M82. Only in the northern edge of the map do we see a lower PAH/[Ne ii] ratio, which might indicate reduced PAH heating efficiency. It is important to note, however, that while the PAH/[Ne ii] ratio does change by factors of a few across the map, this is much less than the variation observed in H₂/PAH. Therefore the latter cannot be solely explained by a change in the strength of the PAH emission alone.

3.3 H₂ masses and temperatures

With the rotational H₂ emission line fluxes in Table 2, we can estimate the gas temperature, column density, and mass distributions, assuming that the ortho/para ratio is in local thermal equilibrium (LTE). In Fig. 10, we plot the excitation diagrams for all regions, except 2 and 9, which did not have sufficient signal-to-noise for fitting the data with these sorts of models. Using the Boltzmann equation, we fit the H₂ lines with two components, each with different temperatures. We then derive the total column density, H₂ temperatures and mass distributions. In fitting the data we set the minimum allowed temperature to 100 K and treated every upper limit as a 3σ detection. The S(1) line is the line with the highest signal-to-noise. Therefore, when we remove the S(0) line as a detection, we obtain essentially the same results.

The diagrams for Regions 1, 3, 4, 5, 7 and 11 look very similar in shape, derived from very similar temperatures for each component. In Region 8, both components have similar temperatures, while in Region 10, there is a large difference in temperature between the two components. We note that Regions 8 and 10 have the most discrepant S(0) upper limits relative to the other fits, making the mass estimates particularly uncertain. In Table 3, we present the properties of the warm H₂ gas in these regions of interest. In all the regions the temperature of the warmer component is above 500 K, especially warm when compared with typical star-forming regions (e.g. Roussel et al. 2007). The column density of the bulk of the warm H₂ gas ranges by a factor of ∼30. The mass ratio between the warmer and colder component varies between 0.0003 and 0.02, with the lowest value observed in Region 8 and the highest in Region 1. The region with the most luminous H₂ emission is Region 8, where the warmer component accounts for ∼95 per cent of the luminosity, whereas the next most luminous region, 11, the colder component accounts for ∼65 per cent of the luminosity.

We also integrated the emission over the entire observed nebula surrounding M82 and measured the line fluxes, which are included in Table 2. The total warm H₂ mass derived from the integrated emission is ∼4.8 × 10⁶ M_⊙. For completeness, we note that Veilleux et al. (2009) finds a mass of H₂ at a temperature T>1000 K (thus hotter than our ‘hot’ gas) of 1.2 × 10⁴ M_⊙, based on NIR observations. The relatively low mass is not surprising since the NIR ro-vibrational lines are probing very high excitation gas which emits very efficiently in these lines. The mass of warm molecular gas in the wind that we are detecting with the IRS is thus more than a factor of 100 more than is estimated from the NIR line maps. We have also estimated the H₂ masses from fitting the H₂ line fluxes with molecular shock models (see Section 4.3).

Considering a cold molecular gas mass of ∼10⁹ M_⊙ (Walter et al. 2002; Salak et al. 2013), we calculate a ratio of warm to cold molecular gas in the outflow of about 0.01. Therefore the warm molecular gas we measure here is a relatively minor component of the total molecular mass in the ISM. Care must be taken when making this comparison as the emission regions for the CO and H₂ do not overlap completely. In fact, the most extended CO emission is only within a few arcminutes of the plane and much less extended (Walter et al. 2002; Salak et al. 2013), although new studies trace the CO-emitting gas to larger distances (Leroy et al., in preparation). However, we note that most of the estimated mass for the warm molecular gas in our observations also comes from regions which lie close to the disc (and likely the base of the wind).

4 DISCUSSION

4.1 H₂ heating processes

In the scenario of thermal excitation of H₂, three mechanisms dominate the heating: (1) UV radiation from starburst, (2) X-rays from the starburst and wind, and/or (3) shocks induced by the outflow. Following Ogle et al. (2010), we use the H₂-to-PAH luminosity ratio to investigate the contribution of UV photons to the total heating of the H₂ gas in the outflow. In starburst galaxies like M82, we expect the H₂ emission to be predominantly powered by UV radiation from young stars but undoubtedly, processes such as low-velocity shocks and turbulent dissipation contribute to the H₂ emission. It should also be noted that H₂/PAH ratios on the order of 0.01 have also been seen from translucent clouds in the cold, diffuse ISM of the Milky Way, which are thought to be excited by UV photons (Ingalls et al. 2011).

There is a clear offset of the wind points from those in the SB disc in terms of H₂, PAH and ionized atomic emission (Figs 7, 8, and 9). The data in regions of the outflowing gas show a large range in both [Ne iii]/[Ne ii] and H₂/PAH, but are clearly offset to higher H₂/PAH at a given [Ne iii]/[Ne ii] ratio compared to regions within the SB disc, suggesting an additional excitation mechanism at work. In fact, models suggest that the large values of the H₂/PAH ratios observed in some regions (>0.1) cannot be explained by photoionization (Guillard et al. 2012).

The H₂ gas could be indirectly excited by X-rays emitted from the expanding bubble of hot and tenuous gas generated by the intense star formation within M82 (Begelman & Cioffi 1989), high-velocity shocks in the entrained gas, X-rays due to inverse Compton scattering of IR photons. In M82, the X-ray emission covers approximately the same region as the optical emission line filaments (Watson, Stanger & Griffiths 1984; Shopbell & Bland-Hawthorn 1998) and radio emission (Seaquist & Odegard 1991). For example, the high H₂/PAH in the northern tip ∼5 kpc from the disc, appears to be associated with a region with X-ray emission (Fig. 11).

We can use the total H₂S(1)–S(3) luminosity and the total X-ray luminosity in the outflow to calculate the H₂-to-X-ray luminosity ratio (Stevens, Read & Bravo-Guerrero 2003). The H to H₂ transition can occur at relatively low column densities, N_H ≈10²¹ cm⁻². At such low columns, even relatively soft X-rays can penetrate. Our estimates for the M82 outflow are L(H₂S(0) − S(3))/L_X(2–10 keV) ∼ 1 (Stevens et al. 2003) and L(H₂S(0) − S(3))/L_X(0.2 − 4keV) ∼ 0.5 (Watson et al. 1984; Strickland & Heckman 2007), over the northern and southern IRS map. Using XDR models (e.g. Maloney, Hollenbach & Tielens 1996) we can estimate whether the H₂ line emission could be powered by X-ray heating, following the procedure described in Ogle et al. (2010). In these models, 30–40 per cent of the absorbed X-ray flux goes into gas heating via photoelectrons, and the cooling by H₂ rotational lines in XDR models is ∼2 per cent of the total gas cooling for a gas temperature of 200 K. At this temperature, the ratio of the first four rotational lines to the total rotational line luminosity is L(H₂ 0–0 S(0)–S(3))/ L(H₂) = 0.6. Combining the above factors, we estimate a maximum H₂ to X-ray luminosity ratio of L (H₂ 0–0 S(0)–S(3))/L_X (2–10 keV) <0.01 (Guillard et al. 2012b). This upper limit on the H₂-to-X-ray luminosity is conservative because it assumes that all the X-ray flux is absorbed by the molecular gas, which is likely not the case (since the medium in which X-rays propagate is inhomogeneous). Therefore, we conclude that X-rays are not the dominant source of heating of the H₂ gas in the M82 wind, since the observed ratios are well above those found in XDR models. This is also seen in some extreme environments, such as powerful radio galaxies, where H₂/PAH ∼0.5–1 and the H₂/X-ray ratio also exceeds that expected from XDRs alone (Ogle et al. 2010).

Ambient gas clouds heated by the mechanical energy of the wind and dense shell fragments from the wind bubble carried along by the wind are the two main sources of X-ray emission in the wind (Strickland, Ponman & Stevens 1997). In Fig. 11, we show the contours of soft X-ray emission from the outflow of M82 overlaid on the H₂/PAH and total (5-17 μm) PAH maps. The X-ray and PAH maps show a remarkable correlation, particularly in the southern part of the bi-conical outflow, although both maps are projections of the true 3D structure and the precise locations of the hot gas and dust grains are unknown. Moreover, there are regions of high surface brightness X-ray emission with relatively little PAH emission (in the north along the outer wind axis). We also see that the middle lobe of X-ray emission in the northern outflow, at about 5 kpc from the disc, appears related to the increase of H₂/PAH7.7 + 8.6 μm surface brightness. The morphology of H₂/PAH closely resembles the morphology of the H₂ (2.12 μm)/PAH maps in Veilleux et al. (2009). The regions where H₂/PAH > 0.1 correspond to regions where H₂ (2.12 μm)/PAH > 0.003.

While there is a reasonable spatial correlation between the surface brightness distribution of the soft X-ray and H₂ emission, a correlation is apparent between soft X-rays surface brightness and H₂/PAH ratio distributions only in high surface brightness X-ray emitting regions. There are several possibilities why this might be so, but the most obvious is that both the X-ray and H₂ emission are fundamentally powered by the same mechanism, namely, the very hot plasma that is acting as the ‘piston’ driving the wind and heating clouds entrained in the wind. However, the warm H₂ and hot X-ray emitting plasma should not be co-spatial on small scales, and instead the H₂ is likely seen in projection against the hot wind, along the edges of the cone in dense knots of gas, heated by shocks (see below). Such a situation would also arise in gas with a wide range of densities or through changes as individual clouds are shredded as they are accelerated by the action of the piston (Cooper et al. 2008). This is also shown by the strong correlation between the X-rays and the PAH emission (Fig. 11b), especially in the Southern wind. As the PAHs and H₂ molecules cannot be expected to survive within the 10⁶–10⁷ K X-ray emitting plasma, it is most likely that they lie along the edges of this phase of the outflow. So while the power source is likely the same, the two emission regions, not surprisingly, represent different physical regions and gas states being impacted by the driving very hot plasma (Strickland & Heckman 2007).

4.2 The velocities and kinetic energy of the molecular gas

Kinematics of the gas can be used to determine if the different gas phases observed in the outflowing wind are related. If the observed kinematics of the warm molecular gas as probed by the H₂ rotational lines is similar to that of the cold molecular gas as probed by the CO lines, the IR atomic neutral and ionized emission lines, or the lines of the warm ionized medium, then this might suggest a direct physical relationship between these phases.

Salak et al. (2013) derived a velocity of 160 km s⁻¹ for the CO(1–0) emitting gas for the northern outflow region. However, the farthest extent of the CO emission is only about 1 arcmin north of the disc of M82, while the IRS map extends to 4 arcmin. The velocity derived by Salak et al. (2013) is not very different from the velocities in the neutral and ionized atomic gas in the infrared within the inner wind region <1 arcmin from the disc (Contursi et al. 2013). If we restrict ourselves to the inner wind region within a few arcminutes of the plane of M82, the velocities in the warm ionized gas are only a few 100 km s⁻¹ or less within a couple of arcminutes of the disc of M82 (Shopbell & Bland-Hawthorn 1998). These velocities are consistent with those observed in IR neutral and ionized atomic lines (Contursi et al. 2013). However, farther out of the plane of M82, much higher outflow velocities (∼600–800 km s⁻¹) are observed in Hα for the outflow in M82 (Shopbell & Bland-Hawthorn 1998). Line splitting is observed in both the CO and Hα with similar velocities (Walter et al. 2002), making the association between the cold molecular gas and the warm ionized gas stronger (cf. Shopbell & Bland-Hawthorn 1998; Salak et al. 2013).

Such a direct relationship would be consistent with a picture where the emission we observe in the extended outflow of M82 is mainly from entrained and swept-up clouds that are predominately shock-heated by a high ram-pressure, hot plasma generated within the starburst region. If we envision the regions of warm H₂ emission as being the areas in the clouds that are shock-heated by the expanding plasma, we would expect these phases to co-exist with the warm ionized gas (Heckman et al. 1990) and the cold/warm molecular gas.

Assuming that the velocity of the warm H₂ is similar to the deprojected velocity of the cold molecular gas derived by Salak et al. (2013), we can calculate its total kinetic energy as |$E_K=0.5\times M_{{\rm H}_2}v_{{\rm outflow}}^2\sim 5.9\times 10^{53}$| erg for |$M_{{\rm H}_2}\sim 4.8\times 10^6 \,\mathrm{M}_{{\odot }}$| and v_outflow=160 km s⁻¹. Our estimates are considerably higher than those of Veilleux et al. (2009) because both their mass estimates are a factor of 500 less and they assumed a velocity of 100 km s⁻¹ for the warm molecular emission. However, they are smaller compared to those estimated for the cold molecular and warm ionized medium. Salak et al. (2013) estimate that the kinetic energy of the wind in the cold molecular gas is ≈1–4 × 10⁵⁶ erg. Walter et al. (2002) estimate a kinetic energy of ≈3 × 10⁵⁵ erg, but they do not map the CO distribution out to large distances from the disc compared to Salak et al. (2013).

4.3 H₂ supersonic turbulent heating and shock models

The hot plasma driving the outflow will likely transfer part of its ram and thermal pressure to swept-up and entrained dense molecular gas by driving shocks and by turbulent mixing between hot and cold gas (see Guillard et al. 2009). To model such shocks, we use a grid of stationary shock models computed with an updated version of the Lesaffre et al. (2013) code to derive the physical shock parameters that explain the H₂ line ratios. Our initial conditions are that the gas is diffuse, magnetized with a pre-shock magnetic field of intensity |$B[\mu G] = \sqrt{{\rm n}_{\rm H}[\rm cm^{-3}]}$|⁠, and irradiated with a radiation field of G₀ = 1 in Habing units. As the ‘Habing field’ would photodissociate the H₂, we assume some self-shielding. In the models, we integrate the extinction along the model from the pre-shock, and we take into account the shielding and self-shielding, which is parametrized in the photodissociation reaction rates of H₂ and CO (Lesaffre et al. 2013). At a density of 10⁴ cm⁻³, the ionization fraction in the pre-shock gas is 6 × 10⁻⁸. This matches approximately our estimate of the ionization fraction of 10⁻⁷ at 10 kpc from the X-ray source. We fit the H₂ line fluxes observed in the different regions with one shock or a combination of two shock velocities depending on how many H₂ lines are detected, for three different pre-shock densities: n_H = 10², 10³ and 10⁴ cm⁻³, [see Guillard et al. (2009) for details of the method]. The results of these fits are shown in Fig. 12, where we present the H₂ excitation diagrams for the best velocity combinations (lowest χ²), and for the intermediate pre-shock density of n_H = 10³ cm⁻³. We present in Table 4 the results of the shock models for each region and the integrated spectrum (labelled ‘Total’). The total gas mass derived from the shock models is |$M_{H_2}=1.7\times 10^7 \,\mathrm{M}_{{\odot }}$|⁠, which results in an estimate of the kinetic energy of E_K ∼ 2 × 10⁵⁴ erg. This is a factor of 3 higher than the estimates calculated in the temperature diagrams of Fig. 10. This difference arises due to the fact that we are integrating the post-shock temperature profile down to 50 K (where the H₂ mid-IR lines start to emit), and because LTE is not assumed in the shock models. Our estimates of the total warm gas mass agree with previous studies using Herschel Spire-FTS (Kamenetzky et al. 2012).

For all regions, the best-fitting shock velocities are in the range 4–35 km s⁻¹ (grid spans 3–40 km s⁻¹; Table 4). These velocities are much lower than the estimated bulk outflow velocity (Lehnert et al. 1999). The regions show a range of best-fitting models, but typically the best-fitting pre-shock densities are n_H = 100–1000 cm⁻³. Our modelling is not intended to be exhaustive, but to simply provide an estimate of the physical parameters in the warm molecular gas of the outflow.

At the range of preshock densities and magnetic field strengths considered, the C-shocks turn into J-shocks when the velocities are above ∼22 km s⁻¹. For some regions the H₂ line fluxes are reproduced by a combination of two C-shocks, in other regions the C-shocks (J-shocks) do not sufficiently populate the highest (lowest) energy levels. In those cases, we fit the H₂ line fluxes by a combination of C- and J-shocks. The shocks travel through the clouds in ∼10⁴–10⁵ yr (for 1 pc clouds), which is much shorter than the dynamical time of the outflow [10⁸ yr; Lehnert et al. (1999)]. Using the estimates of the warm molecular gas mass and the typical cooling time of the shocked H₂ gas (t_cool ∼ 10⁴ yr for n_H = 10³ cm⁻³), we can estimate the mass flow of shocked gas, M_warm/t_cool ∼ 7.5 × 10⁶ M_⊙/10⁴ yr = 750 M_⊙ yr⁻¹, assuming that all the warm H₂ is shock-heated. The dynamical time-scale of the H₂-rich outflow is ∼4.5 kpc/200 km s⁻¹ = 2.5 × 10⁷ yr implying a total mass of shocked gas ≈5 × 10⁹ M_⊙. This estimate is likely to be an upper limit since the H₂ gas may not extend over the whole outflow and because this assumes that the gas is only shocked once, which is highly unlikely. The high mass flow rate of shocked gas may require the reformation of H₂ in the post-shocked gas (see Guillard et al. 2009, 2010) to prolong the survival time of the H₂ gas in the outflow (e.g. Fragile et al. 2004). Clouds may also survive longer if the radiative cooling time is shorter than the time-scale necessary for the shocks to cross the cloud and the growth rate of Kelvin–Helmholtz instabilities (Cooper et al. 2008). Such processes would lower the total molecular mass required to explain the molecular emission from the outflow, although for the C-shocks considered here and even for many of the J-shocks as well, there is little destruction of H₂. At any rate, the total gas mass that is shock heated is similar to that of the cold molecular gas. Therefore, the reservoir for the shocked gas is likely the cold molecular gas as probed through the CO emission from the wind.

4.4 PAH destruction

Since we see that the PAHs grains are on average larger and more ionized in the wind, it is natural to ask whether or not some of the enhanced H₂/PAH could be done to a reduction in the PAH emission, rather than an enhancement of H₂. The decrease in PAH flux is probably just due to a decrease of the radiation field, but there is also a change in the grain properties associated with the wind. The gradient of decrease in H₂ emission is significantly lower than that of the PAH emission, suggesting an additional heating source for the H₂ (see Section 3.2).

One of the ways PAH emission can decrease relative to H₂ is PAH destruction. PAH grains can be destroyed either by hard radiation field, X-rays, or shocks. The EWs of the PAH features in starburst galaxies are observed to decrease sharply when [Ne iii]15.6 μm/[Ne ii]12.8 μm>2, where [Ne iii]/[Ne ii] traces the hardness of the radiation field. This has been interpreted as evidence of PAH destruction by intense radiation fields (Beirão et al. 2006; Wu et al. 2006). In M82, despite an increase of [Ne iii]15.6 μm/[Ne ii]12.8 μm in the outflow, in most of the outflow [Ne iii]15.6 μm/[Ne ii]12.8 μm remains less than 1. [Ne iii]15.6 μm/[Ne ii]12.8 μm∼1 indicates that there is hard (>40.96 eV) radiation present, but not enough to have an effect on the average PAH properties. Furthermore, only a small fraction of the gas in the M82 outflow is photoionized, and PAHs may survive in H i and H₂ gas. PAH destruction is also very sensitive to the clumping factor of the gas. If the gas is dense and clumpy with relatively small volume filling factor and large covering fraction, the PAH grains will be shielded from the hard photons in the radiation field. So PAH destruction by hard photons is not likely to be the cause of the increase of H₂/PAH in the wind.

PAH grains are rapidly destroyed in hot, X-ray emitting gas (Micelotta, Jones & Tielens 2010b). However, the presence of PAHs at large distances (10 kpc from the disc) in the outflow shows that they are protected in cooler phases of the gas, mostly in molecular gas. It is possible, however, that X-rays from the hot outflow fluid could destroy the PAH carriers, but the correlation between the X-ray emitting gas and the PAH maps, particularly in the Southern outflow, suggests that grain destruction by X-rays in the wind is not dominant.

Thermal collisions with very small and big grains can have a significant impact on the cooling of the hot gas ejected into the halo (Dwek, Rephaeli & Mather 1990). However, the impact of this cooling will generally be short lived since the cooling efficiency drops as the grain sputtering occurs in the hot gas (Guillard et al. 2009), unless the hot halo gas is constantly replenished with dusty gas. A physical dust cooling model is needed to fully assess whether thermal emission from collisionally heated dust in the hot plasma is a significant contribution to the gas cooling in the halo of M82, and whether this contribution changes as a function of the radius.

Another possible mechanism for PAH destruction is shocks. As seen in Fig. 12, the H₂ emission that we see can be fitted with shock models having speeds of up to 40 km s⁻¹. PAH destruction likely only becomes relevant for shock velocities >100 km s⁻¹ (Micelotta, Jones & Tielens 2010a). It is unlikely that the low-velocity shocks that heat the molecular gas and enhance the H2 emission are also destroying the PAH grains, although we clearly see a change in the grain population. Similarly, Contursi et al. (2013) find that shocks do not dominate the [C ii] and [O i] emission in the wind. Overall, the enhancement in H₂/PAH is likely caused by H₂ excitation by shocks, although we cannot rule out a relatively small effect due to PAH destruction.

In dense molecular clouds, where the PAH destruction by shocks is less efficient, cosmic rays are the most efficient source of destruction of PAH molecules, mostly because of their bombardment by energetic ions (Micelotta, Jones & Tielens 2011). This is likely to play an important role in the M82 wind, because the intense star formation in M82 is likely to enhance the cosmic ray intensity by a factor of ≈5 compared to the Milky Way (Robitaille & Whitney 2010). Theoretical studies suggest that CRs can destroy PAHs in dense molecular clouds, but in these environments the heating will also be very inefficient, causing a reduction of PAH emission. Micelotta et al. (2011) showed that only PAHs with a number of carbon atoms greater than ≈250 have a survival lifetime greater than the circulation time-scale of the outflow (≈2 × 10⁸ yr). Thus, there may be a size-selection effect of PAHs in the outflow, which is roughly consistent with the size inferred from the band ratios (we mostly see PAH 6.2/7.7 < 0.3, meaning an average size greater than 250 carbon atoms, see Fig. 6), and cosmic rays are likely to play a role.

The calculations of Micelotta et al. (2010a) show that any PAH shocked at >100 km s⁻¹ should be destroyed. To be re-formed in the post-shock gas, an efficient and fast re-formation process is needed. PAHs are a product of high-temperature chemistry involving abundant carbon bearing precursors such as CH₄ and C₂H₂. In the cooled post-shock gas these precursor species are likely not very abundant. In the cold neutral medium, the global PAH coagulation time-scale is of the order of 100 Myr, much longer than the outflow dynamical time (Seok, Hirashita & Asano 2014). AKARI results from Yamagishi et al. (2012) suggest that very small grains still exist in the halo. Therefore, it is reasonable to assume that PAH are produced in the disc, carried up in the outflow, and ‘protected’ from fast shocks in the denser material. If molecular clouds are formed in the post-shock gas in the outflow, then there may be PAH coagulation over time-scales comparable to the outflow time (coagulation time-scales are of the order of a few Myr in molecular clouds), but it should be very small fraction of the volume, contrary to what is observed.

In conclusion, the enhancement of the H₂/PAH is likely caused by H₂ shock excitation, although we cannot rule out a relatively small effect due to PAH destruction.

5 CONCLUSION

• We have mapped the northern and southern superwind in the nearby starburst galaxy M82 using the Spitzer-IRS. We find strong H₂, PAH and Ne line emission that extends over nearly 5 arcmin above and below the plane of the galaxy, and have measured the variation of the H₂/PAH, [Ne iii]/[Ne ii] and PAH feature ratios throughout the wind.

• We estimate a total warm (100–500 K) molecular gas mass to be between ∼5 and 17 × 10⁶ M_⊙ and the total kinetic energy of the warm molecular outflow to be between ≈6 and 20 × 10⁵³ erg assuming an outflow velocity v_outflow = 160 km s⁻¹. The warm molecular gas provides about 0.5–5 per cent of the kinetic energy – comparable to that observed in the atomic gas.

• There are clear variations in the PAH band ratios throughout the wind. The 6.2/7.7 and 11.3/7.7 μm PAH band flux ratios, when compared to models, suggest a higher fraction of larger and/or more ionized PAHs in the wind compared to the starburst disc, especially in the northern part of the wind.

• The warm molecular gas to PAH ratio (H₂/PAH) is enhanced in the outflow by factors of 10–100 as compared to the starburst disc. This enhancement in the H₂/PAH ratio does not seem to follow the trends with the ionization of the atomic gas (as measured with the [Ne iii]/[Ne ii] line flux ratio) seen in the starburst. Our results indicate that UV and X-ray emission are probably not the dominant excitation mechanisms of H₂. Models suggest the observed values could be reproduced with slow C- and J-shocks (v < 40 km s⁻¹) driven into dense molecular gas entrained in the outflow. Using a dynamical time-scale of the H₂-rich outflow of about 10⁷ yr results in an upper limit of 5 × 10⁹ M_⊙ of molecular gas heated and shocked over the lifetime of the wind. Given the rate at which the gas must be shocked, the reservoir for this shocked gas is the cold molecular gas as observed in CO.

• The emission from all of these phases is likely due to clouds entrained in and accelerated by the flow. They are being shock heated by the hot plasma with high ram pressure generated by the collective thermalization of stellar winds and supernova in the starburst region of M82. These clouds are multiphase and the variation in properties are likely due to a range of shock velocities which depend on the density of the pre-shocked gas which depend on the density of the pre-shocked gas.

We would like to thank an anonymous referee for a number of helpful comments which helped improve the quality of this paper. This work is based on observations made with the Spitzer Space Telescope, which was operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. Support for this work was provided by NASA through an award issued by JPL/Caltech. This work is dedicated to the memory of our friend and colleague, CE.

REFERENCES