PROBING THE PHYSICS OF NARROW-LINE REGIONS IN ACTIVE GALAXIES. III. ACCRETION AND COCOON SHOCKS IN THE LINER NGC 1052

NGC 1052 was observed on 2013 November 2–3 using the Wide Field Spectrograph (WiFeS; Dopita et al. 2007, 2010) at the ANU 2.3 m telescope at Siding Spring Observatory. WiFeS is an optical integral field spectrograph providing a field of view of 25'' × 38'' via a total of 25 38'' × 1'' slitlets. The instrument is a double-beam spectrograph, with arms for both the blue and the red part of the spectrum. Thus, a wide wavelength range in the optical is covered simultaneously, and with adequate spectral overlap between the red and blue spectra.

The observations for NGC 1052 consisted of 3 × 1000 s exposures at R.A. = 02:41:04.9, decl. = −08:15:21.1 (J2000) and 3 × 1000 s exposures on the nearby sky at the position R.A. = 02:41:10.5, decl. = −08:11:41.4 (J2000). All observations were obtained using the R_S = 3000 grating in the blue (B3000 grating) and the R_S = 7000 grating in the red arm (R7000 grating). The position of the WiFeS data cube is shown on an image of NGC 1052 in Figure 2 with the physical scale at the assumed distance of 19.5 kpc. Absolute photometric calibration was made using the Space Telescope Imaging Spectrograph (STIS) spectrophotometric standard stars HD 009051 and HD 031128.¹ In addition, two B-type telluric standards were observed, HIP 108975 and HIP 18926. During the observation, the seeing averaged 1.3–1.5'', reasonably well matched to the 1.0'' pixels of the spectrograph.

Zoom In Zoom Out Reset image size

2.1.1. WiFeS Data Reduction

2.1.2. Results from the Data Cube

The WiFeS data cube furnishes a number of novel and interesting results. The emission-line maps reveal that the degree of excitation is in no sense uniform. Indeed, in a restricted X-shaped zone extending nearly 10'' from the nucleus, we find [O iii] lines to be strongly enhanced, approaching the values seen in Seyfert narrow-line regions; see Figure 3. We regard this as probably delineating an ionization cone photoionized by the central engine. In the nuclear region an enhanced region of [O iii] emission is seen to be associated with the arc section scale radio jet (Wrobel 1984), while Hα is enhanced along the axis of rotation of the stars (the minor axis). The [N ii] emission is relatively strongest along the major axis of the galaxy, within the X-shaped [O iii] ionization cone, in which the [O iii]/Hβ ratio rises to 3–5. The spatial extent and shape of this region are very similar to the region within which the thermal X-rays are enhanced (see Kadler et al. 2004b).

Zoom In Zoom Out Reset image size

The dynamical picture—shown in Figure 4—is equally interesting. The stars rotate smoothly around the photometric minor axis and are characterized by a velocity dispersion of ≈200 km s⁻¹, not unusual for an elliptical galaxy (Gerhard et al. 2001; Thomas et al. 2007; Cappellari et al. 2013). There is a sharp cusp in the velocity dispersion close to the nucleus, presumably within the zone of influence of the central black hole. The dynamics of the gas is completely different from the stellar dynamics. The gas velocity component analysis reveals two bubbles of gas expanding along the minor axis of the galaxy, separated by a region of very high velocity dispersion. These bubbles are about 15'' in diameter (∼1.5 kpc) and are well aligned with the Hα emission "funnel" seen in Figure 3. The bubble on the NE side is approaching, so this is located on the near side of the galaxy. This is consistent with the fact that the H₂O maser emission associated with amplification of jet radio continuum by the inner accretion disk is on the opposite side of the nucleus, also suggesting that the E-NE jet is directed somewhat toward us (Claussen et al. 1998; see Figure 1). The velocity of approach in the NE bubble and the recession velocity in the SW bubble are greatest along the axis but decrease strongly toward the boundary of the bubbles. In addition, the internal velocity field is irregular and quite unlike a rotation curve. Furthermore, both bubbles show an internal velocity dispersion of about 120 km s⁻¹, much greater than the instrumental resolution in the red (∼45 km s⁻¹). All these dynamical properties are consistent with the buoyant expansion of turbulent bubbles filled with hot or nonthermal plasma in the ISM of the galaxy.

Zoom In Zoom Out Reset image size

The region between the two expanding Hα-emitting bubbles displays a very high velocity dispersion, up to 300 km s⁻¹. In both extent and morphology it coincides with the region within which the thermal X-rays are enhanced (Kadler et al. 2004b) and the region in which the [N II] emission is unusually strong with respect to Hα, referred to above. The juxtaposition of high velocity dispersion, enhanced [N II]/Hα ratios, and the presence of thermal X-ray emission suggests the presence of strong radiative shocks in a highly turbulent medium (see Ho et al. 2014). Its inner region also corresponds to the zone of enhanced obscuration seen in Figure 1.

