Discovering AGN-driven winds through their infrared emission: II. Mass outflow rate and energetics

TL;DR

This study analyzes a large sample of AGN-driven outflows using infrared and emission line data to accurately determine their mass and energy, revealing higher densities and lower outflow rates than previously thought.

Contribution

It introduces a new method to estimate electron density in outflows, leading to revised, lower outflow mass rates and energetics, and highlights the presence of significant neutral atomic gas.

Findings

Electron densities are about 10^4.5 cm^-3, much higher than previous estimates.

Mass outflow rates are around 10^-2 solar masses per year, lower than earlier estimates.

Neutral atomic gas mass exceeds ionized gas mass, affecting outflow energetics.

Abstract

The global influence of AGN-driven outflows remains uncertain, due to a lack of large samples with accurately-determined outflow properties. In the second paper of this series, we determine the mass and energetics of ionized outflows is 234 type II AGN, the largest such sample to date, by combining the infrared emission of the dust in the wind (paper I) with the emission line properties. We provide new general expressions for the properties of the outflowing gas, which depend on the ionization state of the gas. We also present a novel method to estimate the electron density in the outflow, based on optical line ratios and on the known location of the wind. The inferred electron densities, $n_{\mathrm{e}} \sim 10^{4.5}\,\mathrm{(cm^{-3})}$, are two orders of magnitude larger than typically found in most other cases of ionized outflows. We argue that the discrepancy is due to the fact that the commonly-used [SII]-based method underestimates the true density by a large factor. As a result, the inferred mass outflow rates and kinetic coupling efficiencies are $\dot{M}_{\mathrm{out}} \sim 10^{-2}\, \mathrm{(M_{\odot}/yr)}$ and $\epsilon = \dot{E}_{\mathrm{kin}}/L_{\mathrm{bol}} \sim 10^{-5}$ respectively, 1--2 orders of magnitude lower than previous estimates. Our analysis suggests the existence of a significant amount of neutral atomic gas at the back of the outflowing ionized gas clouds, with mass that is a factor of a few larger than the observed ionized gas mass. This has significant implications for the estimated mass and energetics of such flows.

Figures

Tables

Table 1. Table 1: Best-fitting numerical coefficients of equation 2 . Additional details are provided in the appendix.

$a_1$	$a_2$	$a_3$	$a_4$	$a_5$
-3.766	0.191	0.778	-0.251	0.342

Table 2. Table 2: Correction coefficients for the [OIII] line emissivity expression given in equation 6 .

$\log U$	-2.0	-2.5	-3.0	-3.5
$f\left(\log U\right)$	0.63	0.60	0.41	0.12

Table 3. Table 3: Properties of the SEDs we consider in the modelling. Column (2): BH mass, column (3): Eddington ratio, column (4): BH spin, and column (5): mean energy of an ionizing photon.

Model	$\log M_{\mathrm{BH}}\left(M_{\odot }\right)$	$L/L_{\mathrm{Edd}}$	$a$	(Ryd)
1	8	0.12	0.7	3.15
2	7	0.14	0	2.56
3	9	0.15	0.998	2.65
4	7	0.36	0.7	4.17

Equations

$${\mathrm{E}(B-V)=\mathrm{2.33\times log\,\Bigg{[}\frac{(H\alpha/H\beta)_{obs}}{2.85}\Bigg{]}}},$$

$$\begin{split}&\log U=a_{1}+a_{2}\Big{[}\log\mathrm{\Big{(}\frac{[OIII]}{H\beta}\Big{)}}\Big{]}+a_{3}\big{[}\log\mathrm{\Big{(}\frac{[OIII]}{H\beta}\Big{)}}\Big{]}^{2}+\\ &+a_{4}\Big{[}\log\mathrm{\Big{(}\frac{[NII]}{H\alpha}\Big{)}}\Big{]}+a_{5}\Big{[}\log\mathrm{\Big{(}\frac{[NII]}{H\alpha}\Big{)}}\Big{]}^{2},\end{split}$$

$$n_{\mathrm{e}}\approx n_{\mathrm{H}}=\frac{Q(\mathrm{Lyman})}{4\pi r^{2}cU},$$

$${n_{\mathrm{e}}\approx 3.2\Big{(}\frac{L_{\mathrm{bol}}}{10^{45}\,\mathrm{erg/sec}}\Big{)}\Big{(}\frac{r}{1\,\mathrm{kpc}}\Big{)}^{-2}\Big{(}\frac{1}{U}\Big{)}\,\mathrm{cm^{-3}}},$$

$$\begin{array}{rl}& {\gamma }_{[\mathrm{OIII}]}=C_{[\mathrm{OIII}]}\times h{\nu }_{[\mathrm{OIII}]}\times \frac{n({\mathrm{O}}^{+2})}{n(\mathrm{O})}\times \frac{n(\mathrm{O})}{n(\mathrm{H})},\end{array}$$

$$\begin{split}&\gamma_{\mathrm{[OIII]}}\,(\mathrm{erg\,cm^{3}\,sec^{-1}})\approx\\ &\begin{cases}2.06\times 10^{-24}\times\Big{(}\frac{n(\mathrm{O})/n(\mathrm{H})}{4.25\times 10^{-4}}\Big{)}\times f(\log U),&\text{dusty gas}\\ 3.45\times 10^{-24}\times\Big{(}\frac{n(\mathrm{O})/n(\mathrm{H})}{4.25\times 10^{-4}}\Big{)}\times f(\log U),&\text{dustless gas}\\ \end{cases},\end{split}$$

$$M_{\mathrm{out}}=\frac{\mu m_{\mathrm{H}}{L}_{\mathrm{broad}\mathrm{line}}}{{\gamma }_{\mathrm{line}}{n}_e},$$

$$\log Q(\mathrm{Lyman})=10.06+\log L_{\mathrm{bol}},$$

Full text

1 Introduction

The discovery of correlations between the masses of super massive black holes (SMBHs) and several properties of their host galaxy (stellar velocity dispersion, bulge mass, and bulge luminosity; e.g. Gebhardt et al. 2000; Ferrarese & Merritt 2000; Tremaine et al. 2002; Gültekin et al. 2009) have led to a suggestion that the growth of the SMBH is linked to the stellar mass growth in its host galaxy (e.g., Silk & Rees 1998; Kauffmann & Haehnelt 2000; Zubovas & Nayakshin 2014). Active galactic nuclei (AGN) feedback, in the form of powerful outflows, is invoked as a way to couple the energy released by the accreting SMBH with the ISM of its host galaxy, providing a possible explanation for the observed correlations (e.g., Silk & Rees 1998; Fabian 1999; King 2003; Di Matteo, Springel & Hernquist 2005; Springel, Di Matteo & Hernquist 2005). A number of such models have successfully reproduced these correlations, by requiring that a significant amount of the accretion energy of the AGN will be mechanically-coupled to the ISM of the host galaxy ($\sim$5–10% $L_{\mathrm{bol}}$; e.g., Fabian 1999; Tremaine et al. 2002; Di Matteo, Springel & Hernquist 2005; Springel, Di Matteo & Hernquist 2005; Kurosawa, Proga & Nagamine 2009).

Powerfull winds are routinely detected in the host galaxies of AGN, spanning a large range of evolutionary stages, luminosities, and outflow gas phases (e.g., Nesvadba et al. 2006; Mullaney et al. 2013; Rupke & Veilleux 2013; Veilleux et al. 2013; Cicone et al. 2014; Harrison et al. 2014; Cheung et al. 2016; Zakamska et al. 2016; Baron et al. 2017; Fiore et al. 2017; Rupke, Gültekin & Veilleux 2017; Baron et al. 2018). In order to assess the effect of such flows on their host galaxy evolution, it is necessary to accurately determine their mass and energetics, and compare them to those required by feedback models (see however Harrison et al. 2018). A major source of uncertainty is related to the relative contributions of AGN and supernovae to the observed outflows (see introduction in Baron et al. 2018). It is normally assumed that the source that photoionizes the stationary and outflowing gas is also the main driver of the observed winds (e.g., Karouzos, Woo & Bae 2016a, b; Förster Schreiber et al. 2018). Thus, outflows that are observed in systems with emission line ratios that are consistent with AGN photoionization, are considered to be driven by the AGN.

