THE MEGASECOND CHANDRA X-RAY VISIONARY PROJECT OBSERVATION OF NGC 3115: WITNESSING THE FLOW OF HOT GAS WITHIN THE BONDI RADIUS

NGC 3115 (Figure 1) was observed eight times with the Chandra in 2012 between January and April (ObsIDs 13817, 13819, 13820, 13821, 13822, 14383, 14384, and 14419) for a total of 998 ks. Unless otherwise specified, we only included all these 2012 observations in the data analysis but do not include the 155 ks observations taken in 2001 (ObsID 2040) and 2010 (ObsIDs 11268 and 12095) because the effective area has changed dramatically even since 2010 due to the increasing Advanced CCD Imaging Spectrometer (ACIS) contamination. We note that including the 2001 and 2010 observations generally introduced a ≲10%–20% systematic error in the spectral analysis but does not improve the statistical uncertainties. In all the observations NGC 3115 was placed near the ACIS-S aimpoint. All the data were reprocessed using the chandra_repro script of CIAO 4.4 and CALDB 4.4.10. The default sub-pixel event-repositioning algorithm "EDSER" was used. After we removed the flares using the CIAO deflare script, the cleaned exposure time was 972 ks.

Zoom In Zoom Out Reset image size

To improve the astrometry between different observations required for the high spatial resolution analysis, we have performed relative astrometry correction for each observation. We first created a sub-pixel resolution image in 0.3–6.0 keV with a binning size of 0123 (0.25 ccd pixel) for a ∼4 × 4 arcmin² region around NGC 3115 for each observation. We then used the CIAO wavdetect script to create an initial source list for each image. This source list was only used for astrometry correction. The longest observation (ObsID 13820) was used as the reference for the relative astrometry corrections. The CIAO reproject_aspect script was then used to create new aspect solutions for all the other observations. All of the data were reprocessed again using chandra_repro with the new aspect solutions to complete the astrometry corrections.

We extracted a local background from a 70''–90'' annular region far enough from the center of NGC 3115 so that the source-removed surface brightness of the X-ray emission is basically flat. Background contributes negligibly to the inner ∼10'' but becomes significant only in the outermost regions (Appendix A). Changing the background level by ±10% introduces systematic uncertainties that are much smaller than the statistical uncertainties within 20''. It only introduces a systematic uncertainty that is larger than the statistical uncertainty in gas temperature beyond 20'', although the gas normalization is still hardly affected.