The spectral data used in this study were obtained under project GO-6532 (investigators: Dopita & Koratkar) using the Faint Object Spectrograph (FOS) aboard the Hubble Space Telescope (HST/FOS) on 1997 January 13 (post-COSTAR). All observations were made using the 0.8'' pair square aperture, and the wavelength coverage was 1200–6800 Å. The optical spectra (G400H and G570H grating) and the near-UV (G270H grating) data were obtained at a resolution of $\lambda /{\Delta }\lambda \sim 1300$ (or an approximate resolution 2.1–4.4 Å). Because of the faintness of the lines, the UV spectrum was obtained with the G160L grating, giving a resolution of $\lambda /{\Delta }\lambda \sim 250$ (or an approximate resolution of 6.4 Å). For NGC 1052, the aperture was well matched in size to the central bright region visible in the Hα image of Figure 1, and so it provides the integrated spectrum of the true nuclear LINER in this galaxy.

The G160L data were taken in three exposures covering two and a half orbits. These were followed by the G270H exposure in the next orbit and the G400H and G570H exposures in the final orbit. Only one guide star was acquired at the outset, so the roll of the telescope was controlled by gyroscope, which generally results in a slow drift of the target across the detector over the course of several orbits. Archival data for the two full-orbit G160L exposures were checked for systematic differences in the bright emission lines that would indicate substantial drift. Figure 5 shows that such differences are small.

Zoom In Zoom Out Reset image size

2.2.1. Extracting the Emission-line Spectrum

Figure 6 shows the FOS spectrum, along with the continuum and line model fits used to extract the line fluxes. We use a modified version of lzifu to fit the spectrum. For the continuum fitting, we adopt the MILES stellar libraries (Vazdekis et al. 2010) at wavelengths greater than λ ≈ 3500 Å. At wavelengths smaller than λ ≈ 3500 Å, we fit the continuum with linear combinations of 10 orders of Legendre polynomial. Channels around emission lines are masked out while fitting the continuum. For the emission-line fits, we split the emission lines into 14 groups based on their wavelengths and then perform separate fits to each group. In each fit, we constrain all the lines to have the same velocity and velocity dispersion. All the emission lines are described by single-component Gaussians, except for the Hα line, where a second Gaussian component is required to properly describe the line profile. Where possible, the flux ratios of the emission lines from the same multiplet of a given ion are constrained based on the relative probabilities of their atomic transitions.

Zoom In Zoom Out Reset image size

2.2.2. The FOS Emission-line Spectrum

The dereddened extracted emission-line intensities are listed in Table 1. The line fluxes have been dereddened by the theoretical function for a physical grain model derived by Fischera & Dopita (2005), which provides a very close fit to the Galactic attenuation law with R = 3.2. In common with Gabel et al. (2000), we find an extinction that is unusually high for the LINER population: A_V= 1.05. This corresponds to E(B−V) = 0.33, somewhat lower than the value of E(B−V) = 0.44 estimated by Gabel et al. (2000). The presence of such a large amount of circumnuclear dust suggests that the nuclear region also contains a large amount of interstellar gas, which we will interpret as due to a dusty accretion flow onto the bright nuclear disk. Indeed, this dusty flow can be seen directly in Figure 4 of Kadler et al. (2004a) as a complex of accretion streams within the central ∼2''.

Table 1.  Extracted Nuclear Line Fluxes for NGC 1052 (Dereddened by A_V = 1.05) Compared with the Models Described in Section 4.5

λ(Å)	Line ID	Obs. Fluxa	Model Ab	Model Bc	Model Cd	Model De
		(Hβ = 100)
1548,51f	C iv	72 ± 12g	174	119	144	123
1640	He ii	74 ± 9	72	69	73	79
1661,66	[O iii]	30 ± 10	67	57	70	60
1883,92	Si iii]	51±7	12	10	13	12
1907,09	C iii]	129 ± 7	69	67	82	92
2139,43	N ii]	38 ± 5	8	8	8	8
2323,35	C ii], Si ii]	107 ± 13	116	134	116	122
2470	[O ii]	41 ± 3	58	54	82	89
2422,24	[Ne iv]	37 ± 4	30	27	28	24
2785,93	Mg ii	67 ± 6	55	66	45	48
3345	[Ne v]	4 ± 2	5	5	5	4
3426	[Ne v]	12 ± 3	13	12	15	12
3726,29	[O ii]	403 ± 8	595	617	665	647
3869	[Ne iii]	76 ± 8	59	65	71	97
3889	He i,H i	26 ± 2	26	28	26	27
3968,70	[Ne iii], H	39 ± 8	34	36	35	43
4069,76	[S ii]	77 ± 5	109	121	84	80
4100	Hδ	19 ± 3	26	26	26	26
4340	Hγ	50 ± 3	46	47	47	47
4363	[O iii]	21 ± 3	16	13	16	15
4686	He ii	...	10	10	10	11
4861	Hβ	100 ± 3	100	100	100	100
4959	[O iii]	85 ± 3	73	64	78	76
5007	[O iii]	241 ± 4	211	186	225	220
5198,200	[N i]	31 ± 4	2	3	2	3
5755	[N ii]	12 ± 5	13	14	15	17
5876	He i]	17 ± 6	20	22	20	21
6300	[O i]	183 ± 5	183	241	121	175
6363	[O i]	54 ± 3	58	77	39	56
6548	[N ii]	111 ± 7	110	110	126	123
6563	Hα	274 ± 20	294	294	294	294
6584	[N ii]	330 ± 21	322	324	372	362
6716	[S ii]	143 ± 5	190	188	184	186
6731	[S ii]	131 ± 5	169	178	195	190