A large number of AGN-driven winds are observed in the warm ionized phase, traced by rest-frame optical emission lines (e.g., Mullaney et al. 2013; Harrison et al. 2014; Karouzos, Woo & Bae 2016a, b; Förster Schreiber et al. 2018). In order to properly quantify the mass and energetics of these outflows, it is necessary to determine their kinematics (e.g., outflow velocity, $v_{out}$), spatial extent ($r_{out}$), and electron density ($n_{e}$), all of which are subjected to various uncertainties and systematics (see e.g., Husemann et al. 2016; Villar-Martín et al. 2016; Harrison et al. 2018; Rose et al. 2018; Tadhunter et al. 2018). For example, the electron densities in such outflows are claimed to be rather low, in the range $\sim 10^{2}-10^{3}\,\mathrm{cm^{-3}}$, based on the [SII] $\lambda\lambda$ 6717,6731Å  doublet ratio (e.g., Harrison et al. 2014; Fiore et al. 2017; Karouzos, Woo & Bae 2016b). However, recent studies suggest that the densities are higher, around $\sim 10^{3}-10^{5}\,\mathrm{cm^{-3}}$, a regime where the [SII]-based method cannot be used (Holt et al. 2011; Rose et al. 2018; Santoro et al. 2018; Spence et al. 2018).

So far, only spatially-resolved spectroscopic observations, either by long-spit observations or using integral field units (IFUs), could be used to determine the outflow extent, $r_{out}$, and thus determine the mass and energetics of outflows (e.g., Sharp & Bland-Hawthorn 2010; Fischer et al. 2011; Liu et al. 2013a, b; Harrison et al. 2014; Karouzos, Woo & Bae 2016a, b; Baron et al. 2018; Perna et al. 2019). While such observations offer a more detailed view of the outflows, they are observationally demanding. Typically, only the most extreme outflow cases are followed-up with IFUs (e.g., Mullaney et al. 2013; Harrison et al. 2014), resulting in biased samples that might not represent more typical outflow cases. Furthermore, IFU-based observations suffer from various uncertainties and systematics, such as projection effects and beam smearing, that can significantly affect the derived kinematics and extents of the outflows (e.g., Husemann et al. 2016; Villar-Martín et al. 2016; Fischer et al. 2018; Tadhunter et al. 2018).

In Baron & Netzer (2019, hereafter Paper I), we argued that the outflowing gas in active galaxies contains dust, which is heated by the central AGN, and emits in mid-infrared wavelengths. We examined the infrared spectral energy distribution (SED) of thousands of type II AGN, and showed that emission by such dust is detected in many systems that host ionized gas outflows. Our SED fitting provides the properties of the dust, including its distance from the central source. Since this dust is mixed with the outflowing gas, our method provides the luminosity-weighted (and therefore the mass-weighted) location of the outflow. This infrared emission is not subjected to various systematics affecting ground-based optical IFU observations, such as beam smearing, projection effects, and dust extinction that affects the receding part of the outflow. The estimated location of the dust, $r_{\mathrm{dust}}$, can be combined with 1D spectroscopic observations to determine the mass and energetics of ionized outflows in hundreds of systems.

In this work we focus on a subset of the sample presented in Paper I, in which ionized outflows are detected in many optical emission lines ([OIII], H$\beta$, [OI], [NII], H$\alpha$, and [SII]), and there is a clear detection of a dusty wind component at mid-infrared wavelengths. This combination allows us to accurately determine key properties of the gas in the wind, such as dust reddening, ionization state, and electron density, and thus constrain the mass and energetics of the winds. We describe our sample in section 2, and derive the spectral properties of the objects in section 3. We present new general expressions to derive the ionization parameter, electron density, and line emissivity of the ionized gas in section 4. We then estimate the outflowing gas mass, mass outflow rate, and kinetic energy of the winds in section 5. We discuss our results in section 6, and conclude in section 7.

2 Sample Selection

The sample analyzed in this paper has been described in detail in paper I. It consists of type II AGN with spectroscopically-detected outflows, for which a dusty wind component is detected in mid-infrared wavelengths. For clarity, we summarize here the main properties of the sample.

We start with the publicly-available catalog: AGN Line Profile And Kinematics Archive (ALPAKA; Mullaney et al. 2013), which provides emission line measurements for a sample of 24 264 optically-selected AGN from SDSS DR7 (Abazajian et al. 2009). Mullaney et al. (2013) performed a multi-component fitting to the emission lines [OIII]$\lambda$5007Å, [NII]$\lambda$6584Å, H$\alpha$, and H$\beta$, using continuum-subtracted spectra. The fit included one to three kinematic components to each of the emission lines, and the presence of an additional broader component in [OIII]$\lambda$5007Å  was interpreted as a spectroscopic signture of an ionized outflow. We selected systems that are classified in the ALPAKA catalog as type II AGN in the redshift range $0.05\leq z\leq 0.15$, with spectroscopically-detected outflows.

We cross-matched this sample with the MPA-JHU catalog, which provides stellar mass and star formation rate (SFR) measurements of SDSS galaxies (Brinchmann et al. 2004; Kauffmann et al. 2003b; Tremonti et al. 2004; Salim et al. 2007). We also used optical and infrared (IR) photometric data as follows. We used the $riz$ optical photometry by the SDSS, the $JHK$ near-infrared photometry from the 2MASS All-Sky Extended Source Catalog (Skrutskie et al. 2006), and the W1-W4 mid-infrared photometry from WISE (Wright et al. 2010). We used the astrometric cross-matches provided by the SDSS111https://skyserver.sdss.org/CasJobs to cross-match between systems observed by SDSS and systems observed by 2MASS and WISE.

In Paper I, we examined the infrared SEDs of 2 377 systems that show spectroscopic signatures of an ionized wind. The best-fitting SED included contributions from direct stellar light, which dominates the SED at optical and near-infrared wavelengths, torus and NLR dust emission, which dominate the SED at mid-infrared wavelengths, and dust in star forming (SF) regions which dominates the far-infrared wavelengths. The IR SED also included a new contribution from dust that is mixed with the outflow, is heated by the AGN, and emits at mid-infrared wavelengths. We found that 2 044 systems require the additional dusty wind component. Out of these, we were able to constrain the dust temperature, covering factor, and location in 1 696 systems. The initial sub-sample addressed in this paper consists of these 1 696 type II AGN.

3 Spectral properties

The goal in paper II is to determine the level of ionization, the density, the mass outflow rate and the kinetic energy of the outflowing gas. To study the ionization state and density of the gas, we must detect the outflow in various emission lines: [OIII], [NII], [SII], [OI], H$\alpha$, and H$\beta$ (see e.g., Baron et al. 2017). The outflowing gas is traced by broader emission lines that are blended with narrower emission lines coming from the stationary NLR. The broader emission lines are usually weaker than their corresponding narrow lines, and their characterization requires a careful emission line decomposition procedure. The ALPAKA catalog provides such a decomposition. However, Mullaney et al. (2013) fit the continuum emission locally using a 5 degree polynomial, which they then subtract to obtain an emission line spectrum. Since they do not fit the spectra with a stellar population synthesis model, stellar absorption lines are not properly accounted for, and may affect the resulting emission line spectrum, particularly around the H$\alpha$ and H$\beta$ emission lines. The broader H$\beta$ line is usually the weakest emission line in the spectrum, and its equivalent width is of the same order as the stellar H$\beta$ absorption.

To improve the accuracy of the broad H$\beta$ emission line, we fit stellar population synthesis models to the 1 696 type II AGN sample, and performed our own emission line decomposition. We use Penalized Pixel-Fitting stellar kinematics extraction code (pPXF; Cappellari 2012), which is a public code for extracting the stellar kinematics and stellar population from absorption-line spectra of galaxies (Cappellari & Emsellem 2004). It uses the MILES library, which contains single stellar population (SSP) synthesis models and covers the full range of the optical spectrum with a resolution of full width at half maximum (FWHM) of 2.3Å  (Vazdekis et al. 2010). We use SSP models with stellar ages that range from 0.03 to 14 Gyr, thus ensuring that we can properly describe spectra of systems with different star formation histories. The wide wavelength range of the SDSS spectra and the SSP models includes the Balmer absorption lines: H$\zeta$, H$\epsilon$, H$\delta$, H$\gamma$, H$\beta$, and H$\alpha$, ensuring that the fitted stellar continuum around the H$\beta$ absorption line is robust. The output of the code includes the relative weight of stars with different ages, the stellar velocity dispersion, the dust reddening towards the stars, and the best-fitting stellar model. Since we are not interested in the stellar properties of the galaxies in our sample, we only use the best-fitting stellar models and subtract them from the spectra.