To analyze the diffuse X-ray emission of the hot gas it is necessary to remove contaminating point sources (or compact structures). Point sources were detected with CIAO wavdetect. To detect as many point sources as possible, we used all the observations except the data taken in 2001 because its optical axis position is significantly different (>1') from the other observations. We created images with 1 ccd pixel binning size in four energy bands (0.3–1.0, 1.0–2.0, 2.0–6.0, and 0.3–6.0 keV) and combined images according to energy bands in these observations. The source regions in different energy bands were then visually inspected and combined. We refined the region sizes of the point sources (or compact structures) detected within 4'' by using sub-pixel images with 0.125 pixel binning size and ran CIAO wavdetect again. With these sub-pixel images, a few more weak sources and also some elongated structures were identified within 3'' (Figure 1). Unless otherwise specified, we have removed all these structures except for the source detected at the galaxy center in our nominal data analysis. Including or removing these structures within 3'' gives essentially the same results for the gas component.

As mentioned above, the central peak was detected with wavdetect. Our analysis shows no strong evidence of a point source and the central peak is clearly extended (Section 3). We have determined the upper limit of the potential point source and found that hot gas measurements are not affected by this potential weak AGN (Section 3 and Appendix B). Therefore, we did not remove the central region and we ignored any potential AGN contamination in our analysis.

We extracted spectra in circular annuli centered on the central peak of the extended X-ray emission (0.3–6.0 keV), which is assumed to be the center of the flow (the supermassive black hole). This peak is within ≲ 005 of the soft (0.3–2.0 keV) emission peak. It is separated by 015 from the optical peak we measured using the archival Hubble Space Telescope (HST) data (below), consistent with the position uncertainty. The diffuse gas distribution is assumed to be spherically symmetric which is justified in Section 4.2. All the spectra were analyzed using the X-ray Spectral Fitting Package⁶ (XSPEC).

The unresolved X-ray emission is mainly contributed by unresolved low-mass X-ray binaries (LMXBs), stellar emission from cataclysmic variables and coronally active binaries (CV/ABs), and the diffuse hot gas component that we are interested in. After the CV/AB component is spectrally subtracted in a statistical manner (described below), the very soft emission from the gas and the very hard emission from the unresolved LMXBs can be reliably separated through spectral fitting. The combined spectra of resolved low-L_X (<10³⁷ erg s⁻¹) LMXBs in the bulge of M31 are very similar to more luminous LMXBs (Irwin et al. 2003). Therefore, we can assume the unresolved LMXB emission to be spectrally modeled as the brighter resolved sources. With our deep Megasecond observation, many more LMXBs were detected than in W11, and the unresolved LMXBs are no longer the dominant component at ≳ 2'' in the 0.5–1.0 keV band (Figure 5 below). We modeled the LMXB component as a power law model and fixed the power-index to Γ_LMXB = 1.6 which is consistent with the value of $1.61^{+0.02}_{-0.02}$ measured from the combined spectrum of all the resolved point sources within D₂₅ of NGC 3115 and the value (1.6) of the summed emission from many resolved X-ray binaries in nearby galaxies (Irwin et al. 2003). Using Γ_LMXB = 1.4 or 1.8 gives essentially the same results for the analysis of the hot gas (Appendix A).

The faint and soft sources similar to those of the Galactic Ridge emission (CV/AB) contribute appreciably to the X-ray flux beyond the Bondi region. Hence, including this component is essential in the analysis. In W11, we used the Two Micron All Sky Survey (2MASS) K-band image to estimate the CV/AB contribution and assumed the X-ray flux of the CV/AB component scales linearly with the K-band luminosity. Because the typical spatial resolution of 2MASS is about 2''–3''  which is poorer than the Chandra resolution, the CV/AB in the central regions can be underestimated.

In this paper, we use a higher resolution HST WFPC2 I-band (F814W filter) image to estimate the CV/AB contribution. We assume the intrinsic K-band surface brightness profile follows the I-band surface brightness profile. The archival HST I-band image has a point-spread function (PSF) FWHM and a pixel size of ∼01 which is much smaller than the 2MASS PSF (2''–3'') and is close to the Chandra PSF (∼05 near the aimpoint). Therefore, we can use the HST I-band surface brightness profile to estimate the intrinsic K-band surface brightness profile. In practice, we first smoothed the HST I-band image with a Gaussian kernel. We then re-scaled the HST I-band image to match the 2MASS flux unit within a radius of 10'' centered at the surface brightness peak. After that, we generated a surface brightness profile with a radial binning size of 1'' within 10'' for the HST I-band image and compared it to the 2MASS surface brightness profile. We found that with a Gaussian kernel of 25 FWHM, the smoothed HST I-band surface brightness profile matches the 2MASS profile moderately well. The deviations between the two profiles are at most 17% at all radii which are smaller than the conservative 50% uncertainty of the CV/AB contribution we considered in this paper (see below). The unsmoothed HST I-band image was then scaled with the same factor as the smoothed HST image and was used to estimate the intrinsic K-band surface brightness. With this correction, the flux within the central 1'' is higher than the 2MASS estimate by about a factor of two.

To model the X-ray spectrum of the CV/AB component, we fitted the unresolved X-ray emission of the dwarf elliptical galaxy M32, which is believed to be hot gas-free (Revnivtsev et al. 2007; Boroson et al. 2011). Using an absorbed thermal + power law [PHABS*(APEC+POWERLAW)] model fitted to archival Chandra data, the best-fit temperature was $T_{\rm CV/AB}=0.76^{+0.07}_{-0.06}$ keV and the power law index was $\Gamma _{\rm CV/AB} = 1.92^{+0.07}_{-0.11}$. The CV/AB normalizations of each annulus were determined by the L_X–L_K scaling relation derived from M32. We investigated the systematic uncertainties by varying the CV/AB normalizations by ±50%, comparable to the galaxy by galaxy variation of this component (Appendix A). This only changed the measured gas properties slightly within ∼10''. All of the systematic uncertainties are smaller than the statistical uncertainties.

The eight observations in 2012 were observed in three different periods (two observations between January 18–January 23, three observations between January 26–February 5, and three observations between April 4–April 7). For each period, the pointings and roll angles of the observations are nearly identical. Therefore, we merged the spectra for observations taken in each period using the CIAO specextract script. The three merged spectra of each extraction region in the 0.5–6.0 keV energy range were fitted jointly to the three component absorbed (PHABS) model—a thermal (APEC) model for the gas, a power law (POWERLAW) with a slope of 1.6 for the unresolved LMXBs, and a combination of thermal + power law (APEC + POWERLAW) model for the CV/ABs. We fixed the absorption at the Galactic value of N_H = 4.32 × 10²⁰ cm⁻² (Dickey & Lockman 1990). The systematic uncertainty introduced by N_H is not significant (Appendix A). We fitted the temperature of the thermal gas component. The normalizations of the LMXB component were allowed to vary in the three merged spectra to account for possible variability. For the thermal gas component, when the three normalizations are untied in the joint fitting, they are within the uncertainties of each other. Moreover, the ccd responses did not change much during the observations and the extended hot gas should not be time varying. Therefore, we tied the hot gas normalizations in the fitting. The metallicity was fixed to the solar value using the wilm abundance table (Wilms et al. 2000) and the systematic uncertainties of this assumption are discussed in Appendix A. In brief, thawing the metallicity only introduces a small systematic bias in hot gas temperature compared to the statistical uncertainty and increases the gas normalizations by a factor of four to five, and the derived gas density only increases by a factor of two without affecting the density slope (Appendix A). None of our conclusions are sensitive to the adopted metallicity. Unless otherwise specified, all the spectra were fitted using the c-statistic.

4.1.1. Single Temperature Model

We model the projected spectra of the diffuse gas component with a single temperature optically thin thermal plasma model (APEC) described in Section 2.3. The projected temperature profile is shown in the upper panel of Figure 3. It is clear that the projected temperature is around 0.3 keV in the outer region and increases sharply to about 0.7 keV within about 5''. Within the inner 2'' or 3'', the new data now strongly suggest that the projected temperature drops at the center, indicating that there is significant soft emission near the center. Previous analysis of the moderate length 2001+2010 Chandra data by W11 suggested an increase in best-fit projected temperature toward the center, although the uncertainty was too large to be conclusive. Such a simple interpretation of a monotonic increase in projected temperature is no longer valid. We also noticed that a single temperature fit to the projected temperature within the inner ∼3'' is no longer adequate. In Section 4.1.2 below, we present two-temperature fitting results to the data.

Zoom In Zoom Out Reset image size

4.1.2. Two-temperature Model

Single temperature fitting to the projected spectra within the inner 3'' suggests a central drop in temperature. However, a single temperature model may not be sufficient to characterize the projected spectra; at the very least the spectra should consist of gas at different temperatures due to projection effects of gas at larger radii. Motivated by the expected central rise in temperature for RIAFs accreting in hot modes (e.g., Narayan & Yi 1994; Fabian & Rees 1995; Brighenti & Mathews 1999; Narayan & McClintock 2008; Guo & Mathews 2013) and evidence of central temperature peaks ≲ 300 pc from the galactic nuclei in a few early type galaxies (Pellegrini et al. 2003, 2012; Humphrey et al. 2008), we searched for possible evidence of a hot thermal component potentially hidden in a multi-temperature structure within the central 3''. We first examined the spectra with different binning sizes including or removing structures within 3''. We also combined all the different observations taken in 2012 to create a single spectrum in each extraction region, so that the combined spectra have enough counts for visual inspection. This also allows us to group the spectra with a minimum of 25 counts per spectral bin to use chi-squared or F-statistics. We found that while a single temperature model generally gives a good enough fit "globally" ($\chi ^2_{\nu } \approx 1$) for most spectra, there is a notable systematic excess of emission at about 1 keV. An example of a spectrum in an annular region of 1''–3'' is shown in Figure 4. The best-fit temperature of a single temperature model is $0.37^{+0.11}_{-0.06}$ keV with χ² = 75.9 and 80 degrees of freedom. When we added one more thermal component (an APEC model with the same abundance and redshift as the first thermal component), this two-temperature model gave best-fit temperatures of $1.23^{+0.25}_{-0.21}$ and $0.29^{+0.05}_{-0.05}$ keV with χ² = 53.6 and 78 degrees of freedom. The lower panel of Figure 3 shows the spectrum of the two-temperature model. A simple F-test with the F-statistic of 16.3 and probability of 1.26 × 10⁻⁶ strongly suggests that the second component is needed. A more formal test for an additional component was also performed by simulating 1000 spectra with the single temperature model and then comparing the likelihood ratio of the single temperature model with respect to the two-temperature model (likelihood ratio test in XSPEC). We found that all 1000 of the simulated likelihood ratios are smaller than the observed ratio, strongly suggesting that the two-temperature model is preferred. We performed a similar likelihood ratio test to determine whether the extra component is a narrow single line emission or a broader thermal component by simulating 1000 spectra with an extra Gaussian model in XSPEC. The line width was fixed to zero. We found that 99.2% of the simulated likelihood ratios of the single temperature + Gaussian model with respect to the two-temperature model are smaller than the observed ratio, again, strongly favoring the two-temperature model. We conclude that there is evidence of a hotter thermal component with temperature ≳ 1 keV within the central 3'' region.

Zoom In Zoom Out Reset image size

It is possible that there is a wider distribution of temperature structure along the line of sight. Unfortunately, the statistics of our data do not allow us to test beyond a two-temperature model. It is also known that a multi-temperature structure is difficult to quantify from the data (e.g., Kaastra et al. 2008; Gayley 2013). Here, we simply characterize the thermal component within 3'' with a two-temperature model. This two-temperature model at least may be able to characterize the rough lower and upper limits of the temperature distribution. Temperature profiles of the two-temperature model are shown in the lower panel of Figure 3. The fittings were performed by joint-fitting different observations as described in Section 2.3 rather than fitting the combined spectra. Different colors represent different spatial binning. It is interesting that the hotter temperature rises all the way toward the center, consistent with most hot accretion models. Considering the hotter component within 3'' and the single temperature profile beyond that, fitting the projected temperature profile to a power law gives $T \propto r^{-[0.44^{+0.29}_{-0.33}]}$ (90% confidence) for r < 5'' and $T \propto r^{-[0.34^{+0.25}_{-0.25}]}$ for 5'' < r < 40''. Ignoring the central 1'' region gives essentially the same results. It is also interesting that the lower temperature of the two-temperature model is consistent with the low temperature of ∼0.3 keV outside the Bondi radius.

Since it is quite certain that there are at least two (or multiple) temperature structures within about 3'', we also tested whether there is any evidence of significant two-temperature structure beyond 3''. Between 3'' and 5'', the best-fit temperatures and normalizations of the hotter component are consistent with a single temperature fit; the best-fit lower temperature is consistent with zero, which is unphysically low and the flux (0.5–2.0 keV) of the soft thermal component is less than 10% of the hot component. Beyond 5'', either the higher and lower temperatures are consistent with each other within 90% confidence, or the best-fit temperatures are more sensitive to systematic uncertainties. The best-fit temperatures (∼0.3 keV) and normalizations of the cooler component are consistent with single temperature fits. The cooler component dominates over the hot component by a factor of about two to four in normalization. Performing additional tests by tying the higher (and lower) temperature of a few radial bins together does not change the conclusion. We conclude that a single temperature model is sufficient to approximate the thermal structure of the hot gas beyond 3''.

The surface brightness profiles for the hot gas, the CV/AB, and the unresolved LMXB components are shown in Figure 5. Unlike the surface brightness profiles in W11 which assume all the emission above 2 keV is contributed by the LMXBs, in this paper, the surface brightness of the hot gas and unresolved LMXB were calculated from the best-fit models spectroscopically. This self-consistently takes into account the uncertainties in the hot gas and LMXB contributions. We use a single temperature model for the hot gas in this plot. In the central region where the gas is multi-temperature, this single temperature model, which has a low temperature of ∼0.3 keV near the center (Section 4.1.1), should be able to take into account most of the gas emission of the hot gas in the 0.5–1 keV band. Using a two-temperature model gives essentially the same result but gives larger uncertainties.

Zoom In Zoom Out Reset image size

The hot gas is robustly detected within ∼20''–40''. Hot gas is the dominant component in this energy band out to ∼10''–20''. The CV/AB component is not significant within a few arcsec but becomes significant in the outermost regions. However, varying the CV/AB contributions by ±50% does not change any of our conclusions qualitatively at all radii and does not significantly change any result quantitatively within a few arcsec (Appendix A).

The distribution of optical light of NGC 3115 is highly elliptical. We have tested whether the thermal X-ray emission deviates from azimuthal symmetry by extracting surface brightness profiles of the hot gas in four 90 degree sectors in NW, NE, SE, and SW directions. We conclude that there is no evidence of azimuthal variation for the single temperature model and the hot component of the two-temperature model, and therefore, spherical approximation is adequate for our analysis. However, there is some weak evidence that the cooler component of the two-temperature model is distributed more along the major axis of the galaxy within 3'' and this implication is discussed in Section 5.4. Nevertheless, we also present systematic tests by assuming the gas is distributed elliptically as the optical light and as a thick disk in Section 5.1.

Figure 6 shows the XSPEC APEC normalization per unit surface area, which is proportional to the emission measure of $\int n_e^2 dl$, where n_e is the electron density and l is the column length along the line of sight. The single temperature model is shown in black in the upper panel and thick gray in the middle and lower panels. The hotter (cooler) component of the two-temperature model within 3'' is shown in color in the middle (lower) panel, with different colors representing different spatial binning. In general, the emission measure of the hotter component is lower than the single temperature model but barely consistent within the error bars while the cooler component is comparable to the single temperature model. This indicates that the hotter gas density may be slightly lower and the cooler gas density may be similar to that determined by a single temperature model (as shown below), but the density profiles should not be too sensitive to these two models (since $n_e \propto \sqrt{\rm emission\,\, measure}$). In the central 15, the uncertainties of the hot component appear to be significantly larger than that of the single temperature model. The large error bar is due to the poor temperature constraint in that region, with a higher temperature upper limit so that the gas normalization is degenerate with the hard emission from, e.g., LMXBs. We also noted that there is a weak AGN in the central 1'' which can contribute at most up to 30% of the total emission and can increase the hot component uncertainties. Note that this upper limit is very conservative (see Appendix B).

Zoom In Zoom Out Reset image size

Within 3'', the soft normalization per unit area decreases with radius (lower panel in Figure 6), suggesting that the cooler component should not be projected emission from a much larger relatively uniform background structure. The centrally peaked soft normalization per unit area suggests that the cooler component should be located physically inside about 150 pc (3'') from the galaxy center (see Section 5 for more detailed discussions).

With the XSPEC APEC normalizations (or emission measure) in each annulus, we can deproject the density profile using the onion peeling method as was done in W11 and also outlined in detail in Kriss et al. (1983) or Wong et al. (2008). In brief, this technique calculates the emission measure of each spherical shell starting from the outermost annulus toward the center, and the emission measure of each subsequent shell is calculated by subtracting the projected emission measure from the outer shells. Unlike W11 who deprojected the density from the surface brightness and assumed a certain spectral shape (or gas temperature), in this work, we directly deprojected the density profile from the spatially resolved emission measure fitted from spectra, and therefore, the uncertainties of the spectral shape (or temperature) have been partially taken into account. Note that the deprojected density determined here is effectively the root-mean-squared of the density. If the filling factor is less than one or if the gas is not homogeneous, our deprojected density is overestimated (see also Shcherbakov et al. 2013). Unfortunately, X-ray observation alone cannot determine the filling factor or clumpiness. Theoretical models may provide constraints to these factors.

The deprojected electron density profiles of the single and two-temperature models are shown in Figure 7. The errors were estimated by running 10⁶ Monte Carlo simulations. The density profile of the single temperature model (upper panel of Figure 7) is easier to interpret by assuming a single-phase plasma. However, the two-temperature model is more difficult to interpret. It is possible that the two components are in two (or more) distinct phases or the particles can be distributed in a broader than Maxwellian distribution. In any case, the distribution of the two-temperature model cannot be constrained without additional assumptions (e.g., filling factor, gas distribution). Motivated by hints that the colder component may be distributed as a (∼3'') disk along the major axis of the galaxy and also by the fact that the hotter component appears to be more spherical (see Sections 4.2 and 5.4), we assume a simple model that the hotter component characterizes the spherically distributed hot gas in projection within 3'' (i.e., assuming a filling factor of 1 for the hotter component). We assume the cold component is concentrated in a small disk-like region that can be ignored when doing the spherical deprojection analysis. The origin of such a cold component is discussed in detail in Section 5. With the two-temperature model within 3'', the emission measure (normalization) cannot be constrained well enough if we thaw the temperatures. In particular, using a narrow spatial binning size of 1'' gives too large statistical uncertainties and also the systematic uncertainties in the central 1'' can be larger. To improve the constraints, we used a larger spatial bin of 15 for deprojection. We fixed all the higher and lower temperatures to their best-fit values and assessed the uncertainties of the normalizations. The deprojected density of the hotter component is shown in the middle panel of Figure 7. As expected, it is slightly lower than the single temperature model. We noted that if we assume the cold component as the spherically distributed gas and ignore the hot component (lower panel of Figure 7), the density profile is closer to the single temperature model. Note that the density profiles under these three different assumptions are remarkably similar, suggesting that these models are measuring similar emission measure that is more sensitive to density than temperature. Note also that if the filling factor is less than one or if the gas is clumpy, the density we measured is biased high.

Zoom In Zoom Out Reset image size

Fitting the density profile of the hotter component of the two-temperature model within 5'' to a power law gives $\rho \propto r^{-[0.89^{+0.35}_{-0.45}]}$ (90% confidence; note that W11 present 1σ confidence). The less physically motivated density profile of the cooler component of the two-temperature models gives $\rho \propto r^{-[1.08^{+0.31}_{-0.24}]}$ in 0''–5''. The single temperature model gives a power law index of $1.05^{+0.25}_{-0.25}$ in 0''–5'', $0.90^{+0.24}_{-0.30}$ in 0''–4'', and $0.62^{+0.26}_{-0.38}$ in 0''–3''. The density profile becomes steeper in the 5''–40'' outer region, with a power law index of $1.34^{+0.20}_{-0.25}$.

Because the characteristic temperature (∼0.3 keV) of the softer emission is very similar to the temperature of the outer region beyond ∼5'', it is possible that the softer emission comes from projected gas of the outer regions. However, the surface brightness (or emission measure per unit area) is quite steep in the central 3'' (lower panel in Figure 6), which does not appear to come from projected gas of a larger outer region. It would be ideal if we could fit a projection model (such as the XSPEC projct, although it has some limitation in the assumed geometry) to test whether projection can account for all the soft emission. Unfortunately, the statistics of the data do not allow us to perform such a test with high confidence. Instead, we performed three conservative tests to quantify the allowed projected soft emission.

We first tested whether the soft emission within 3'' can be projected from a spherical distribution of ∼0.3 keV gas beyond 5''. To maximize the projection effect, we assumed that all gas beyond 5'' has a temperature of 0.3 keV. This may overestimate the cool gas contribution as the hot gas temperature between ∼5''–10'' is slightly hotter and most of the projection should come from this region. We then fitted the gas normalizations in each annulus beyond 5'' with the gas temperature fixed to 0.3 keV. By using the onion peeling method (Section 4.3; Kriss et al. 1983; Wong et al. 2008), the projected gas (flux) contribution to the inner 3'' region can be calculated. In the 1''–3'' region, we fitted a two-temperature model with the lower temperature fixed to 0.3 keV to estimate the total flux of the cooler component. The inner 1'' is ignored due to the potential contamination from a weak AGN. We found that projected gas can only account for 11% of the soft emission in the 1''–3'' annulus. Even if the soft emission is at its lower limit of the 90% confidence interval (95% one-sided limit), projected gas can only account for 16% of the soft emission. The statistical uncertainty of the projected gas normalization is of the order of 20%–30%, and therefore, the uncertainty in projected gas cannot account for the difference.

We then constructed an oblate spheroid model of the gas halo with constant ellipticity which roughly follows the optical light (minor radius/major radius = 0.6; Kormendy & Richstone 1992). Here, we extracted spectra in elliptical annular regions with the radial binning sizes along the major axis equal to our circular annulus sizes. By doing a similar analysis as the spherical model above, we found that projected gas can still at most account for about 22% of the soft emission in the 1''–3'' annulus.

Finally, we assume a thick circular disk of uniform gas with thickness of 6'' and an outer radius of 40'' aligned along the optical major axis. The rotation axis is assumed to be parallel to the plane of the sky along the optical minor axis. We extracted a spectrum from a 6 × 80 arcsec rectangular region aligned along the major axis, with a 3'' circular region at the center and point sources excluded. Again, we fitted the APEC normalization of a 0.3 keV gas in this region. The projected gas within the 3'' region is proportional to the projected volume. We found that projected gas can only account for 9%–13% of the soft emission within the central 1''–3'' annulus. This disk model accounts for a smaller amount of projected gas compared to the ellipsoid model. This may be due to the very bright optical disk in the disk region that overestimates the CV/AB contribution in X-ray and hence underestimates the projected gas. Another reason may be that the gas emission is not disk-like as assumed.

In summary, most of the soft emission within the central 3'' cannot be explained by projected emission from a spherical, an oblate spheroid, or a thick disk distribution of cooler gas. Projection may work, e.g., if the outer cooler gas is preferentially distributed toward the line-of-sight of the supermassive black hole, which seems unlikely.

It has been suggested that tidally stripped cores of giant stars can have very soft spectra (E ≲ keV) with relatively large luminosities (>100 L_☉) which can last for 10³–10⁶ yr around supermassive black holes (e.g., Di Stefano et al. 2001; Davies & King 2005) and could conceivably account for the soft emission inside 3''. The tidal radius of a billion solar mass black hole is R_t ≈ 1.2 R_S(M_*/M_☉)^−1/3(M_BH/10⁹ M_☉)^−2/3(R_*/5 R_☉), where R_S is the Schwarzschild radius, and M_* and R_* are the stellar mass and radius, respectively (MacLeod et al. 2012). Therefore, most main-sequence stars with R_* ≲ 5 R_☉ are directly swallowed by the billion solar mass black hole, but the envelopes of giant stars can be tidally stripped.

First, for the nearest two point-like sources or extended regions located near the major axis at 1'', spectral analysis suggests that their spectra are perfectly consistent with a typical LMXB spectrum with a power-law index of 1.65 ± 0.13. These weak sources do not affect the spectral analysis of the hot gas and do not provide excess soft emission in the central 3'' region. We have removed all the detected point sources within the 1''–3'' annular region, and there is no indication of any other point source. If the soft emission comes from these soft stripped cores, the luminosity of each source has to be lower than ∼10³⁶ erg s⁻¹ to remain undetected.

Second, the rate of all stars passing through the corresponding tidal radius (including those swallowed and disrupted) of NGC 3115 has been estimated to be about ${\dot{N}} \sim 5 \times 10^{-5}$ yr⁻¹ (Wang & Merritt 2004), which is roughly consistent with other estimations of giant galaxies with supermassive black holes (e.g., Magorrian & Tremaine 1999; Syer & Ulmer 1999). The fraction of giant stars (R_* > 5–10 R_☉) passing through the tidal radius is about 5%–10% (e.g., MacLeod et al. 2012). Taking the upper limit of 10% gives the tidal stripping rate of giant stars to be ${\dot{N}_G} \sim 5 \times 10^{-6}$ yr⁻¹. Even if a stripped core can maintain its luminosity for as long as 10⁶ yr, there are only about five stripped cores expected to be luminous at one time around the center of NGC 3115. The total unabsorbed luminosity of the soft component of the two-temperature model in the 1''–3'' annular region is L_0.5-2 keV = 2 ± 0.4 × 10³⁷ erg s⁻¹, implying that at least 20 low-luminosity stripped cores are required to account for the soft emission. This is a factor of >4 greater than the expected number of stripped cores at any one time.

Third, from the spatial scale of ∼100 pc of the soft emission and a typical velocity dispersion of 250 km s⁻¹, these cores would have traveled for 0.4 Myr. If most of the soft emission comes from these stripped cores, the average lifetime of the emission has to be close to the upper limit of the expectation.

In summary, to account for all the soft X-rays with stripped giant cores, each core needs to be less luminous than ∼10³⁶ erg s⁻¹, the tidal stripping rate of giant stars needs to be higher than the predicted value of a few × 10⁻⁶ yr⁻¹, and the average lifetime of the emission should be longer than a fraction of a Myr. These set very tight and challenging constraints on the properties of the stripped cores to meet in order to explain most of the soft emission seen at the galaxy center.

Given the difficulties of the projected gas and stripped cores scenarios, it is likely that most of the soft emission is physically located within the Bondi region. In fact, it is expected that the interstellar medium (ISM) can be in a multi-temperature phase⁷ (e.g., McKee & Ostriker 1977; Muno et al. 2004; Randall et al. 2006). Simulations have shown that hot gas accreted toward supermassive black holes can be chaotic and has a wide range of temperature within the Bondi radius (Barai et al. 2012; Gaspari et al. 2013; Das & Sharma 2013). Numerical simulations also suggest that thermal instabilities of non-rotating cooling gas occurs when t_cool/t_ff ≲ 10, where t_cool and t_ff are the cooling time and free fall timescales (Sharma et al. 2012; Gaspari et al. 2013). For NGC 3115, t_cool/t_ff ≈ 100 (Shcherbakov et al. 2013), so cooling may not be important to induce thermal instability. However, it may be that there are some regions in NGC 3115 where t_cool/t_ff ≲ 10 locally and therefore cooling can become important. This cooling gas may also be cooling out of the X-ray emitting hot phase (T ≳ 0.1 keV) and therefore only the hotter phase with longer t_cool is detected (see Figure 15 in Gaspari et al. 2013). It may also be that the gas is not free falling (accretion rate greatly suppressed by, e.g., rotation), and therefore, the relevant timescale is the accretion time rather than the free fall time. The accretion rate of any cooling gas out of the X-ray band should, however, not be too high to trigger a powerful AGN at the moment. From a theoretical point of view, such a multi-temperature phase, however, is less likely to be clumpy as clumpiness or fragmentation in accretion flow is more likely to occur with higher accretion rate ($\eta \dot{M} c^2/L_{\rm Edd} \gtrsim 0.02$, where η is the radiative efficiency; Wang et al. 2012) or larger L_{X, AGN}/L_Edd ∼ 0.01, where L_{X, AGN} is the AGN luminosity (Barai et al. 2012).

Hot gas in an early type galaxy with significant stellar rotational velocity is likely to be rotating to some degree, in particular, if a significant faction of the hot gas comes from stellar mass loss. For NGC 3115, the total hot gas within 10 R_B is about 5 × 10⁶ M_☉. Assuming the specific stellar mass lost rate of 1.5 × 10⁻¹² yr⁻¹ estimated by Mathews (1989) and a total stellar mass of 5 × 10¹⁰ M_☉ within 10 R_B (Kormendy & Richstone 1992) implies that only 70 Myr is needed to build up the hot gas from stellar mass loss, which is much shorter than the expected age of an early type galaxy (on the order of 10 Gyr). Therefore, we expect the angular momentum of the hot gas in NGC 3115 to be comparable to the stellar component. In fact, X-ray observations and numerical simulations suggest that hot gas can be rotating collectively with some rapidly rotating early type galaxies (e.g., NGC 4649; Brighenti et al. 2009). The implications of rotation of hot gas are discussed briefly in Section 7 below.

As mentioned in Section 4.2, there is some weak evidence that the cooler gas component of the two-temperature model is preferentially located on the major axis while the hot component is more spherically symmetric within the 3'' region. Interestingly, there is also a very small and distinct optical thin disk with a radius of about 3'' along the major axis of the galaxy as shown with HST (Figure 4 in Kormendy et al. 1996 and Figure A8 in Ledo et al. 2010). It is possible that some of the gas can circularize and have enough time to cool toward a small disk region. If this is the case, the cooler gas can be dynamically uncoupled (or weakly coupled) with the hotter gas halo/flow. Such a geometry can also allow us to deproject the density profile of the more spherical gas by considering only the hot gas component within 3''. A rigorous test of this small cooler disk model is beyond the scope of this paper.

Since the dynamical timescale (sound crossing time) is much shorter than the heating, cooling, or conduction timescales near the Bondi radius (Shcherbakov et al. 2013), in the absence of a (supermassive) black hole, the hot gas should be in hydrostatic equilibrium (HSE) with the galactic potential. We tested whether the measured gas density profile is consistent with hot gas in HSE with the galactic potential without a black hole. We model the hydrostatic gas density profile with the total stellar mass profile described in Shcherbakov et al. (2013). The gas mass is neglected as it is much smaller than the stellar mass, and dark matter is also neglected near the Bondi scale in which we are interested. Figure 8 shows the density predicted by the adiabatic model (thick dashed green) and the isothermal model (dot-dot-dashed red). Note that the measured temperature profile of the cooler component of the two-temperature model is perhaps fairly isothermal and the entropy of the hotter component is close to adiabatic⁸ justify the comparison to isothermal and adiabatic HSE models. The HSE profiles, as well as the measured data points, are normalized at 5''. A single temperature model of the hot gas, which roughly corresponds to the cooler thermal component of the two-temperature model, is not consistent with the two HSE models, even accounting for the possible normalization (density) uncertainty at 5'' (thin dashed green). The hot component of the two-temperature model is not consistent with the isothermal HSE model, as expected from the rising temperature toward the center. It is also not very consistent with the adiabatic HSE model. If we account for the normalization (density) uncertainty by shifting up the adiabatic HSE model by the density uncertainty at 5'' (thin dashed green), the inner region is more consistent with the adiabatic HSE model but the density in the region between 5'' and 10'' deviates more from the model.

Zoom In Zoom Out Reset image size

The classical Bondi flow model without a galactic potential is shown as a dot-dashed orange line in Figure 8. Compared to the HSE models, the inner most region fits better but is still inconsistent with the data. There is a large discrepancy beyond ∼8'', due to the ignorance of the galactic potential in this model. In fact, the total enclosed stellar mass at about 1''–2'' is reaching ∼10⁹ M_☉, and therefore, the galactic potential has to be taken into account self-consistently. We demonstrate the effect of the galactic potential with a Bondi-like modeling including the galactic potential with the flow rate fixed to be the classical Bondi accretion rate (solid blue line in Figure 8). A more realistic model with stellar feedback and conduction is presented in Shcherbakov et al. (2013), although the existence of a billion solar mass black hole has to be an assumed prior in that model. Our model here is more consistent with the data in a larger radial range (≲ 10'') compared to the HSE and the classical Bondi models. In particular, the model agrees quite well with the hot component of the two-temperature, although the uncertainties of the gas density are underestimated because the uncertainties in temperature of the hotter component were not taken into account (Section 4.3). There is still a large discrepancy between the single temperature model (or the cooler component of the two-temperature model).

Our data suggest that adiabatic or isothermal HSE with the absence of a black hole is ruled out. Rotation of hot gas would give an even flatter density profile. The short dynamical timescale argues against non-HSE with the absence of a black hole. The X-ray data alone suggest that the rise in density toward the Bondi radius is more likely due to the gravitational influence of the supermassive black hole, in which the existence is supported by optical observations (Kormendy & Richstone 1992; Kormendy et al. 1996; Emsellem et al. 1999). Note that a similar technique has also been used to detect a massive black hole in the giant elliptical galaxies NGC 4649 with its Bondi radius of ∼1'' using X-ray data alone (Humphrey et al. 2008). The gravitational influence of the black hole should, however, compress and heat the hot gas at the center. It is puzzling that the single temperature profile shows a decrease in temperature toward the center. As discussed in Sections 4.2, 4.3, and 5.4, the softer component might be located more in a small disk region and the hotter component of the two-temperature model might be the more spherical accretion/outflow component.

The X-ray data alone support that we are witnessing the onset of an accretion (out)flow due to the gravitational influence of the billion solar mass supermassive black hole. The hot gas is likely to be in transition from the ambient gas in the galactic potential near the Bondi radius. Within ∼1''–2'', the galactic potential becomes negligible. Combining detailed theoretical modeling (e.g., Shcherbakov et al. 2013) together with deeper radio observations to constrain the plasma properties in the vicinity of the black hole should allow us to distinguish among different accretion models.

Our results also suggest that in all other hot accretion flows, the transition regions (around the Bondi scale) connecting the ambient gas and the asymptotic flow near the black hole should contribute significant X-ray emission. Theoretical models should take into account such a transition region in order to probe the inner most accretion region (e.g., Quataert & Narayan 2000).

Recently, Yuan et al. (2012) performed simulations of hot accretion flows which span a much larger dynamical range compared to many previous simulations. They found that all their numerical simulations show single power law self-similar profiles spanning close to the Bondi radius down to 1–10 R_S. Such results agree with many other previous simulations. In particular, they found that all radial profiles scale with similar power law indexes regardless of viscosity, magnetic field, or initial conditions. For instance, their simulated density profile scales as ρ∝r^{−3/2 + p}, with p = 0.65–0.85 from their simulations,⁹ which is consistent with our measured density profile of $\rho \propto r^{-[0.62^{+0.26}_{-0.38}]}$ within 3'' (141 pc) for our single temperature model and close to ρ∝r⁻¹ for a wider range of radii and for the hotter component of the two-temperature model. Hot gas near the Bondi radius may still be in transition from the ambient ISM to the accretion flow, and therefore the density slope near the Bondi radius may not reach the asymptotic value of the accretion flow (e.g., Bondi 1952; Quataert 2002). For accretion flow in a galactic potential, the transition may be smoother than the accretion flow in a uniform ambient ISM (e.g., Quataert & Narayan 2000). If the asymptotic density profile toward the black hole is close to the predictions by Yuan et al. (2012), the density profile of NGC 3115 may be more smoothly in transition from the region around the Bondi radius all the way toward the event horizon, and therefore a single power-law in density between ∼R_S and ∼R_B may be applicable.

For comparison, Wang et al. (2013) recently estimated a density power law index of one for Sgr A*, which generally agrees with NGC 3115. Since the spectrum in Sgr A* is not spatially resolved, they had to assume a temperature profile of the form T∝r⁻¹ to estimate the density profile. Similar to NGC 3115, the flow in Sgr A* may also be in transition near its Bondi radius so that the asymptotic T∝r⁻¹ behavior may not be applicable. A flatter temperature profile would give a flatter density profile, which is more consistent with the density profile of our single temperature model in NGC 3115.

At about 1'' (47 pc), the electron density of the single temperature model is 0.15 cm⁻³. While there are larger uncertainties in the two-temperature model and the metallicity, the density should not be off by more than an order of magnitude. The accretion rate at 1'' (47 pc) is then estimated to be $\dot{M}_{\rm acc}(47\, {\rm pc}) = 4 \pi \lambda R^2 \rho c_s = 9 \times 10^{-3} \,M_{\odot }\, {\rm yr}^{-1}$, where λ = 0.25 and γ = 5/3 for an adiabatic process, $c_s = \sqrt{\gamma k_B T/\mu {\rm m}_p}$ is the adiabatic sound speed, μ = 0.63 is the mean molecular weight, and T = 0.3 keV is assumed. The uncertainty of T also will not introduce an error in the accretion rate by more than an order of magnitude.

The mass accretion rate estimated at 1'' (47 pc) from the black hole is a factor of a few smaller than the accretion rate of (2–4) × 10⁻² M_☉ yr⁻¹ estimated around a larger radius of 4''–5''. This highlights the systematic uncertainty in estimating the accretion rate near or beyond the Bondi radius in other galaxies with unresolved/underresolved Bondi radii. The estimation near 1'' (47 pc) is closer to the upper limit of 2 × 10⁻³ M_☉ yr⁻¹ we estimated in a more self-consistent model (Shcherbakov et al. 2013).

The upper limit of the X-ray luminosity of the central AGN is 4.4 × 10³⁷ erg s⁻¹, which is about six orders of magnitude smaller than the accretion luminosity (5 × 10⁴³ erg s⁻¹) at 1'' if we assume a 10% radiative efficiency. As discussed in W11, this discrepancy can be explained if the accretion rate near the black hole is suppressed as predicted by the ADIOS or CDAF models with the scaling relation $\dot{M} \propto r^p$. Yuan et al. (2012) found in their simulations that the accretion rate is constant within about 10 R_S and the accretion rate scales as $\dot{M} \propto r^p$ beyond 10 R_S, as expected in the ADIOS model. Assuming the accretion rate at 10 R_S is suppressed by six orders of magnitude, the scaling relation gives p ≈ 1.3 for a billion solar mass black hole. Such a high value of p is outside the theoretical upper limit of 1 and also larger than most of the simulated results of 0.5–0.7 (Yuan et al. 2012). However, observational uncertainties may bring p close to one, which is consistent with the latest version of the ADIOS model (Begelman 2012). Such a high p = 1 value suggests a very flat density profile of ρ∝r^−0.5, consistent with the power law slope of $0.62^{+0.26}_{-0.38}$ (90% confidence) within 3'' (141 pc) for the single temperature model.

It is also likely that other factors are needed to explain the large discrepancy in X-ray luminosity and accretion rate. For example, accretion may be highly suppressed by rotation close to the event horizon of the black hole (Proga & Begelman 2003; Li et al. 2013) and significant outflow can also be generated during the rotational accretion process (Blandford & Begelman 1999; Li et al. 2013), although the suppression can be less effective if the gas viscosity is sufficiently high (Narayan & Fabian 2011). Numerical simulations show that rotation can flatten density and temperature profiles (Brighenti et al. 2009). Stellar feedback should also suppress the accretion as discussed in Hillel & Soker (2013) & Shcherbakov et al. (2013), with the latter authors also including conduction as a possible suppression mechanism. It is also possible that the radiation efficiency can be lower than the 10% canonical value we assumed (Ho 2008).

The discussion above was based on steady or quasi-steady state flows without feedback. However, the dynamics of hot gas within a Bondi radius is not only governed by the black hole (and the galactic) potential alone. If there is (was) a sudden release of feedback energy from the black hole (recently) as we see in other giant elliptical galaxies, dynamical disturbance will be important (e.g., Fabian 2012; McNamara & Nulsen 2012). Since there is no strong evidence of AGN feedback in NGC 3115 (e.g., strong radio source, jet, or X-ray bubble), we do not consider strong dynamical flow in this paper. For weak AGN such as NGC 3115 or Sgr A*, other feedback mechanisms such as stellar feedback or conduction can also play important roles (Hillel & Soker 2013; Soker et al. 2013; Shcherbakov et al. 2013). Therefore, resolving the gas profiles within the Bondi radius is only a step further toward the understanding of black hole accretion. Realistic theoretical modeling and simulations should be performed to match observations, and we present our effort by including conduction and stellar feedback in NGC 3115 in our companion paper (Shcherbakov et al. 2013).

We modeled the CV/AB contribution in X-ray by assuming the L_X–L_K scaling relation determined from M32. The galaxy by galaxy variation is about 30% (rms) for 11 nearby early type galaxies (Bogdán & Gilfanov 2011) and about 30%–40% for 3 nearby early type galaxies studied by Revnivtsev et al. (2008). We varied the CV/AB normalizations by ±50% which is comparable to these galaxy by galaxy variations and also larger than the statistical uncertainties of about 30% (90% confidence level) in the spectral fitting of M32. Varying the CV/AB normalizations only introduces very small systematic biases to the temperatures or gas normalizations (densities) within ∼10'' and larger ones beyond that (left column in Figure 9). All of the systematic uncertainties are within the statistical uncertainties (90% confidence).

Zoom In Zoom Out Reset image size

The power-law index of the combined spectrum of all the resolved point sources within D₂₅ is $\Gamma _{\rm LMXB} = 1.61^{+0.02}_{-0.02}$, consistent with the value of 1.6 determined in the bulge of M31 (Irwin et al. 2003). We assessed the systematic uncertainties of the power-law index of unresolved LMXBs by changing Γ_LMXB to 1.4 and 1.8 (middle column in Figure 9). Changing the LMXB power-law index to 1.8 only introduces negligible systematic uncertainties of less than 3% in temperature. A power-law index of 1.4 changes the temperature by less than 2% in general, but with a large change of 12% in the 4''–5'' bin. These are all well within the 90% confidence of the statistical uncertainties. The normalizations are generally changed by 10%–20% and at most by about 40%. This generally leads to density biases of less than 10% and at most by ∼20%–30%. All of these are comparable to the statistical uncertainties.

We also examined the systematic uncertainties of the normalizations of the unresolved LMXBs by fixing them to their upper and lower limits (90% confidence) of the nominal fitting (right column in Figure 9). Fixing the LMXB normalizations to their lower limits (solid red) generally lowers the gas temperatures and normalizations (densities), but within the statistical uncertainties. When we fixed the LMXB normalizations to their lower limits (solid green), the temperatures are biased higher, particularly in the outer regions. We noticed that in the seventh and the ninth bins, there are local minima with deviations of c-statistics from the global minimum ΔC < 1. The profiles at these local minimums (dashed green lines) are much closer to the nominal profiles. These larger biases in temperatures in the outer regions may be due to the underestimation of the hard photons from LMXBs and therefore the thermal component in the model tries to compensate for the excess hard photons with higher temperatures. The existence of the second minimum indicates the significance of a thermal component that is closer to the nominal model. Nevertheless, the gas normalization and density profiles are hardly biased by the systematic uncertainties of the unresolved LMXB normalizations.

Ideally, the absorption should be fitted from the spectrum, but thawing the absorption gives an unphysically low n_H of zero, suggesting a degeneracy between n_H and the hot gas component. Given the generally low n_H intrinsic to early type galaxies and also the insignificant amount of total H i (less than a few 10⁷ M_☉) in NGC 3115 (Roberts et al. 1991; Karachentsev et al. 2004; Sage & Welch 2006), we fixed n_H = 4.32 × 10²⁰ cm⁻² to the Galactic value (Dickey & Lockman 1990). Fitting the n_H with bright resolved LMXB in NGC 3115 generally gives a higher n_H by at most 60% but consistent with the Galactic value within 1σ–3σ in uncertainties. We assessed the systematic uncertainty in absorption by fixing n_H = 10 × 10²⁰ cm⁻² (left column in Figure 10). Even if the absorption is this high, the temperatures are only slightly biased lower and the gas normalizations (densities) are biased higher to the upper limit of the uncertainties of the nominal model. This does not change our results qualitatively.

Zoom In Zoom Out Reset image size

Fitting the metallicity of the hot gas from X-ray spectra gives a very low abundance of 0.14 solar. However, this sub-solar abundance is not expected since the hot gas should be contributed by stellar feedback which should give solar or super-solar abundance (see Shcherbakov et al. 2013 for a more detailed discussion). It is also known that metallicity determination can be biased low, particularly due to low spectral resolution of ccds and multi-temperature structure of the hot gas (Buote 2000; Su & Irwin 2013). Motivated by a more realistic physical situation, we fixed the abundance to the solar value. At temperatures below 1 keV, the emission is dominated by line emission which is proportional to metallicity. This introduces a degeneracy between metallicity and gas density (because emission is also proportional to density squared). Thawing the abundance only introduces a systematic uncertainty in temperature that is smaller than the statistical uncertainty and increases the gas normalizations by a factor of four to five (middle column in Figure 10). The derived density only increases by a factor of two without affecting the density slope. All our conclusions remain the same.

The background only contributes from less than about 1(2)% of the 0.5–2.0 (0.5–6.0 keV) emission at the central 3'' up to 6(10)% at 8''. The background increases to about 50% (60%–70%) of the emission in the 0.5–2.0 (0.5–6.0) keV energy band between 20''–40''. We changed the background level by ±10%. This generally introduces less than 2% systematic uncertainties in temperature within ∼10'', but larger beyond that (right column in Figure 10). The temperature beyond 20'' can be significantly biased. Nevertheless, the gas normalizations is biased by less than 1% within ∼10'' and at most 25% in the outermost bin (but still smaller than its statistical uncertainty). The bias in the density profile is even smaller. Note that the outermost regions only contribute a small fraction of emission to the projected emission of the inner regions. Therefore, the deprojected density profile of the inner regions is not sensitive to the precise value of the density in the outermost regions.