^aAbsolute dereddened flux in aperture: $F({{H}_{\beta }})=9.04\times {{10}^{-14}}$ (erg cm⁻²s⁻¹). ^bF_X = 0.07 erg cm⁻², accretion shock velocity 150 km s⁻¹, cocoon shock velocity 200 km s⁻¹. ^cF_X = 0.07 erg cm⁻², accretion shock velocity 150 km s⁻¹, cocoon shock velocity 300 km s⁻¹. ^dF_X = 0.22 erg cm⁻², accretion shock velocity 150 km s⁻¹, cocoon shock velocity 200 km s⁻¹. ^eF_X = 0.22 erg cm⁻², accretion shock velocity 150 km s⁻¹, cocoon shock velocity 300 km s⁻¹. ^fLine doublets are indicated as follows: 1548,51 = λ1548 + λ1551. ^gErrors given are measurement errors only. Errors due to reddening correction are not computed.

Download table as:  ASCIITypeset image

The line fluxes of Table 1 cannot be fit by a simple one-zone photoionization model. As realized by Gabel et al. (2000), the line ratios are inconsistent with such a model. While the [S ii] λ6717/λ6731 line ratio and the presence of a strong [O ii] λλ3726, 3729 doublet suggest a density of order $300\lesssim {{n}_{e}}\lesssim 1000$ cm⁻³, ratios such as [S ii] $\lambda 4068,76/6717,31$, [N ii] λ5755/λ6584, [O ii] λ2470/λλ3726, 3729, and [O ii] λ4363/λ5007 are all indicative of densities in the range 10⁶ $\lesssim$ n_e $\lesssim$ 10⁷ cm⁻³. Gabel et al. (2000) model these two components as physically independent entities photoionized by the same power-law spectrum and derive appreciably different ionization parameters for the two components. Surprisingly, the denser component has the higher ionization parameter, so if the model is correct, it would have to be emitted much closer to the central engine. If so, why then does the dense component not show the broader line widths expected with the Keplerian motion around the central engine? In what follows we will present a physically motivated accretion model with photoionization that appears to provide a natural answer to this question.

The WiFeS observations of NGC 1052 clearly reveal the presence of a high velocity dispersion gas along the major axis of the galaxy, associated with enhanced [N ii] λ6584 emission. In addition, these data also provide evidence for a jet-induced outflow with embedded ionization cone along the minor axis of the galaxy. The high velocity dispersion along the major axis may be due to a variety of causes. First, it could be simply line-of-sight overlap between the two expanding bubbles. Second, it could be due to a physical interaction between the ambient medium and the radio jets. Third, it could be due to intrinsic velocity dispersion and shocks in a turbulent accretion flow. In all of these scenarios, shocks play a major role, and the presence of thermal X-ray emission coextensive with the region of high velocity dispersion tends to support this picture. The temperature of this X-ray plasma (kT = 0.4−0.5 keV), if shock-heated, would require shocks with a velocity of ∼400 km s⁻¹.

Figure 1 reveals the presence of dusty (accretion?) streams, as well as both clumpy and filamentary emission-line gas in the region occupied by the thermal X-ray plasma. Figure 3 shows that the [N ii] λ6584 emission is particularly strong in this region, which is bounded by the X-shaped region of high [O iii] λ5007 emission, which coincides with the inner boundary of the expanding large-scale bubbles of gas seen in Figure 4.

We are thus led to a model that has both inflows and outflows on the large scale. An accretion flow is falling into a small, few arc section diameter circumnuclear region. The accretion flow is mostly confined within 45° of the major axis of the galaxy. At the same time, a hot plasma consisting of both relativistic gas from the nucleus and a hot X-ray plasma with kT ∼ 0.4–0.5 keV is rising buoyantly at velocities of at least 150 km s⁻¹ through the large-scale accretion flow. This outflow is goverened by the potential of the galaxy such that the bubble rises toward the minor axis of the galaxy. The interaction between these rising bubbles of hot plasma and the large-scale accretion flow is the fundamental cause of the regions of large velocity dispersion seen in the WiFeS data.