The subtracted emission line spectra show varios emission lines: [OIII] $\lambda\lambda$ 4959,5007Å, $\mathrm{H\alpha}$ $\lambda$ 6563Å, [NII] $\lambda\lambda$ 6548,6584Å, [OII] $\lambda\lambda$ 3725,3727Å, $\mathrm{H\beta}$ $\lambda$ 4861Å, [OI] $\lambda\lambda$ 6300,6363Å, and [SII] $\lambda\lambda$ 6717,6731Å  (hereafter $\mathrm{[OIII]}$, $\mathrm{H\alpha}$, $\mathrm{[NII]}$, $\mathrm{[OII]}$, $\mathrm{H\beta}$, $\mathrm{[OI]}$, and $\mathrm{[SII]}$). We model each emission line as a sum of two Gaussians - one which represents the narrow component and one that represents the broader, outflowing, component. We perform a joint fit to all the emission lines under several constrains: (1) we force the intensity ratio of the emission lines [OIII] $\lambda\lambda$ 4959,5007Å  and [NII] $\lambda\lambda$ 6548,6584Å  to the theoretical ratio of 1:3, (2) we tie the central wavelengths of all the narrow and broad lines so that the gas has the same systemic velocity, and (3) we force the widths of all the narrow lines to show the same velocity dispersion, and we do the same for the broader lines.

Out of our initial sample, we detect a broad component with at least 2$\sigma$ in [OIII], [NII], [SII], H$\alpha$, and H$\beta$, in 234 systems. We focus on these 234 systems throughout the paper. In figure 1 we show the distribution of redshift, SFR, stellar mass, and dust-corrected narrow [OIII] luminosity for this subsample. The distributions are roughly similar to the distributions of these properties in the parent sample of 1 696 type II AGN. The two-sample KS test for the distributions give p-values of 0.53, 0.04, 0.50, and 0.03, for the redshift, SFR, stellar mass, and [OIII] luminosity respectively.

In figure 2 we show the narrow and broad emission lines on a line-diagnostic diagram (Baldwin, Phillips & Terlevich 1981; Veilleux & Osterbrock 1987). We show three separating criteria. The first is a theoretical upper limit which separates starbursts and AGN-dominated galaxies (Kewley et al. 2001, Ke01; black line). The second is a modified criterion which includes composite galaxies showing contributions from both SF and AGN (Kauffmann et al. 2003a, Ka03; grey line). The third is a line from Cid Fernandes et al. (2010) to separate between LINERs and Seyferts (CF10; dashed black). Clearly, both the stationary (narrow lines) and the outflowing (broad lines) gas in most of the 234 systems are dominated by AGN photoionization, with the majority of the systems classified as Seyferts and some as LINERs.

In the left panel of figure 3 we show the FWHM of the broad emission lines. Since the fit requires all the broad lines to have the same velocity dispersion, the FWHM is the same for all the broad lines. In the right panel of figure 3 we show the velocity shift of the broad emission lines with respect to the systemic velocity, defined by the centroid of the narrow emission lines. As in the left panel, this velocity shift is the same in all the lines.

Next, we use the measured H$\alpha$/H$\beta$ flux ratios to estimate the dust reddening towards the two kinematic components. Assuming case-B recombination, a gas temperature of $10^{4}$ K, a dusty screen, and the Cardelli, Clayton & Mathis (1989, CCM) extinction law, the colour excess is given by:

[TABLE]

where $\mathrm{(H\alpha/H\beta)_{obs}}$ is the observed line ratio. In the left panel of figure 4 we show the distribution of the colour excess for the narrow and broad lines respectively. For comparison, we show the distribution of the colour excess for the narrow emission lines in the ALPAKA catalog. Our narrow-line-based distribution of $\mathrm{E}(B-V)$, and the one from the ALPAKA catalog are roughly consistent, with somewhat larger number of sources with higher $\mathrm{E}(B-V)$ in the latter distribution. We attribute this difference to the different procedures used to fit and subtract the stellar continuum emission. For this reason, we do not compare the colour excess distributions in the broad emission lines.

In the right panel of figure 4 we show the colour excess towards the broad lines versus the colour excess towards the narrow lines. Surprisingly, we find no correlation between the two. This is consistent with the findings of Rose et al. (2018), who analysed a sample of AGN-driven outflows in ultra luminous infrared galaxies. The figure clearly shows that the narrow and the broad lines suffer from different amounts of extinction (see also Baron et al. 2017, 2018). Various earlier studies assume the same amount of dust extinction, and use the estimated colour excess in the narrow lines to correct the broad line luminosities for dust extinction. Figure 4 suggests that such a procedure is not justified.

We derived dust-corrected line luminosities using the measured $\mathrm{E}(B-V)$. We propagate the uncertainties on the best-fitting parameters to obtain the uncertainties of the line luminosities, which are typically in the range 0.05–0.1 dex. We use the narrow emission line luminosities, which are based only on the narrow components from the multi-Gaussian fits, to estimate the bolometric luminosity of the AGN with two different methods. The first uses the narrow H$\beta$ and [OIII] luminosities (equation 1 in Netzer 2009), and the second uses the narrow [OI] and [OIII] luminosities (equation 3 in Netzer 2009). In figure 5 we show the distribution of AGN bolometric luminosity in our sample. The two methods give consistent estimates. One can see that our sample spans about two orders of magnitude in AGN bolometric luminosity, roughly from $\log L_{\mathrm{bol}}=43.5$ to $\log L_{\mathrm{bol}}=45.5\,\mathrm{erg/sec}$.

4 Gas properties from photoionization modelling

In the next section we estimate the mass outflow rates and energetics of the observed outflows, which depend on various properties of the emission line gas, such as its electron density and level of ionization. Studies usually assume constant fiducial values for the electron densities and line emissivities in all the objects in the sample (see e.g., Harrison et al. 2014; Karouzos, Woo & Bae 2016a; Fiore et al. 2017), which correspond to specific gas properties that might not be appropriate in a general case. In this section we provide more general expressions for the electron densities and line emissivities, which depend on the observed properties of individual sources.

Many of the gas properties depend on the ionization parameter, $U=Q(\mathrm{Lyman})/4\pi r^{2}n_{\mathrm{H}}c$, where $Q(\mathrm{Lyman})$ is the number of the hydrogen-ionizing photons, $n_{\mathrm{H}}$ is the hydrogen density, and $c$ is the speed of light. In section 4.1 we show that the ionization parameter can be constrained using observed emission lines ratios. We then present a novel method to estimate the electron density in the gas, which is based on the estimated ionization parameter and on the known location of the gas, in section 4.2. Finally, we provide more general expressions of the H$\alpha$ and [OIII] line emissivities, where the latter depends on the estimated ionization parameter, in section 4.3.

4.1 Ionization parameter estimation

The expressions we provide are based on extensive photoionization modelling, which details are listed in appendix A. We use a ”single cloud” approach, where we calculate the properties of spherical shells of gas, with different density, pressure, column density, and distance from the central ionizing source. These models do not intend to replace more complex realistic models, that account for the geometrical distribution of clouds, with a range of densities and locations (see Netzer 2013 for a more complete description of multi-component photoionization models). Nevertheless, they provide useful expressions that can be used in the more realistic scenarios.

We focus on models with optically-thick geometrically-thin shells of gas, and do not model cases in which the gas is optically-thin and extends over regions that are as large, or even larger, than the distance from the central source (e.g., see Baron et al. 2018). The latter models require information on the gas spatial distribution, which is not available in this case. In addition, the outflowing gas can be a combination of matter-bounded (optically-thin) and radiation-bounded (optically-thick) clouds (see e.g., Binette, Wilson & Storchi-Bergmann 1996). We did not consider such a combination since (1) the optically-thick models we examine provide a general and consistent description of the four emission lines studied here ([OIII], [NII], H$\alpha$, and H$\beta$), and (2) the contribution of matter-bounded clouds to Balmer line luminosities, and to the outflowing dust emission, must be small unless the covering factor of this component is considerably larger than the one found here for radiation-bounded gas. This is inconsistent with our measurements, and thus the presence of such a component must have a small effect on the expressions derived here.

We examined a large range in gas metallicity (from 0.5$\mathrm{Z_{\odot}}$ to 2$\mathrm{Z_{\odot}}$), ionization parameter (from $\log U=-3.8$ to $\log U=-2$), and several different slopes of the ionizing continuum (with mean energy of an ionizing photon of 2.56 Ryd to 4.17 Ryd). We compared the predictions of constant-density versus constant-pressure models, and found consistent results within the range of properties we considered. Our analysis suggests that the expressions we provide are general, and do not depend on specific assumptions made during the modelling.