In our model, the hot plasma acts as a source of viscosity in the accreting gas, allowing the infalling gas to settle at sub-Keplerian velocity. In these circumstances, an accretion shock will form on the two faces of the dense central accretion-disk-like structure. The strongly shock-enhanced emission in the accretion shock region produces the central bright disk so evident in Figure 1 and mapped in detail in the subarcsecond resolution integral field observations of Sugai et al. (2005). These authors find that this region is characterized by a relatively slow rotation with a velocity amplitude ∼100–120 km s⁻¹ at a radial distance of ∼150 pc. The inner outflow axis defined by the radio axis is strongly inclined to the rotation axis of this accretion disk by about 50°.

Shocks may also arise internally in the central accretion disk if the orbits of the accreting material are somewhat chaotic. Such a model was suggested for the case of M87 by Dopita et al. (1997). In such a case, radiative internal shocks serve to increase the binding energy of the accepted material. In other words, these shocks act as a source of viscosity in the accreted material, helping to convey gas toward the central engine.

The misalignment of the radio jet axis and the accretion disk rotation axis has important consequences resulting from the interaction of the radio plasma with the (relatively thick) accretion disk. In the Sugai et al. (2005) data, this interaction shows two [O iii] λ5007 emission ridges associated with the Kadler et al. (2004a) radio jets. The approaching (E-NE) jet produces a ridge of [O iii] emission blueshifted with respect to the systemic velocity by −400 km s⁻¹ out to a radial distance of about 75 pc, while the receding radio jet gives a redshifted [O iii] emission ridge with a radial velocity of 200 km s⁻¹ out to a radial distance of ∼100 pc. Beyond this radius, the mean radial velocity of both ridges falls rapidly toward the systemic velocity. At the same time, the [O iii]/Hβ ratio rises rapidly.

Here we interpret these observations in terms of the cocoon shock model for the interaction of a radio jet with an ambient ISM (Bicknell et al. 1997, 1998). In this model, the radio jet passes though a termination shock and becomes surrounded by a cocoon of hot plasma that is a mixture of relativistic plasma and thermal plasma entrained or swept up by the jet material. It is this cocoon that gives rise to diffuse radio emission and possibly some of the diffuse X-rays. The optical emission arises not from this component but from a strong radiative wall shock propagating into the surrounding medium as a result of the overpressure of the cocoon material. These radiative shocks may also give rise to thermal X-ray emission if the shock velocity is high enough (≳200 km s⁻¹).

In the case of NGC 1052, the misalignment of the jet and accretion disk axes ensures that these radiative wall shocks sit in the already shocked material of the accretion disk and are much brighter on the side of the jet closest to the midplane of the accretion disk. From the observed radial velocities from Sugai et al. (2005), we infer that the cocoon shock velocity is of order ∼300 km s⁻¹. If the direction of jet propagation is almost in the plane of the sky, then the way in which the radial velocity falls to the systemic velocity at the head of the radio jets can be explained in terms of projection. On the other hand, the increase of the [O iii]/Hβ ratio can be ascribed to lower-density, faster, and/or partially radiative shocks at the head of the cocoon.

The geometry of this model is schematically illustrated for clarity in Figure 7, in which we choose a viewing angle close to the midplane of the galaxy. The misalignment of the inner accretion disk and the resulting one-sided nature of the cocoon shock is evident, as is the escape of the relativistic plasma into the expanding buoyant bubbles of gas surrounded by weak optical emission lines seen in the WiFeS observations (Figure 4).

Zoom In Zoom Out Reset image size

In our cocoon shock model, the radiative wall shock is propagating into a medium that has already been shocked to high density by the accretion flow. Thus, the initial density is orders of magnitude greater than the initial density in the accretion shocks. This two-shock model then would provide a natural explanation for the high- and low-density components observed to be present in the FOS observations.

We now proceed to develop this model in more detail. For this purpose we use the photoionization and shock modeling code MAPPINGS v4.2 (Sutherland et al., in preparation). This code now incorporates a global revision of the input atomic data for all elements up to and including zinc, as well as a full revision of the nonequilibrium ionization/cooling solver. As such, it is more comprehensive from an atomic physics point of view, as well as being faster and more stable than previous versions of the code.

All our modeling is done with a dusty plasma having a three times solar abundance set and depletion factors as defined in Dopita et al. (2013). Such high abundances are appropriate for the centers of massive and evolved galaxies in the local universe (Zahid et al. 2013). We set up a medium in which the component of magnetic pressure support and the gas pressure are equal (equipartition). The component of the magnetic field transverse to the direction of propagation of the shocks limits the maximum compression of the postshock gas that can be achieved. Defining the Alfvén Mach number ${{\mathcal{M}}^{{\rm A}}}$ as the ratio of the shock velocity to the Alfvén velocity ${{v}_{{\rm A}}}={{B}^{2}}/4\pi {{\rho }_{0}}$, where B is the transverse component of the magnetic field and ρ₀ is the preshock density, the maximum compression factor in the shock, ρ₁/ρ₀, is given by

$$\begin{eqnarray}&&\frac{{{\rho }_{1}}}{{{\rho }_{2}}}={{2}^{1/2}}{{\mathcal{M}}_{{\rm A}}}.\end{eqnarray} \tag{ 1 }$$

An accretion shock produced in a medium that is initially in equipartition is already magnetically dominated in terms of its internal pressure support by the time it has cooled enough to emit in the [S ii] λλ6717, 6731 lines. Thus, the ISM pressure that provides the best fit to the density-sensitive [S ii] λ6717/λ6731 ratio provides a good estimate of the initial pressure in the infalling gas for a given accretion shock velocity. On the basis of the observed velocity dispersion, we estimate the accretion shock velocity to be about 150 km s⁻¹, but in order to attempt to constrain it further, we investigated the effect of changing the accretion shock velocity from 80 to 260 km s⁻¹. With an accretion shock velocity of 150 km s⁻¹, the parameters of the preshock plasma are found to be T = 10⁴ K, n_H = 24 cm⁻³, and B = 36 μG, ${\rm log} P/k=6.2$. Such a pressure is not unreasonable for an elliptical galaxy. For example, in a detailed study of five nearby relaxed giant ellipticals Werner et al. (2012) find central pressures (in the hot X-ray plasma) averaging ${\rm log} P/k=6.5$.

For the atomic phase, we assume this to be preionized by the diffuse X-ray spectrum (the molecular gas is too dense for this to be a significant source of preionization in this phase). Kadler et al. (2004b) give the absorbed luminosity of L = (1.4) × 10⁴¹ erg s⁻¹ and an (intrinsic) luminosity in the X-rays of L = (1.7−2.0) × 10⁴¹ erg s⁻¹. In addition, they find a very flat power-law photon index Γ = 0.2−0.3. The measured hydrogen absorption column density of N_H ∼ (6−7) × 10²¹ cm⁻² is consistent with a nuclear extinction of at least A_V ∼ 3 mag. For our input radiation field we have adopted a theoretical spectrum for a black hole of M_BH = 10⁹ M $_{\odot }$ accreting at an Eddington factor of 10⁻³ taken from the library of Done et al. (2012) and Jin et al. (2012a, 2012b, 2012c). The form of this is rather similar to that observed by Kadler et al. (2004b).

A constraint on the radiation field intensity is provided by the [O iii] lines. As the energy density in the radiation field comes to dominate over the energy density in the shocks and the ionization parameter exceeds about ${\rm log} U=-3.0$, the [O iii] /Hβ ratio increases to eventually approach the values typical of Seyfert galaxies. At the same time, as the ionization parameter increases, the temperature in the [O iii] zone decreases so that the [O iii] λ4363/λ5007 ratio decreases as the shock component of the spectrum decreases in relative importance. This is illustrated in Figure 8 for the case v_s = 170 km s⁻¹. In this figure we have used the mixing fraction of the cocoon shock defined by the observations (see Section 4.4). The best fit to observations is obtained with F_X = 0.22 erg cm⁻².

Zoom In Zoom Out Reset image size

Since the ionization parameter characterizing the X-ray photon field is determined by the geometry of the diffuse medium with respect to the AGN—which we cannot exactly determine a priori—we treat the dilution factor of the X-ray field as a (somewhat) free parameter. In the full family of models we have used three values of the X-ray field corresponding to F_X = 0.06, 0.22, and 0.66 erg cm⁻², which is achieved at radial distances of 2.8 × 10²⁰, 2.2 × 10²⁰, and 1.3 × 10²⁰ cm from the central engine, respectively. This should be compared with the size of the central bright emission region in the Hα emission, which is roughly 3 × 10²⁰ cm in radius at the distance of NGC 1052 (Sugai et al. 2005).

The observations constrain the velocity of the cocoon shock to be between 200 and 300 km s⁻¹. Since the cocoon shock is propagating into the accretion disk material that has already been processed through the accretion shock, we set the preshock pressure equal to postshock pressure in the accretion shock and once again assume equipartition. We also set the ionization state in equilibrium with the EUV radiation field. As a result of the low ionization parameter of the local radiation field, the preshock gas is only partly ionized: H⁺/H = 0.7. This gives a set of preshock parameters: T_e = 3000 K, n_H = 28,000 cm⁻³, B = 600 μG, and ${\rm log} P/k=8.2$ cm⁻³ K. This very high pressure would have to be comparable to the pressure in the radio plasma that is driving the cocoon shock. As a result of the high preshock density, the electron density in the cooling zones of the shock exceeds 10⁶ cm⁻³, and consequently both the cooling time and cooling length in the shock are very short, of the order of 2 years and 3 × 10¹⁴ cm, respectively.