The [OIII]/H$\beta$ and [NII]/H$\alpha$ emission line ratios are based on the strongest lines in the optical spectra of AGN. These ratios are insensitive to dust reddening and are widely used to distinguish between AGN and SF as the main source of photoionization (Baldwin, Phillips & Terlevich 1981; Veilleux & Osterbrock 1987; Kewley et al. 2001; Kauffmann et al. 2003a; Cid Fernandes et al. 2010). The [OIII]/H$\beta$ ratio is mostly sensitive to the ionization parameter in the gas, but it also depends on the gas metallicity. The [NII]/H$\alpha$ ratio is mostly sensitive to the gas metallicity, but it also depends on the ionization parameter. For AGN-dominated systems, both ratios depend slightly on the shape of the ionizing continuum (e.g., see figure 10 and Groves, Dopita & Sutherland 2004a, b; Kewley et al. 2013).

We looked for an expression that ties the ionization parameter in the gas with the observed [OIII]/H$\beta$ and [NII]/H$\alpha$ line ratios. It is given by:

[TABLE]

where $a_{1}$, $a_{2}$, $a_{3}$, $a_{4}$, and $a_{5}$ are constants which are determined by a likelihood method. In appendix A we provide details of this fitting process and derive best fitting values for the five constants. These are listed in table 1. Equation 2 is only valid for ionization parameters in the range $\log U=-2$ to $\log U=-3.8$. The typical deviation in $U$, given the range of parameters listed in the appendix, and the two types of models (constant density and constant total pressure) is 0.11 dex. This includes a large range of metallicities and several SED shapes.

4.2 Electron density estimation

We use an independent method to estimate the electron density in the gas using optical emission lines. The method relies on the ionization parameter of the ionized gas, and the known location of the gas with respect to the central source. Using the definition of the ionization parameter, $U$, the electron density in the gas is given by:

[TABLE]

where $Q(\mathrm{Lyman})$ can be estimated from the AGN bolometric luminosity and the assumed SED (see equation 8 and related text in appendix A). For the typical AGN in our sample, and the assumed SEDs, equation 3 can be written as:

[TABLE]

where $L_{\mathrm{bol}}$ is the AGN bolometric luminosity and $r$ is the known distance of the gas from the central source. The ionization parameter is estimated using the observed emission line ratios (equation 2). Equation 4 can be used to estimate the electron density of both stationary and outflowing gas, as long as the AGN is the main source of ionizing radiation.

4.3 Line emissivities estimation

In this section we provide more general expressions for the line emissivities in ionized gas. We focus on the H$\alpha$ and [OIII] emissivities, $\gamma_{\mathrm{H\alpha}}$ and $\gamma_{\mathrm{[OIII]}}$, as these emission lines are typically used to estimate mass outflow rates of ionized outflows. Both emissivities depend on the electron temperature in the gas, with $\gamma_{\mathrm{H\alpha}}\propto T_{e}^{-1}$ and $\gamma_{\mathrm{[OIII]}}\propto e^{-1/T_{e}}/\sqrt{T_{e}}$ (e.g., Draine 2011). The electron temperature depends on the gas density, metallicity, and the ionization parameter. Our photoionization models suggest that the electron temperature in the ionized gas is in the range 10 000 K to 14 000 K, consistent with observations (e.g., Perna et al. 2019). We therefore select $T_{e}$ = 12 000 K, resulting in H$\alpha$ emissivity of $\gamma=3\times 10^{-25}\,\mathrm{erg\,cm^{3}\,sec^{-1}}$.

The [OIII] line emissivity is given by:

[TABLE]

where $n_{\mathrm{e}}C_{\mathrm{[OIII]}}$ is the collisional excitation rate of [OIII], which depends on the electron temperature exponentially ($C_{\mathrm{[OIII]}}\propto e^{-1/T_{e}}/\sqrt{T_{e}}$; Draine 2011), $h\nu_{\mathrm{[OIII]}}$ is the photon energy, $n(\mathrm{O^{+2}})/n(\mathrm{O})$ is the fraction of $\mathrm{O^{+2}}$ ions, and $n(\mathrm{O})/n(\mathrm{H})$ is the oxygen abundance relative to hydrogen. The emissivity depends on the temperature, primarily through $C_{\mathrm{[OIII]}}$, on the gas metallicity, which sets $n(\mathrm{O})/n(\mathrm{H})$, and on the ionization parameter, which affects $n(\mathrm{O^{+2}})/n(\mathrm{O})$. We note that different ionizing SEDs result in somewhat different $n(\mathrm{O^{+2}})/n(\mathrm{O})$ ratios, but this effect is of secondary importance.

Cano-Díaz et al. (2012) provided an estimate for the [OIII] line emissivity which is based on equation 5. They assumed that $n(\mathrm{O^{+2}})/n(\mathrm{O})\sim 1$, which is a reasonable assumption for their specific system, but does not apply in the general case. In particular, this assumption is not valid for most of the objects in our sample. Furthermore, they estimated the emissivity assuming solar abundance of oxygen and did not take into account that a non-negligible fraction of the oxygen atoms are depleted onto dust grains. Assuming electron temperature of 12 000 K, solar abundance of oxygen, $n(\mathrm{O})/n(\mathrm{H})=4.25\times 10^{-4}$, and a mean energy of an ionizing photon of 2.56 Ryd (model 2 SED; see appendix), the [OIII] line emissivity is given by:

[TABLE]

where $f(\log U)$ represents the dependence on the ionization parameter, and is listed in table 2. The ionization parameter is estimated from observed line ratios in equation 2. The [OIII] emissivity for the dusty gas case is a factor of $\sim 0.6$ smaller than the emissivity for the dustless gas case. This is a direct result of the depletion of oxygen onto dust grains, where only 60% of the oxygen atoms remain in the gas phase. In the next section, we use the [OIII] emissivity for the dusty gas case. Finally, we note that the constants in equation 6 represent the average [OIII] emissivity in the models we considered, with the assumed ionizing SED and gas metallicity. While these give a general description of the ionized gas, they cannot replace full photoionization models, where the emissivities can be different by a factor of $\sim$1.5.

5 Mass outflow rate and energetics

Having separated the stationary and outflowing gas components, we now turn to estimate the outflowing gas mass and its kinetic energy. The ionized gas mass in the outflow can be estimated either from the broad H$\alpha$ luminosity, or from the broad [OIII] luminosity. The gas mass is given by (see e.g., Baron et al. 2017; Fiore et al. 2017):

[TABLE]

where $\mu$ is the mass per hydrogen atom, which we fix at 1.4, $m_{\mathrm{H}}$ is the hydrogen mass, $L_{\mathrm{broad\,line}}$ is the extinction-corrected broad line luminosity ($L_{\mathrm{H\alpha}}$ or $L_{\mathrm{[OIII]}}$), $n_{\mathrm{e}}$ is the electron density in the outflowing gas, and $\gamma_{\mathrm{line}}$ is the effective line emissivity ($\gamma_{\mathrm{H\alpha}}$ or $\gamma_{\mathrm{[OIII]}}$), which is given in section 4.3. The [OIII] line emissivity depends on the gas metallicity, which is typically unknown. Lacking additional information, we assume solar metallicity.

We also estimate the mass outflow rate, the kinetic power of the wind, and the kinetic coupling efficiency. The mass outflow rate is given by $\dot{M}_{\mathrm{out}}=M_{\mathrm{out}}/t_{\mathrm{out}}$, where $t_{\mathrm{out}}=r_{\mathrm{out}}/v_{\mathrm{out}}$. The kinetic power in the wind is given by $\dot{E}_{\mathrm{kin,\,out}}=\frac{1}{2}\dot{M}_{\mathrm{out}}v_{\mathrm{out}}^{2}$. The kinetic coupling efficiency relates the wind kinetic power to the AGN bolometric luminosity: $\epsilon=\dot{E}_{\mathrm{kin,\,out}}/L_{\mathrm{bol}}$ (see e.g., Harrison et al. 2018).

The above estimates require knowledge of the location of the wind, the electron density in the wind, and the effective outflow velocity. Since we use spatially-integrated (1D) spectroscopy and photometry, our estimates are based on luminosity-averaged properties, in particular the mass-weighted average properties of the wind. We describe in section 5.1 our estimate of the wind location, which is based on the method presented in Paper I. We then describe various methods to estimate the electron density in the wind in section 5.2. Finally, in section 5.3 we discuss our definition of outflow velocity.

5.1 Effective wind location