In our model, lines such as [S ii] λλ6717, 6731, [O ii] λλ3727, 3729, [N ii] λλ6548, 6584, and to a lesser extent [O iii] λλ4959, 5007 have their origin predominantly in the accretion shock. On the other hand, lines produced in dense gas such as [S ii] λλ4069, 4076, [O ii] λ2470, [N ii] λ5755, and [O iii] λ4363, as well as other lines with a high critical density such as [O i] λλ6300, 6363, are produced in the cocoon shock component. This difference provides the means of determining the relative emission fraction, ϕ, of the cocoon shocks in the combined spectrum as seen in the FOS aperture. The total flux is defined in terms of the separate fluxes from the cocoon shock and the accretion shock as ${{F}_{{\rm tot}.}}=\phi {{F}_{{\rm cocoon}}}+(1-\phi ){{F}_{{\rm accretion}}}$.

This difference provides a very direct way of constraining ϕ. In Figure 9 we show the sensitivity of five ratios to the shock velocity in the accretion flow and mixing fraction of the cocoon shock component. The mixing fraction is very tightly constrained and found to lie between $0.28\leqslant \phi \leqslant 0.35$, regardless of the shock velocity in the accretion shock, which is here varied from 80 up to 200 km s⁻¹. We cannot directly estimate the mixing fraction from the observations of Sugai et al. (2005), but these are at least consistent with the fact that the cocoon shock does not dominate within the FOS aperture.

Zoom In Zoom Out Reset image size

The best-fit model is derived from a global fit to the emission-line spectrum. For this purpose we compare the modulus of the logarithmic difference between theory and the observations of NGC 1052 for the predicted line intensities (relative to Hβ), separately for the UV and the optical lines, and attempt to minimize this difference by varying the shock velocity while keeping the molecular fraction, UV radiation field, and preshock gas pressure constant. This procedure is known mathematically as minimizing the L1 norm in the logarithmic space,²  and it allows us to fit both bright lines and faint lines simultaneously to similar goodness of fit.

The result of a run with fixed cocoon shock velocity (300 km s⁻¹) and varying accretion shock velocity is shown in Figure 10. Marginally, we obtain the best-fit model for an accretion shock velocity of 150 km s⁻¹. However, what is most striking about this figure is the insensitivity of the solution to accretion shock velocity—essentially no value in the computed range is definitely excluded.

Zoom In Zoom Out Reset image size

The result is also insensitive to the assumed cocoon shock velocity. All that we require is that the ram pressure of the shock is sufficient to compress the gas up to the range implied by the observed line ratios, which are sensitive to high density. For illustration of this point, we include in Table 1 four models that provide an almost equally good fit to the observations. Models A and B both have an EUV field F_X = 0.07 erg cm⁻² and an accretion shock velocity of 150 km s⁻¹, but they differ in their cocoon shock velocities (200 and 300 km s⁻¹, respectively). Models C and D have a stronger EUV field, 0.22 erg cm⁻², but have the same accretion velocity (150 km s⁻¹). They likewise differ in their cocoon shock velocities (200 and 300 km s⁻¹, respectively).

Note that all of these models provide a good description to the full range of ionization states observed, ranging from O i (IP = 13.61 eV) up to Ne v (IP = 126.2 eV). A very good fit is obtained in the optical, with the exception of the [N i] line—which we find to be an exceptionally difficult line to model correctly as it arises in a thin zone where hydrogen transitions from the ionized to neutral state. It is worth remarking that the intensity of the He ii λ4686 line predicted by the models and given in Table 1 also looks to be in good agreement with the actual measurement of this line by Sugai et al. (2005, their Figure 1).

Comparing our results with the Gabel et al. (2000) two-phase dusty photoionization model, we find that our model gives an appreciably improved fit to the UV lines. This is because the high-ionization zones of our models are appreciably hotter than in photoionization models, thanks to the influence of shock heating. In general, the UV lines are enhanced in our model. Likewise, the [Ne v] lines in the near-UV are much stronger in our model. In the red, Gabel et al. (2000) note that their models underestimate the [O i] λλ6300, 6363, [N ii] λλ6548, 6584, and [S ii] λλ6717, 6731 lines. Doubtless the fit would have been improved had they adopted the higher chemical abundances of these elements such as we have used in our model. Clearly a reasonable fit to the observed spectrum can be had with a two-component photoionization model, but it would still need to be made consistent with the observed spatial and dynamical structure of NGC 1052.

The model developed in the previous section and the observations of Sugai et al. (2005) can provide interesting insights on both the black hole mass and the power in the radio jet of NGC 1052. Let us first consider what constraint can be placed on the mass of the black hole in NGC 1052 from the dynamical data. Sugai et al. (2005) quote the amplitude of the rotation velocity as 100–130 km s⁻¹ at a radius of 145 pc (1.5''). The inclination of the accretion disk at this radius can only be crudely estimated as ∼30°–50° from the shape of the line-free continuum distribution shown in their Figure 2, or the HST Hα observation of Figure 1. Assuming that the gravitational potential this close to the nucleus is dominated by the black hole, we (very crudely) estimate ${{M}_{{\rm BH}}}=(4-14)\times {{10}^{8}}\;{{M}_{\odot }}$.

An alternative approach is offered by the measurement of the central velocity dispersion in the WiFeS data of Figure 4: 230 km s⁻¹. Using the relationship between black hole mass and velocity dispersion initially derived by Tremaine et al. (2002) and refined by Graham & Scott (2013) would then imply a black hole mass in the range ${{M}_{{\rm BH}}}=(2-7)\times {{10}^{8}}\;{{M}_{\odot }}$. This is probably the more reliable estimate of the two, although both are very crude. Nonetheless, they at least give an idea of the size of the central engine.

The power of the radio jet can be estimated from the jet–cocoon theory developed by Bicknell et al. (1997), in which the evolution is driven by the jet parameters, the jet energy flux $\dot{{{E}_{j}}}$, and the relativistic β = v/c. Define the r-axis to be in the direction along the jets and the z-axis to be the perpendicular direction. If the mean pressure in the cocoon is P, then toward the head of the cocoon (in the r direction) the "head" pressure in the direction of propagation is higher than the mean pressure in the cocoon. This is caused by the relatively small area, A, over which the jet termination shock is jittering and the stagnation of the cocoon flow in the direction of propagation. Here we can take the pressure as a factor ζ times the average lobe pressure, where typically ζ ∼ 2–10. However, given that the jet appears to be more directly impinging on the cocoon of NGC 1052, we will assume ζ ∼1 for the purpose of calculation. The velocity of advance of the jet and the velocity of the wall shocks surrounding the cavity excavated by the jet are given by, respectively,

$$\begin{eqnarray}&&\frac{dr}{dt}={{\left[ \frac{\beta \dot{{{E}_{j}}}}{\rho cA} \right]}^{1/2}}\sim {{\zeta }^{1/2}}{{\left[ \frac{P}{\rho } \right]}^{1/2}}\end{eqnarray} \tag{ 2 }$$

and

$$\begin{eqnarray}&&\frac{dz}{dt}={{\left[ \frac{P}{\rho } \right]}^{1/2}}\end{eqnarray} \tag{ 3 }$$

where ρ is the density in the surrounding preshock medium and $\dot{{{E}_{j}}},\beta$ and A have been defined above. For the 300 km s⁻¹ cocoon shock model, the postshock pressure P is found to be 5.2 × 10⁻⁵ dyne cm⁻², and the preshock density is ρ = 7.8 × 10⁻²⁰ gm cm⁻³. The relativistic β = 0.38 has been measured from the proper motion of hotspots in the jets at radio wavelengths in the MOJAVE survey of Lister et al. (2013). From these observations, evidence for precession in the jet is inferred. All that remains is to estimate A. Noting that the jet is highly collimated, we assume a jet opening angle of $\lesssim$ 5°, which implies an area A $\lesssim$ 10³⁸ cm². Substituting these into Equation 2 gives us $\dot{{{E}_{j}}}\lesssim 4\times {{10}^{44}}$ erg s⁻¹. This jet energy flux is comparable to those inferred for the GPS and CSS radio sources by Bicknell et al. (1997).

Let us now test whether the accretion shock arises solely at an accretion shock into the outer accretion disk or has in addition contributions from internal dissipative shocks in the accretion structure itself, as suggested for the case of M87 by Dopita et al. (1997). The key parameter that the model should be capable of reproducing is the observed absolute Hβ flux, which is a direct measure of the rate of increase of binding energy of the accreted gas. From the dereddened flux given in Table 1 and the distance of 19.5 Mpc, we estimate L_{H β} ∼ 4 × 10³⁹ erg s⁻¹. The theoretical surface flux in Hβ in our best-fit model is ${{S}_{{\rm H}\beta }}=(2-4)\times {{10}^{-3}}$ erg cm⁻²s⁻¹. The radius of the bright accretion structure from Sugai et al. (2005) is about 100 pc. Therefore, the total area of the accretion disk (top and bottom) is A ∼ 5 × 10⁴¹ cm², and we predict a luminosity of L_β ∼ (1 − 2) × 10³⁹, in fairly good agreement with what is observed. We conclude that internal dissipative shocks within the accretion disk itself cannot contribute more than 50% to the observed luminosity.

In the infall model the accretion rate is $\dot{M}=A\rho {{v}_{{\rm s}}}$, where ρ is the mean density of the infalling gas, ∼5 × 10²³ gm cm⁻³, and v_s ∼ 150 km s⁻¹ is the accretion shock velocity. This gives $\dot{M}\sim 6$ M $_{\odot }\;y{{r}^{-1}}$. This accretion rate seems remarkably high. However, we can check this by estimating the column density in the accreting material, using this to estimate the extinction it would produce, and comparing this with what is observed. For simplicity we assume that the accretion velocity remains constant, so that the accretion flow density falls off as r⁻². The hydrogen column density is therefore (approximately) the radius of the accretion disk times the preshock density n_H ∼ 24 cm⁻³. This implies a column density of ∼7 × 10²¹ cm⁻². If the gas-to-dust ratio is similar to that of the Galaxy (Bohlin et al. 1978), then A_V∼ 3, which is a little higher than the A_V∼ 1.1–1.4 inferred from the emission-line spectrum of the central disk (see Section 2.2.2) but is probably consistent with the observations when the clumpy nature of the accretion flow is taken into account.