In Paper I we argued that the outflowing gas in active galaxies contains dust, which is heated by the central AGN, and emits at mid-infrared wavelengths. We have shown that this emission component is detected in many type II AGN that host ionized gas outflows. This new component offers novel constrains on the outflow properties. Specifically, we used the dust temperature obtained from SED fitting and the AGN bolometric luminosity to estimate the distance of this dust from the central source. Our method provides luminosity-weighted, and therefore mass-weighted distance from the center of the galaxy. Since this dust is mixed with the outflowing gas, the method provides an estimate of the mass-weighted location of the outflow. In figure 6 we show the distribution of outflow locations, traced by the dust infrared emission, for our sample. The distribution is similar to the distribution of our parent sample in Paper I. As discussed in Paper I, the uncertainty on the location is roughly 0.25 dex.

5.2 Electron density

The electron density of the outflowing gas is a major source of uncertainty. The most robust estimates are obtained from full photoionization modelling of the outflowing gas (e.g., Baron et al. 2018; Revalski et al. 2018). Such models require detailed, spatially-resolved spectroscopic observations, which are available only for a small number of objects, and hence they cannot be used to explore the entire range of population properties. Most earlier studies either assume a constant value for the electron density for all objects in the sample (e.g., Liu et al. 2013b; Harrison et al. 2014; Husemann et al. 2016; Fiore et al. 2017), or use the ratio of the weak [SII] emission lines to put limits on the density (e.g., Nesvadba et al. 2006, 2008; Karouzos, Woo & Bae 2016b). Typical [SII]-based estimates are in the range 200–1000 $\mathrm{cm^{-3}}$ (e.g, Fiore et al. 2017). Here we discuss the merits and the limitations of this and other methods used to determine the gas density in the outflow

5.2.1 [SII] line ratio

In the left panel of figure 7 we show the [SII] density diagnostic diagram. The [SII]-based estimator is sensitive to electron densities which are close to the critical density of the [SII] $\lambda$ 6717,6731Å  lines: about 1600 and $1.5\times 10^{4}\,\mathrm{cm^{-3}}$ for a gas temperature of $10^{4}$ K (Draine 2011). The [SII] ratio saturates at the low ($n_{e}<10^{2}\,\mathrm{cm^{-3}}$) and high ($n_{e}>10^{4}\,\mathrm{cm^{-3}}$) density regimes, and cannot be used as a density estimator in these ranges. We use the broad [SII] lines measured in section 3 to estimate the electron density in the outflow, and show 20 randomly-selected measurements on the diagram. Due to the significant uncertainty in the emission line decomposition, the ratios are consistent with the entire range of possible ratios, from roughly 0.5 to 1.5. In particular, most measurements are consistent with the high density limit of about $10^{4}\,\mathrm{cm^{-3}}$.

In the right panel of figure 7 we show two distributions of the electron density based on the [SII] lines. The first is using the initial guess provided to the fitting function and without taking account the large uncertainties. This distribution peaks at roughly $10^{2.5}$ $\mathrm{cm^{-3}}$, similar to the estimates by others (e.g. Fiore et al. 2017). Moreover, the exact location of the fit depends on the initial guess. We also show the distribution of [SII]-based densities after accounting for the measurement uncertainties. This, and the information provided in the left panel of the diagram, demonstrate that electron densities above about $10^{3.5}\,\mathrm{cm^{-3}}$ cannot be constrained using the [SII] line ratio.

The [SII]-based electron density is subjected to an additional systematic. For the range of ionization parameters considered here, the [OIII] and H$\alpha$ emission lines are emitted throughout most of the ionized cloud, while the [SII] lines are primarily emitted close to the ionization front in the cloud, where the electron density decreases with increasing depth into the cloud. Thus, the [OIII] and H$\alpha$ lines tend to trace higher electron density regions, compared to the [SII] lines. The difference between the electron density traced by [OIII] and H$\alpha$ and that traced by the [SII] becomes more significant for increasing ionization parameter, and for increasing hydrogen density in the cloud. Since the mass of the outflowing gas is estimated using the [OIII] or the H$\alpha$ emission line luminosities (equation 7), the combination of [OIII] or H$\alpha$ line luminosity with [SII]-based electron density may result in inconsistent gas mass estimates, with significant deviations in cases with high ionization parameters and large hydrogen densities.

5.2.2 [OIII]/H$\beta$, [NII]/H$\alpha$, and $r_{\mathrm{dust}}$

We now use the novel method presented in section 4.2 to estimate the electron density in the outflowing gas. Given the coefficients in table 1, we can use equations 2 and 4 to estimate the electron density, where we define the distance of the outflowing gas to be $r_{\mathrm{dust}}$, the mass-weighted location of the dust in the wind (see section 5.1). We propagate the uncertainties of the line ratios, the effective location of the dust, and the analytical expression of the ionization parameter, to obtain the uncertainty on the electron density, which is roughly 0.6 dex. In the right panel of figure 7 we show the distribution of the electron densities in our sample using this method. The distribution is centered around $n_{\mathrm{e}}\sim 10^{4.5}\,\mathrm{cm^{-3}}$, suggesting that the [SII]-based method, shown on the left side of the diagram, underestimates the true electron density in the outflowing gas by roughly 2 orders of magnitude.

5.2.3 Auroral and transauroral [OII] and [SII] emission line ratios

Holt et al. (2011) discussed the uncertainties involved in estimating the electron density using the [SII] $\lambda$ 6717,6731Å  emission lines. They pointed out that this estimator is not sensitive to $n_{\mathrm{e}}\gtrsim 10^{4}\,\mathrm{cm^{-3}}$. To estimate the electron density, they used the transauroral and auroral emission lines [SII] $\lambda$4068,4076Å  and [OII] $\lambda$7318,7319,7330,7331Å, which have critical densities around roughly $10^{6}\,\mathrm{cm^{-3}}$ (e.g., Draine 2011; Holt et al. 2011). Their method is based on the summed flux of each doublet, while the traditional [SII] method is based on the doublet flux ratio. As such, their estimated electron density does not suffer from the uncertainty involved in emission line decomposition, which is significant when using the traditional [SII] method. Their method is sensitive to a broad range of electron densities, from roughly $10^{2}\,\mathrm{cm^{-3}}$ to $10^{6}\,\mathrm{cm^{-3}}$. However, it is sensitive to reddening corrections.

Rose et al. (2018) applied this method to study the AGN-driven outflows in 9 ultra luminous infrared galaxies (ULIRGs) observed with VLT/Xshooter. The exceptional quality of the VLT/Xshooter observations allowed them to perform emission line decomposition of the weak emission lines [SII] $\lambda$4068,4076Å  and [OII] $\lambda$7318,7319,7330,7331Å, and thus provide accurate estimates of the electron density in the outflowing gas. They found relatively high electron densities in the wind, from roughly $\log n_{e}(\mathrm{cm^{-3}})=3.4$ to $\log n_{e}(\mathrm{cm^{-3}})=4.8$ (see also Spence et al. 2018). Santoro et al. (2018) applied the method to an additional AGN, finding a remarkably high density of $\log n_{e}(\mathrm{cm^{-3}})\simeq 5.5$. These estimates exceed the densities typically assumed or measured for ionized outflows in active galaxies. Our independent method confirms this in a much larger sample of non-ULIRG type II AGN. The densities we infer are somewhat larger than those found by Rose et al. (2018). A detailed comparison between the samples should be carried out with caution, since the two have very different size and consist of different types of objects, i.e. low luminosity type II AGN versus ULIRGs.

The two methods are sensitive to a large range of electron densities, and are expected to give consistent results. There are several differences between the method presented in Holt et al. (2011) and our method. First, we assume that the distance of the gas from the central source is known, which is not required in the auroral and transauroral emission lines method. Furthermore, the method presented in Holt et al. (2011) can be used to put constrains on the dust reddening in the gas, while we need to estimate it using the observed H$\alpha$/H$\beta$ line ratio. However, the auroral and transauroral emission lines are very weak, and can only be detected in high-quality observations. In addition, this method requires a larger wavelength coverage, roughly from 3700Å  and 8000Å. Therefore, while less general, our method is expected to be more practical in many cases.

5.3 Velocity

There are several uncertainties associated with the estimate of the bulk velocity of the outflow (e.g., Liu et al. 2013b; Harrison et al. 2014; Karouzos, Woo & Bae 2016b; Fiore et al. 2017; Baron et al. 2018). First, since we work with spatially-integrated spectra, the wind geometry is unknown and we cannot correct for projection effects. Second, dust-extinction may affect the observed emission line profiles, resulting in underestimation of the outflow extent and kinematics.