Although the accretion rate is high, only a small fraction of this accreted matter can find its way down to feed the central black hole. The radio jet parameters of NGC 1052 estimated in Section 5.2 imply (for a black hole mass of $(2-7)\times {{10}^{8}}{{M}_{\odot }}$) an Eddington fraction of $\lesssim (0.025-0.1)$ and an accretion rate onto the black hole of $\lesssim$ 10⁻² M $_{\odot }\;y{{r}^{-1}}$. Thus, either the radio jet is very effective in ejecting the accreted gas, or we are seeing NGC 1052 in a special epoch in which the accretion rate into the outer accretion disk is very high, while the accretion from the inner accretion disk to the black hole is very low.

We have shown that the optical emission luminosity in NGC 1052 most likely results from accretion and mechanical luminosity in the jet. However, we should ask whether this scenario is valid also for other LINERs, particularly those hosted by spiral galaxies, where radio luminosities and radio jets are not so common as in elliptical galaxies, and where the accretion rates could be lower.

We have WiFeS data on two other galaxies, NGC 6812 and NGC 7213, taken from the S7 first data release (Dopita et al. 2015). These LINER galaxies have many features in common with NGC 1052 and cover a range of types. While NGC 1052 itself is an E4 galaxy, NGC 3031 is of SA(s)ab type, NGC 6812 is classified as an S0, and NGC 7213 is SA(s)a type in the de Vaucouleurs classification scheme. The remarkable similarity of their spectra (extracted from the nuclear region with an aperture of 4'' in diameter) is shown in Figure 11. NGC 7213 has been classified as a LINER or a Seyfert 1.5 (Evans et al. 1996) based on the broad Hα visible in Figure 11. However, it is clear that a broad-line LINER classification is more appropriate.

Zoom In Zoom Out Reset image size

All of these galaxies—like NGC 1052—are flat-spectrum radio sources (Healey et al. 2007). These exhibit significant polarization and are believed to show appreciable Doppler boosting related to a strong jet component. In particular, NGC 7213 contains a compact and variable radio core (Blank et al. 2005). This galaxy, like NGC 1052, also possesses an amorphous core emission region ∼130 pc in diameter.

One of the best-studied LINERs in a spiral galaxy host is NGC 3031, better known as M81. This also contains a flat-spectrum core, with clear evidence for a precessing radio jet (Martí-Vidal et al. 2011), which would set up ideal conditions for the type of jet/accretion disk interaction proposed here. M81 contains a black hole with a mass of about $2\times {{10}^{7}}\;{{M}_{\odot }}$ (Martí-Vidal et al. 2011 and references therein), which displays variable hard X-ray activity (Pian et al. 2010). In the soft X-ray component there is clear evidence for thermal plasma components with two well-constrained temperatures of 0.18 ± 0.04 keV and 0.64 ± 0.04 keV (Page et al. 2003). If produced by shocks, this would require shock velocities of 250 and 470 km s⁻¹, respectively.

M81, like NGC 1052, contains a compact (∼140pc diameter) disk of gas bright in [N ii] λ6584 and Hα imaged by Devereux et al. (1997)—see also Pogge et al. (2000). Leading into this disk is a spiral-like dust lane. The UV spectrum of this galaxy has been studied by Ho (1996) and Maoz et al. (1998). These spectra show many of the same features as NGC 1052, but with more evidence of an underlying BLR seen in Mg ii λ2800 and in Hα and Hβ and in the strong featureless continuum. We also obtained high signal-to-noise optical–UV spectroscopic data under project GO-6532 (investigators: Dopita & Koratkar) using the FOS. These data are presented in Figure 12 in the same format as Figure 6. The similarity of the two spectra is remarkable. In particular, we see the same evidence of a two-component density structure. For example, the S ii λ6731/λ6717 ratio in NGC 3031 is indicative of an electron density n_e ∼ 1200, while the S ii λλ4069, 4076/λλ6717, 6731 ratio is characteristic of densities ∼3 × 10⁵.

Zoom In Zoom Out Reset image size

Regrettably, the presence of the broad-line component makes the determination of the reddening effectively impossible in this relatively low resolution FOS spectrum, as we cannot reliably ascribe a Balmer decrement to the LINER gas. However, on the assumption that the dust extinction is negligible, we can estimate the density-sensitive line ratios, given in Table 2. These indicate densities of roughly a factor of two greater than NGC 1052 in both components of density. In the two-shock model advocated here, this is most likely accommodated by an increase in the mixing fraction of the cocoon shock (cocoon shock emission to total emission) to 0.5–0.6 (see Figure 9), presumably thanks to the lower attenuation of this component. Increasing the dust extinction would drive this estimate to even higher values. The key point, however, is that M81—like NGC 1052—requires a two-density structure that is naturally explained by the two-shock model proposed here.