Various studies define the bulk velocity in different ways (e.g., Liu et al. 2013b; Karouzos, Woo & Bae 2016b; Fiore et al. 2017; Harrison et al. 2018). We follow the definition by Karouzos, Woo & Bae (2016b), since the properties of their sample, such as wind velocity, AGN bolometric luminosity, and redshift, closely match the properties of our sample. According to Karouzos et al. (2016b), $v_{\mathrm{out}}=\sqrt{v_{\mathrm{shift}}^{2}+\sigma^{2}}$, where $v_{\mathrm{shift}}$ is the velocity shift of the broad emission line centroid with respect to the narrow lines (see figure 3), and $\sigma$ is the velocity dispersion of the broad emission lines. According to Karouzos et al. (2016b), using different definitions can result in mass outflow rates and kinetic powers which differ from the one used here by a factor of 1–3.

Having estimated the wind location and velocity, and the electron density in the outflow, we can now estimate the gas mass, mass outflow rate, and kinetic power of the wind. In figure 8 we show the distribution of ionized gas mass, mass outflow rate, kinetic power, and the kinetic coupling efficiency for the objects in our sample. We show two estimates for these properties, one based on the dust-corrected broad H$\alpha$ luminosity (black histogram), and one based on the dust-corrected broad [OIII] luminosity (grey histogram). The [OIII]-based estimates are consistent with the H$\alpha$-based estimates. This suggests that our assumption of solar metallicity is justified.

Based on the broad H$\alpha$, we find that the median ionized gas mass in the outflow is roughly $M_{\mathrm{out}}\sim 10^{4}\,\mathrm{M_{\odot}}$, the median mass outflow rate is $\dot{M}_{\mathrm{out}}\sim 10^{-2}\,\mathrm{M_{\odot}/yr}$, and the median kinetic coupling efficiency is roughly $\epsilon=\dot{E}_{\mathrm{kin}}/L_{\mathrm{bol}}\sim 10^{-5}$. These values are 1–2 orders of magnitude lower than typical estimates in type II AGN with similar bolometric luminosity (see e.g., Karouzos, Woo & Bae 2016b; Fiore et al. 2017; Harrison et al. 2018). These differences are the result of the much higher ionized gas density found in this work; roughly 2 orders of magnitude higher than in most previous studies.

6 Discussion

In paper I of this series we have shown that the IR SED of the dusty wind component in a large number of type-II AGN can be used to place constraints on the wind location. Here, in paper II, we analysed a sub-sample of 234 type II AGN with signatures of ionized gas outflows in their optical spectra, and with a detection of a dusty wind component in mid-infrared wavelengths, with well determined temperature, covering factor, and location. The detection of the outflow in multiple emission lines (H$\beta$, [OIII], [OI], [NII], H$\alpha$, and [SII]) allowed us to accurately determine key-properties of the ionized flow. In particular, we introduced a novel method to estimate the ionization parameter of the ionized gas. The combination of known ionization parameter and location allow us to use a novel method to estimate the electron density in the ionized wind, and thus to estimate its mass outflow rate and kinetic power.

Our analysis of spatially-integrated spectroscopy and photometry provides the mass-weighted average properties of the ionized outflows. These estimates are not as accurate as estimates that are based on spatially-resolved IFU observations, but they provide the largest compilation of ionized outflow properties to date, and can be used to study the statistical properties of winds in type II AGN. Below, we discuss our results in the broader context of AGN feedback (section 6.1). We then elaborate on the neutral atomic gas which is associated with the observed ionized winds (section 6.2).

6.1 AGN feedback

AGN feedback, in the form of powerful outflows, is invoked as a way to couple the energy released by the accreting SMBH with the ISM of its host galaxy, thus providing a possible explanation for the observed correlations between the masses of SMBHs and their host bulges (e.g., Silk & Rees 1998; Fabian 1999; King 2003; Di Matteo, Springel & Hernquist 2005; Springel, Di Matteo & Hernquist 2005). A number of such models have successfully reproduced the observed correlations, by requiring that a significant amount of the accretion energy of the AGN will be mechanically-coupled to the ISM of the host galaxy ($\sim$5–10% $L_{\mathrm{bol}}$; e.g., Fabian 1999; Tremaine et al. 2002; Di Matteo, Springel & Hernquist 2005; Springel, Di Matteo & Hernquist 2005; Kurosawa, Proga & Nagamine 2009). According to Hopkins & Elvis (2010), in a two-stage feedback scenario, the initial energy requirement may be 10 times lower.

The coupling efficiency used in hydrodynamic cosmological simulations is often defined as the fraction of the total energy carried out by the nuclear wind close to its launching location. It remains unclear how much of this energy is transferred to the galactic-scale winds, due to various energy losses such as shocks, radiation, and work that is done against the gravitational potential of the host galaxy (e.g., Veilleux et al. 2017; Richings & Faucher-Giguère 2018; see Harrison et al. 2018 for a comprehensive review). Furthermore, galactic-scale winds can be multi-phased, and it is often unclear what is the expected coupling efficiency of the different gas phases. Therefore, the coupling efficiencies required by numerical feedback models are not necessarily similar to kinetic coupling efficiencies derived from observations (Harrison et al. 2018).

The typical mass outflow rate of the objects in our sample is roughly $10^{-2}\,\mathrm{M_{\odot}/yr}$, corresponding to kinetic coupling efficiencies around $\epsilon=\dot{E}_{\mathrm{kin}}/L_{\mathrm{bol}}\sim 10^{-5}$, 3–4 orders of magnitude lower than the typical requirement in most numerical simulations. Our estimates are about 1–2 orders of magnitude lower than typical estimates in AGN host galaxies (e.g., Harrison et al. 2014; Karouzos, Woo & Bae 2016b; Fiore et al. 2017; Harrison et al. 2018), which are based on electron densities in the range $n_{e}\sim 10^{2}-10^{3}\,\mathrm{cm^{-3}}$. Indeed, when we scale the mass outflow rates and coupling efficiencies found by Fiore et al. (2017) to have electron densities of $n_{e}\sim 10^{4.5}\,\mathrm{cm^{-3}}$, we find consistent results for AGN with the same bolometric luminosity. Furthermore, some of the earlier studies are based on IFU observations that target systems which show the most extreme outflow signatures. Our sample is less biased towards the most extreme outflow cases.

In the left panel of figure 9 we show the mass outflow rate versus the AGN bolometric luminosity for the objects in our sample. For comparison, we show the measurements by Fiore et al. (2017) for ionized outflows which, for the same $L_{\mathrm{bol}}$, are much larger. Unlike Fiore et al. (2017), we find no correlation between the mass outflow rate and the AGN bolometric luminosity in our sample. The middle panel shows the kinetic energy of the outflow versus $L_{\mathrm{bol}}$, where we find no correlation as well. In the right panel of figure 9 we show the mass outflow rate versus the star formation rate (SFR) in the host galaxy, where the latter is measured using the Dn4000 index (see details in Paper I). Surprisingly, the strongest ionized outflows in our sample are found in galaxies with high SFRs. This dependence might suggest that star formation is the main driver of the observed winds (see e.g., Wild, Heckman & Charlot 2010; Cicone, Maiolino & Marconi 2016). Since both the narrow and broad kinematic components in all the objects in our sample are classified as AGN using line diagnostic diagrams, it is possible that the process that photoionizes the stationary (NLR) and outflowing gas is not necessarily the process that produces the observed outflows. While the AGN in our sample populate only 3 out of the 4 quadrants indicated in the diagram, which might suggest that there are no active galaxies with low SFRs and strong winds, the small correlation coefficient ($\rho=0.16$) and the large uncertainty on the mass outflow rates, prevent us from drawing more definite conclusions.

The lack of correlation between $\dot{M}_{\mathrm{out}}$ and $L_{\mathrm{bol}}$ might be related to the uncertainty in our mass outflow rate and bolometric luminosity measurements. Combining all the uncertainties discussed in section 5, we find the uncertainty on the mass outflow rate to be roughly 1 dex. The uncertainty on $L_{\mathrm{bol}}$ is 0.3–0.4 dex (Netzer 2009). The combination of limited dynamical range (2 dex in $L_{\mathrm{bol}}$) and the above uncertainties is not enough to explain the complete lack of correlation between $\dot{M}_{\mathrm{out}}$ and $L_{\mathrm{bol}}$.

Finally, the lack of correlation with the AGN bolometric luminosity might be related to projection effects. If the AGN drives the observed winds through radiation pressure, we expect the gas to be accelerated within the AGN ionization cones, which are perpendicular to the line of sight in our type II AGN. Therefore, the measured wind velocity is underestimated by a factor of $\sin\alpha$, where $\alpha$ roughly represents the opening angle. Assuming that the scatter in the observed mass outflow rate is driven by scatter in the opening angle, such that $\log\dot{M}_{\mathrm{obs}}=\log\dot{M}_{\mathrm{outflow}}+\log(\sin\alpha)$, the scatter in the observed mass outflow rate is given by: $\Delta(\log\dot{M}_{\mathrm{obs}})=\sqrt{\Delta(\log\dot{M}_{\mathrm{outflow}})^{2}+\Delta(\log(\sin\alpha))^{2}}$. Under these assumptions, the lack of correlation between $\dot{M}_{\mathrm{out}}$ and $L_{\mathrm{bol}}$ over 2 orders of magnitude implies that the scatter in $\sin\alpha$ is roughly 1-2 dex, depending on the exact shape of the $\sin\alpha$ distribution, and assuming $\log\frac{\dot{M}_{\mathrm{out}}}{1\,M_{\mathrm{\odot}}/\mathrm{yr}}=1.3\log\frac{L_{\mathrm{bol}}}{10^{45}\,\mathrm{erg/sec}}$ (Fiore et al. 2017). We find this very large range in $\sin\alpha$ unlikely.

Given the above considerations, we suggest that the lack of correlation between $\dot{M}_{\mathrm{out}}$ and $L_{\mathrm{bol}}$ is not entirely due to measurement uncertainties and/or projection effects, and is partially driven by the scatter in the intrinsic properties of the winds. This is at odds with the recent results by Fiore et al. (2017), who found such a correlation assuming a constant value electron density in all the objects in their sample. Our new estimates of the electron density for individual sources both decrease the measured mass outflow rates and introduce a scatter. The lack of intrinsic correlation can be attributed to several different factors. First, different outflow histories may result in different spatial extents of the observed winds. Furthermore, since outflows are most likely to go through the path of least resistance, the observed winds may have very different geometries, depending on the mass distribution in their host galaxy.

6.2 Neutral atomic gas in the outflow

So far, we have focused on the warm ionized phase of the outflow, which is traced by strong optical emission lines. However, this gas can have large enough column density to be optically thick to the Lyman continuum radiation. For the ionization parameters found here, this requires a column in the range $N_{\mathrm{H}}=10^{19.5}\,\mathrm{cm^{-2}}$ to $N_{\mathrm{H}}=10^{20.5}\,\mathrm{cm^{-2}}$. The neutral gas at the back of such clouds is largely unconstrained by the observed optical lines, yet it may represent a significant fraction of the outflowing gas mass. The lower-ionization optical lines, such as [OI] $\lambda\lambda$ 6300,6363Å, [NII] $\lambda\lambda$ 6548,6584Å, and [SII] $\lambda\lambda$ 6717,6731Å, peak close to the ionization front. However, their emission drops significantly beyond it, and they cannot be used to place constraints on the fraction of neutral gas in the outflow.

The NaID absorption line is commonly used to constrain neutral gas outflows in active galaxies (e.g., Rupke & Veilleux 2013; Santoro et al. 2018). We examined the NaID absorption profile as a function of different outflow properties in our sample, and found no variation. However, the lack of variation does not necessarily imply that there is no neutral gas in the outflow, since the NaID absorption does not trace the same regions as the optical emission lines.

The dust mid-infrared emission can be used to constrain the neutral fraction of the outflows because the ratio of its infrared luminosity to the dust-corrected H$\alpha$ luminosity, L(IR)/L(H$\alpha$), is sensitive to the ionization parameter and to the total column density (ionized and neutral). Since the dust is mixed with the outflowing gas, L(IR)/L(H$\alpha$) should not depend on the gas covering factor. We find that the observed L(IR)/L(H$\alpha$) ratios in our sample are broadly consistent with those expected from photoionization models, given the empirically-estimated ionization parameters, for hydrogen column densities of $N_{\mathrm{H}}\sim 10^{20.5}-10^{21.5}\,\mathrm{cm^{-2}}$. This suggests significant amount of neutral atomic gas in most cases. It also supports the earlier conclusion (section 4 that most of the outflowing gas is radiation bounded. The inferred hydrogen column density implies that the mass of the neutral gas is a factor of a few larger than the observed ionized phase, with important implications to the mass and energetics of the winds.

Unfortunately, the dust infrared luminosity, L(IR), was estimated in Paper I though SED fitting with significant uncertainties of $\sim$ 0.5 dex. This, and the uncertainty on $U$ mentioned earlier, prevents us from reaching more definite conclusions about the exact amount of neutral atomic gas in individual cases. Furthermore, for a conical outflow, it is possible that the H$\alpha$ emission originating in the receding side of outflow is completely obscured by dust in the ISM of the host galaxy. While the ISM of the host galaxy is optically-thick to optical radiation, it is optically-thin to mid-infrared radiation, and the observed L(IR) traces the full extent of the dust in the outflow. Therefore, the measured L(IR)/L(H$\alpha$) ratio is probably overestimated by a factor of $\sim 2$.

7 Summary and conclusions

In the first paper of this series (Baron & Netzer 2019), we argued that the outflowing gas in active galaxies contains dust. The dust is heated by the central AGN and emits at mid-infrared wavelengths. We analysed the infrared spectral energy distribution of thousands of type II AGN and showed that this dust emission is detected in many systems that host ionized gas outflows. This infrared component offers novel constraints on the outflow properties, in particular, its location. In this work we focused on a subset of 234 galaxies from this sample, in which ionized outflows are detected in many optical emission lines ([OIII], H$\beta$, [OI], [NII], H$\alpha$, and [SII]), and there is a clear detection of a dusty wind component in mid-infrared wavelengths. The combination of the two allowed us to determine key properties of the gas in the wind, such as dust reddening, ionization state, and electron density, and thus constrain the mass and energetics of the winds. While our sample is based on spatially-integrated spectroscopy and photometry, and is thus less accurate than samples that are based on spatially-resolved observations, it forms the largest compilation of ionized outflow properties to date, and can be used to study the statistical properties of winds in type II AGN. Our results can be summarized as follows:

• •

We presented a novel method to estimate the ionization parameter in the outflowing gas. We also provided more general expressions for the H$\alpha$ and [OIII] line emissivities, where the latter depends on the ionization state of the outflowing gas, and thus on the ionization parameter. The new line emissivities, which are then used to estimate the mass and energy of the ionized outflows, are expected to be more accurate than previously-used estimates.

• •

We presented a novel method to estimate the electron density in ionized gas, based on the strong optical line ratios [OIII]/H$\beta$ and [NII]/H$\alpha$, and on the known location of the gas. We applied this method to estimate the electron density in the observed outflows and found electron densities of $n_{\mathrm{e}}\sim 10^{4.5}\,\mathrm{cm^{-3}}$, which are about two orders of magnitude higher than most published estimates. We argue that the commonly-used method to estimate the electron density in the wind, which is based on the [SII] emission lines, underestimates the true densities by 1–2 orders of magnitude. Our estimates are more consistent with recent studies that are based on the density-sensitive auroral and transauroral [SII] $\lambda$4068,4076Å  and [OII] $\lambda$7318,7319,7330,7331Å  emission lines.

• •

We found that the median ionized gas mass in the outflow is roughly $M_{\mathrm{out}}\sim 10^{4}\,\mathrm{M_{\odot}}$, the median mass outflow rate is $\dot{M}_{\mathrm{out}}\sim 10^{-2}\,\mathrm{M_{\odot}/yr}$, and the median kinetic coupling efficiency is roughly $\epsilon=\dot{E}_{\mathrm{kin}}/L_{\mathrm{bol}}\sim 10^{-5}$. These values are 1–2 orders of magnitude lower than typical estimates in type II AGN with similar bolometric luminosity. This is a direct consequence of the high electron densities found in this study compared with the [SII]-based method.

• •

We found no correlation between the mass outflow rate of the ionized outflow and the AGN bolometric luminosity in our sample. This lack of correlation can be attributed to outflows having different locations, geometris, and different fractions of ionized gas in the outflowing clouds. Surprisingly, we found that the strongest outflows in our sample tend to occur in galaxies with high star formation rates, which might suggest that these are driven by supernovae in the host galaxy. Thus, the source that photoionizes the outflowing and stationary gas is not necessarily the source that drives the observed winds.

• •

Our study suggests the existence of a significant amount of neutral atomic gas at the back of the outflowing ionized gas clouds. We made an attempt to estimate the neutral gas fraction by combining the dust infrared luminosity with the H$\alpha$ luminosity. We suggest that the neutral gas mass at the back of the outflowing clouds is a factor of a few larger than the observed ionized gas mass, which has significant implications for the estimated mass and energetics of such flows. The large uncertainties on the dust infrared luminosity prevent us from estimating the neutral gas mass in individual sources.

Acknowledgments

We thank the referee for useful comments and suggestions that helped improve this manuscript. This work is supported in part by the Israel Science Foundation grant 284/13 and by DFG/DIP grant STE 1869/2-1 GE 625/17-1. The spectroscopic analysis was made using IPython (Pérez & Granger 2007). We also used the following Python package: astropy222www.astropy.org/.

This work made use of SDSS-III333www.sdss3.org data. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. SDSS-III is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS-III Collaboration including the University of Arizona, the Brazilian Participation Group, Brookhaven National Laboratory, Carnegie Mellon University, University of Florida, the French Participation Group, the German Participation Group, Harvard University, the Instituto de Astrofisica de Canarias, the Michigan State/Notre Dame/JINA Participation Group, Johns Hopkins University, Lawrence Berkeley National Laboratory, Max Planck Institute for Astrophysics, Max Planck Institute for Extraterrestrial Physics, New Mexico State University, New York University, Ohio State University, Pennsylvania State University, University of Portsmouth, Princeton University, the Spanish Participation Group, University of Tokyo, University of Utah, Vanderbilt University, University of Virginia, University of Washington, and Yale University.

Appendix A Photoionization models

We model the central source using standard assumptions about AGN SEDs. The SED consists of a combination of an optical-UV continuum emitted by an optically-thick geometrically-thin accretion disk, and an additional X-ray power-law source that extends to 50 keV with a photon index of $\Gamma=1.9$. The normalisation of the UV (2500Å) to X-ray (2 keV) flux is defined by $\alpha_{OX}$, which we take to be 1.37.

The grid of models consists of a geometrically-thin, optically-thick, shells of dusty gas, with ISM-type grains and a density of $n_{\mathrm{H}}=10^{4}\,\mathrm{cm^{-3}}$. We used a range of gas metallicities: 0.5, 0.65, 0.8, 1, 1.5, and 2$\mathrm{Z_{\odot}}$, and ionization parameter in the range $\log U=-3.8$ to $\log U=-2$. We consider four different shapes for the ionizing continuum SED, and we list their properties in table 3. Specifically, the SEDs have different slopes, with mean energies of ionizing photon of 2.56, 2.65, 3.15, and 4.17 Ryd. The ranges of metallicities, ionization parameters, and ionizing SEDs were chosen to fully cover the range of possible gas and AGN properties, and thus to provide a general description of the sources in our sample. For the typical AGN in our sample, we use SED 2 (mean energy of an ionizing photon of 2.56 Ryd), where the number of ionizing photons is given by:

[TABLE]

where $L_{\mathrm{bol}}$ is the AGN bolometric luminosity, which can be estimated using narrow emission line luminosities (see Netzer 2009). Using the other SEDs (models 1, 3, and 4) result in a value which is different by roughly 10% from the value given in equation 8. Substituting equation 8 into equation 3 results in the expression given in equation 4 in the main text.

We run a grid of 192 models with version 17.00 of cloudy (Ferland et al. 2017). In figure 10 we show the predicted [OIII]/H$\beta$ versus [NII]/H$\alpha$ by the different models. The colour represents the ionization parameter of the gas, and the size of the markers represents the metallicity, where the smallest markers correspond to 0.5$\mathrm{Z_{\odot}}$, and the largest markers to 2$\mathrm{Z_{\odot}}$. The shapes of the markers represent the ionizing SED, where the squares, stars, circles, and triangles correspond to SEDs with a mean energy of ionizing photon of 2.56, 2.65, 3.15, and 4.17 Ryd respectively. The predictions of these models are in agreement with predictions by other studies (e.g., Groves, Dopita & Sutherland 2004a, b).

The models presented in figure 10 are calculated under the assumption of constant density in the dusty gas. We also examined a similar grid of models, where we assume a constant total pressure (gas and radiation) in the gas. We find that the two sets of models give consistent optical line ratios for a given ionization parameter, and thus our conclusions do not depend on the assumption of constant density versus constant pressure in the gas. This is expected since we model gas clouds with ionization parameters of $\log U\leq-2$, where gas pressure is larger than radiation pressure. Therefore, these constant-pressure models do not represent a case of radiation pressure confinement (RPC), and are effectively constant-density models. RPC takes place when the radiation pressure is significantly larger than the gas pressure (Dopita et al. 2002; Stern, Laor & Baskin 2013), which is achieved for gas clouds with $\log U\geq-2$. The calculations of Stern, Laor & Baskin (2013) show that such models can reproduce the spectrum of Seyfert galaxies with the most extreme [OIII]/H$\beta$ line ratios (e.g., figure 10 in Stern, Laor & Baskin 2013), but cannot describe the bulk of the population. In particular, they do not cover a large region of the BPT diagram where many of our sources reside (figure 2). In such models, a range of ionization parameters (from $\log U\sim 3$ to $\log U\sim-1$) produce similar optical line ratios, and thus the line ratios cannot be used to place constrains on the average ionization parameter in the gas. In the rest of the section, we examine only constant-density models, since constant-pressure models give consistent results for our selected gas properties.

Figure 10 shows that the metallicity and the shape of the ionizing continuum are degenerate. That is, different combinations of gas metallicity and SED shape result in similar [OIII]/H$\beta$ and [NII]/H$\alpha$ line ratios. Therefore, these two properties cannot be determined from the optical line ratios alone, and require additional observables. On the other hand, the colour of the markers varies smoothly across the diagram, suggesting that the ionization parameter can be determined from the [OIII]/H$\beta$ and [NII]/H$\alpha$ line ratios, almost regardless of the gas metallicity or the shape of the ionizing continuum.

The ionization parameter expression derived from our modeling is shown in equation 2. We find the coefficients $a_{1}$, $a_{2}$, $a_{3}$, $a_{4}$, and $a_{5}$, that minimize the residuals between the ionization parameter predicted by the photoionization models and the ionization parameter derived using equation 2. We list the best-fitting coefficients in table 1. In figure 11 we show the predicted [OIII]/H$\beta$ versus [NII]/H$\alpha$ by the different models (same as figure 10), where the background colour represents the ionization parameter that is derived according to equation 2. Overall, there is a good agreement between the ionization parameters predicted by the models and those derived using equation 2. In figure 12 we show the residuals between the two, which span the range of -0.2 to 0.2 dex. The standard deviation of the residuals is $\sigma=0.11$ dex, which we adopt as the uncertainty of the analytic expression.

To test the robustness of our analytical expression, and in particular the dependence of line ratios on the gas density, we considered a more complicated model, consisting of four geometrically-thin, optically thick, gas clouds. All gas clouds have solar metallicity and ISM-type grains. The central source is characterized by a mean energy of an ionizing photon of 2.56 Ryd. We varied the distance of each gas cloud and its hydrogen density, such that all models have the same ionization parameter of $\log U=-2.75$. We set the hydrogen density in the clouds to be $\log n_{\mathrm{H}}(\mathrm{cm^{-3}})=5$, $\log n_{\mathrm{H}}(\mathrm{cm^{-3}})=4$, $\log n_{\mathrm{H}}(\mathrm{cm^{-3}})=3$, and $\log n_{\mathrm{H}}(\mathrm{cm^{-3}})=2$.

As a first test, we examine the emission line ratios [OIII]/H$\beta$ and [NII]/H$\alpha$ of each individual cloud. Since the different clouds were defined to have the same ionization parameter, we expect their line ratios to be similar. The resulting line ratios differ from each other by less than 0.1 dex, supporting our suggestions that the hydrogen density does not affect the observed line ratios, for a constant ionization parameter. We used these line ratios and equation 2 to predict the ionization parameter, and found $\log U$ of -2.68, -2.78, -2.87, -2.96 respectively. These values differ from the input ionization parameter, $\log U=-2.75$, by at most $\sim$0.2 dex, consistent with the residuals shown in figure 12.

Next, we consider the case of a combination of clouds with different densities and distances. We normalize the four models to have a similar covering factor. This is similar to using the mean of the four previous ionization parameters. We extracted the predicted emissivities of [OIII], [NII], H$\alpha$, and H$\beta$ in the four different clouds, and summed them to obtain the total emissivity in each of the lines for the model. We then used the total [OIII]/H$\beta$ and [NII]/H$\alpha$ line ratios, and estimated the ionization parameter using the analytical expression from equation 2. The predicted ionization parameter is $\log U=-2.83$. The difference between the input ionization parameter and the predicted is 0.08 dex, which is smaller than our adopted uncertainty. This suggests that the simple analytical form in equation 2 is expected to hold for more complex scenarios, where the outflow is stratified, with a large range of density and location.