ALMA Multiple-transition Molecular Line Observations of the Ultraluminous Infrared Galaxy IRAS 20551–4250: Different HCN, HCO⁺, and HNC Excitation, and Implications for Infrared Radiative Pumping

Continuum and molecular line emission are dominated by a nuclear compact component; however, in the brightest CO J = 3–2 integrated-intensity (moment 0) map in Figure 3, a spatially extended structure is seen southeast of the nucleus. A similar extended structure is seen in the stellar emission probed in the near-infrared K-band (2.2 μm) image (Duc et al. 1997), suggesting that this extended CO J = 3–2 emission originates in the host galaxy.

In the intensity-weighted mean velocity (moment 1) maps of CO J = 3–2, J = 1–0, 3–2, 4–3, and 8–7 of HCN, HCO⁺, HNC, and CS J = 5–4 in Figures 5–8 and Imanishi et al. (2016b), as well as CO J = 1–0 in Ueda et al. (2014), the overall dynamics are dominated by rotational motion in such a way that the northeastern part is redshifted and the southwestern part is blueshifted, relative to the nucleus. However, in the moment 1 map of the brightest CO J = 3–2 emission line, a dynamically decoupled component from the overall rotation is seen at the southwesternmost region (the clump at the lower-right edge in Figure 5). A plausible explanation is that some type of merger event happened previously, which is quite reasonable given that IRAS 20551−4250 is a ULIRG and ULIRGs are usually driven by gas-rich galaxy mergers (Sanders & Mirabel 1996).

In Figure 4, although we fit emission lines with single Gaussian component, when we examine the line profiles of the very bright CO J = 3–2 and HCO⁺ J = 4–3 emission lines in more detail, skew patterns are recognizable. To investigate the skewed asymmetric line profiles in more detail, we fit these bright emission lines with a Gaussian component, only using data points at the red part of the emission peak (Figures 9(a) and (b)) and at the blue part of the emission peak (Figures 9(c) and (d)). When we fit the red component, the data in the blue part indicates an excess, compared to the best Gaussian fit. Conversely, when we fit with a Gaussian, using data at the blue part of the emission peak only, the extrapolation of the best-fit Gaussian to the redder part is higher than the actual data, particularly for CO J = 3–2. The full-width at half maximum (FWHM) values of the best fit Gaussian are 260 km s⁻¹ and 165 km s⁻¹ for the blue and red components of the CO J = 3–2 emission line, respectively. These are 205 km s⁻¹ and 160 km s⁻¹ for the blue and red components of the HCO⁺ J = 4–3 emission line, respectively. We thus quantitatively confirm that the line widths are larger for the blue components than the red components both for CO J = 3–2 and HCO⁺ J = 4–3. Figure 10 displays the contours of the blue and red components of the CO J = 3–2 and HCO⁺ J = 4–3 emission lines. The red component is slightly offset to the northeastern direction from the blue component; this could be explained by the overall rotational motion of IRAS 20551−4250. A natural interpretation for this skewed profile is that turbulence is stronger at the blueshifted molecular gas-emitting region, which broadens the line width at the blue side. In the CO J = 3–2 moment 1 map in Figure 5, the signature of a merger-induced distinct emission component is seen at the blueshifted southwestern part of the nucleus. This possible merger-induced turbulence component may contribute to the observed broader emission line profile for the blue component.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In Figure 11(a), we show a zoom-in of the bottom part of the very bright CO J = 3–2 emission line, with the best Gaussian fit of this line (Figure 4 and Table 4) overplotted. There are clear excesses at both the blue and red sides of the CO J = 3–2 emission line, although at the red side, possible contamination from the H¹³CN J = 4–3 (${\nu }_{\mathrm{rest}}=345.340$ GHz) emission line makes quantitative discussion of the broad CO J = 3–2 emission component difficult. This kind of profile is typically interpreted to be due to the broad emission component by outflow activity (Feruglio et al. 2010; Alatalo et al. 2011; Aalto et al. 2012; Maiolino et al. 2012; Cicone et al. 2014; Garcia-Burillo et al. 2015). Figure 11(b) displays the contours of the emission at the blue and red parts of the broad emission component defined in Figure 11(a). The peak position of the red broad component is slightly shifted to the southeast, compared to that of the blue broad component (1 pix left and 1 pix lower with the pixel scale of 01 pix⁻¹). Considering the peak positional accuracy of the blue and red broad components with (beam-size)/(S/N), the significance of this positional displacement is marginal. However, this pattern is different from the global rotation of IRAS 20551−4250. We may be witnessing CO J = 3–2 outflow toward (away from) us being ejected in the northwestern (southeastern) direction from the nucleus.

Zoom In Zoom Out Reset image size

In Figure 11(a), we added a Gaussian with FWHM = 540 km s⁻¹ and velocity peak at 12,886 km s⁻¹ for the broad outflow component, and another Gaussian with FWHM = 190 km s⁻¹ to incorporate the H¹³CN J = 4–3 emission line at the red part of the bright CO J = 3–2 emission line. The Gaussian fit fluxes of the broad CO J = 3–2 emission line component and the H¹³CN J = 4–3 emission line at the tail of the very bright CO J = 3–2 emission line are estimated to be ∼7.7 (Jy km s⁻¹) and ∼0.9 (Jy km s⁻¹), respectively.

We here estimate a molecular outflow rate from our CO J = 3–2 emission line data. The peak flux value of the blue broad emission component in Figure 11(b) is ∼0.4 (Jy km s⁻¹). We adopt the value of the blue broad component because the red broad component is likely to be contaminated by the H¹³CN J = 4–3 emission line. Assuming that the outflow-origin red broad component is comparable to the blue broad component, we obtain ∼0.8 (Jy km s⁻¹) for outflow-origin CO J = 3–2 emission. This is a factor of ∼10 smaller than the Gaussian fit flux of the broad CO J = 3–2 emission line component, but we adopt the value of ∼0.8 (Jy km s⁻¹) to obtain the conservative estimate of a molecular outflow rate. We obtain the CO J = 3–2 luminosity with $L{{\prime} }_{\mathrm{CO}J=3-2}\,\sim 7.3\times {10}^{6}$ (K km s⁻¹ pc²) (Solomon & Vanden Bout 2005). Assuming that CO J = 3–2 emission is optically thick and thermalized, and adopting the ULIRG-like CO luminosity to molecular mass (${M}_{{\rm{H}}2}$) conversion factor with $\sim 0.8\,{M}_{\odot }$ (K km s⁻¹ pc²)⁻¹ (Cicone et al. 2014), we obtain a molecular outflow mass of ${M}_{\mathrm{outf}}\sim 5.8\,\times \,{10}^{6}\,{M}_{\odot }$. In Figure 11(b), the peak position difference between the blue broad and red broad emission component is ∼014 or ∼120 (pc). Adopting the outflow peak position offset with R ∼ 60 (pc) from the nucleus and outflow velocity with V = 500 km s⁻¹ (in Figure 11(a), the broad wing component extends to approximately ± 500 km s⁻¹ with respect to the systemic velocity of ${v}_{\mathrm{opt}}={\rm{12,900}}$ km s⁻¹), we obtain a molecular outflow rate with ${\dot{M}}_{\mathrm{outf}}\sim 150$ (${M}_{\odot }$ yr⁻¹), where we adopt the relation of ${\dot{M}}_{\mathrm{outf}}=3\times {M}_{\mathrm{outf}}\times V/R$ (Maiolino et al. 2012; Cicone et al. 2014). Assuming ∼30% AGN contribution to the infrared luminosity of IRAS 20551−4250 (i.e., ${L}_{\mathrm{AGN}}\,\sim 1.1\times {10}^{45}$ erg s⁻¹), the derived molecular outflow rate with ${\dot{M}}_{\mathrm{outf}}\sim 150$ (${M}_{\odot }$ yr⁻¹) agrees within a factor of ∼2 with the relation seen in other ULIRGs (Cicone et al. 2014). The molecular outflow kinetic power is estimated to be ${P}_{\mathrm{outf}}\,\equiv 0.5\times {\dot{M}}_{\mathrm{outf}}\times {V}^{2}\sim 1.2\times {10}^{43}$ (erg s⁻¹), which is ∼1% of the AGN luminosity. The molecular outflow momentum rate is ${\dot{P}}_{\mathrm{outf}}\equiv \,{\dot{M}}_{\mathrm{outf}}$ × V ∼ 4.6 × 10³⁰ (kg m s⁻²), which is ∼12 × L_AGN/c. These values are comparable to those observed in other ULIRGs with detectable molecular outflow activity (Cicone et al. 2014). Note that all of molecular outflow mass (${M}_{\mathrm{outf}}$), molecular outflow rate (${\dot{M}}_{\mathrm{outf}}$), molecular outflow kinetic power (${P}_{\mathrm{outf}}$), and molecular outflow momentum rate (${\dot{P}}_{\mathrm{outf}}$) could increase by an order of magnitude, if we use the Gaussian fit flux of the broad CO J = 3–2 emission line component.

From Table 4 and Imanishi et al. (2016b), we obtained the ratios of HCN-to-H¹³CN J = 3–2 flux in (Jy km s⁻¹) to be ∼22 ± 3, HCO⁺-to-H¹³CO⁺ J = 3–2 flux in (Jy km s⁻¹) to be ∼60 ± 9, and HNC-to-HN¹³C J = 3–2 flux in (Jy km s⁻¹) to be ∼45 ± 15. We must note that the detection significance of HN¹³C J = 3–2 is only ∼3σ in both Gaussian fit in the spectrum and moment 0 map. Thus, discussion of HNC could be more uncertain than HCN and HCO⁺. Adopting the ¹²C-to-¹³C abundance ratios in ULIRGs with 50–100 (Henkel & Mauersberger 1993; Henkel et al. 1993; Martin et al. 2010; Henkel et al. 2014), we find that the flux attenuation by line opacity for HCN J = 3–2 is estimated to be a factor of 2–5, while those of HCO⁺ J = 3–2 and HNC J = 3–2 are a factor of ∼1–1.5 and ∼1–2, respectively, where we assume that isotopologue emission lines are optically thin (Jimenez-Donaire et al. 2017).

In a geometry where radiation sources and molecular gas are spatially well mixed, the flux attenuation and optical depth (τ) are related as $\tfrac{\tau }{1-\exp (-\tau )}$. From this relationship, we obtain τ = 2–5 for HCN J = 3–2. Using the estimated H¹³CN J = 4–3 flux with ∼0.9 (Jy km s⁻¹) (Section 4.1), we obtain the ratio of HCN-to-H¹³CN J = 4–3 flux in (Jy km s⁻¹) to be ∼13. The flux attenuation by line opacity for HCN J = 4–3 is estimated to be 4–7 or τ = 4–7, which is comparable to that for HCN J = 3–2 (τ = 2–5). In summary, our isotopologue observations suggest that the HCN J = 3–2 and J = 4–3 emission lines are considerably flux-attenuated by line opacity.

In contrast, the relatively small flux attenuation with <2 suggests that the line opacity for HCO⁺ J = 3–2 and HNC J = 3–2 is τ < 1.5. Based on the observed HCN, HCO⁺, and HNC J = 3–2 fluxes at v₂ = 1f and v = 0, and calculation of vibrational excitation by infrared radiative pumping, using the available 5–30 μm spectrum of IRAS 20551−4250, Imanishi et al. (2016b) argued that HCN abundance is higher than those of HCO⁺ and HNC. This higher HCN abundance scenario is a plausible explanation for the higher HCN line opacity than HCO⁺ and HNC at J = 3–2.

Based on the fluxes of vibrationally excited (v₂ = 1f) and vibrational ground (v = 0) HCN, HCO⁺, and HNC emission lines at J = 3–2, it has been estimated that the infrared radiative pumping mechanism plays a role in rotational excitation at v = 0 in IRAS 20551−4250, at least for HCN and HNC (Imanishi et al. 2016b). Namely, once HCN and HNC are vibrationally excited to v₂ = 1 by absorbing infrared ∼14 μm and ∼22 μm photons, respectively, they are decayed back to v = 0 and their rotational J-transition fluxes at v = 0 can be higher than those for collisional excitation alone (Rangwala et al. 2011). We here attempt to confirm this result, based on our newly collected ALMA data. The J = 1–0 line is not present at the vibrationally excited level (v₂ = 1). For J = 8–7, v₂ = 1f emission lines are not clearly detected, mainly due to large background noise at this high ALMA frequency range. We thus focus on J = 4–3 data, and follow the same logical steps as those made by Imanishi et al. (2016b) for J = 3–2 data.

Infrared radiative pumping is estimated to play a significant role if the following condition is met:

$$\begin{eqnarray}&&{T}_{\mathrm{ex} \mbox{-} \mathrm{vib}}\gt \displaystyle \frac{{T}_{0}}{\mathrm{ln}\tfrac{{A}_{v2=1-0,\mathrm{vib}}}{{A}_{J=4-3,\mathrm{rot}}}},\end{eqnarray} \tag{ 1 }$$

where ${T}_{\mathrm{ex} \mbox{-} \mathrm{vib}}$ is the v₂ = 1f vibrational excitation temperature, T₀ is the energy level at v₂ = 1f J = 4, ${A}_{v2=1-0,\mathrm{vib}}$ is the Einstein A coefficient from v₂ = 1f to v = 0, and ${A}_{J=4-3,\mathrm{rot}}$ is the Einstein A coefficient from J = 4 to J = 3 at v = 0 (Carrol & Goldsmith 1981; Sakamoto et al. 2010; Mills et al. 2013).

Adopting the values in Table 7, and ${A}_{v2=1-0,{vib}}$ = ∼1.7 s⁻¹, ∼3.0 s⁻¹, and ∼5.2 s⁻¹ for HCN, HCO⁺, and HNC, respectively (Deguchi et al. 1986; Aalto et al. 2007; Mauclaire et al. 1995), we obtain ∼160 (K), ∼185 (K), and ∼210 (K) in the right-hand column of Equation (1) for HCN, HCO⁺, and HNC, respectively. The observed ${T}_{\mathrm{ex} \mbox{-} \mathrm{vib}}$ values in the left-hand column of Equation (1)—derived based on, e.g., Equation (2) of Imanishi et al. (2016b)—are ∼360 (K), <280 (K), and ∼360 (K) for HCN, HCO⁺, and HNC, respectively (Table 8). For HCN and HNC, the condition is fulfilled. If the HNC v₂ = 1f J = 4–3 flux is approximately half that shown in Table 4 (Section 3), the left-hand column of Equation (1) for HNC is ∼260 (K), which is still higher than the right-hand column. For HCO⁺, the upper limit of ${T}_{\mathrm{ex} \mbox{-} \mathrm{vib}}$ still allows the condition to be met.

Table 7.  Parameters for Molecular Transition Lines

Line	Frequency	${E}_{{\rm{u}}}/{k}_{{\rm{B}}}$	${A}_{\mathrm{ul}}$	Flux
	(GHz)	(K)	(10−3 s−1)	(Jy km s−1)
(1)	(2)	(3)	(4)	(5)
HCN J = 1–0	88.632	4.3	0.024	1.4 ± 0.1
HCN J = 3–2	265.886	25.5	0.84	5.9 ± 0.1
HCN J = 4–3	354.505	42.5	2.1	11.6 ± 0.2
HCN J = 8–7	708.877	153.1	17.4	13.7 ± 1.3
HCN J = 3–2, v2 = 1f	267.199	1050.0	0.73	0.25 ± 0.07
HCN J = 4–3, v2 = 1f	356.256	1067.1	1.9	0.60 ± 0.13
HCO+ J = 1–0	89.189	4.3	0.042	2.2 ± 0.1
HCO+ J = 3–2	267.558	25.7	1.5	8.4 ± 0.1
HCO+ J = 4–3	356.734	42.8	3.6	16.9 ± 0.2
HCO+ J = 8–7	713.341	154.1	30.2	12.5 ± 1.0
HCO+ J = 3–2, v2 = 1f	268.689	1217.4	1.3	<0.088
HCO+ J = 4–3, v2 = 1f	358.242	1234.6	3.4	<0.20
HNC J = 1–0	90.664	4.4	0.027	0.37 ± 0.04
HNC J = 3–2	271.981	26.1	0.93	3.2 ± 0.1
HNC J = 4–3	362.630	43.5	2.3	6.6 ± 0.2
HNC J = 8–7	725.107	156.6	19.4	7.2 ± 1.3
HNC J = 3–2, v2 = 1f	273.870	692.0	0.85	0.20 ± 0.07
HNC J = 4–3, v2 = 1f	365.147	709.6	2.2	0.98 ± 0.19

Note. Column (1): Molecular transition line. Column (2): Rest-frame frequency in (GHz). Column (3): Upper energy level in (K). Column (4): Einstein A coefficient for spontaneous emission in (10⁻³ s⁻¹). Values in Columns (3) and (4) are from the Cologne Database of Molecular Spectroscopy (CDMS) (Müller et al. 2005) via Splatalogue (http://www.splatalogue.net). Column (5): Flux in (Jy km s⁻¹) estimated based on Gaussian fit (Table 4, column 9 of this paper; Imanishi & Nakanishi 2013b; Imanishi et al. 2016b). For the undetected HCO⁺ v₂ = 1f J = 3–2 and J = 4–3 emission lines, upper limits based on the integrated intensity (moment 0) maps are adopted.

Download table as:  ASCIITypeset image

Table 8.  Excitation Temperature for IRAS 20551−4250

Molecule	(v2,J; v, J)	${T}_{\mathrm{ex}}$ (K)
(1)	(2)	(3)
HCN	(1f, 4; 0, 4)	360
HCO+	(1f, 4; 0, 4)	<280
HNC	(1f, 4; 0, 4)	360
HCN	(0, 3; 0, 1)	7
	(0, 4; 0, 3)	36
	(0, 8; 0, 4)	42
HCO+	(0, 3; 0, 1)	7
	(0, 4; 0, 3)	38
	(0, 8; 0, 4)	36
HNC	(0, 3; 0, 1)	10
	(0, 4; 0, 3)	41
	(0, 8; 0, 4)	42

Notes. Column (1): Molecule. Column (2): Transition. Column (3): Excitation temperature (${T}_{\mathrm{ex}}$) in (K), derived from observed fluxes. Different beam sizes among different J-transition lines are not considered, in order to be consistent with our previous estimates (Imanishi et al. 2016b).

Download table as:  ASCIITypeset image

The left-hand column of Equation (1) will decrease if flux attenuation by line opacity is significant for some v = 0 lines and the intrinsic line opacity-corrected flux at v = 0 and v₂ = 1f is used for the calculation. From our isotopologue observations, we see that HCN J = 4–3 flux is estimated to be attenuated by line opacity with a factor of 4–7. Adopting a factor of 5 flux attenuation for HCN J = 4–3, the ${T}_{\mathrm{ex}}$ value of HCN decreases only to ∼230 (K). Equation (1) is still fulfilled for HCN, not changing the above conclusion. For HCO⁺ and HNC, flux attenuation by line opacity is not significant at J = 3–2 (Imanishi et al. 2016b). Assuming that this is also the case for J = 4–3, the condition in Equation (1) is still valid for HNC J = 4–3. Thus, it is confirmed from our new ALMA J = 4–3 data that the role of infrared radiative pumping is significant in IRAS 20551−4250, at least for HCN and HNC.

As discussed in Imanishi et al. (2016b), if the column density at v = 0 is the same for HCN, HCO⁺, and HNC, and if all molecular lines are emitted from the same regions, the ratio of the infrared radiative pumping rate among HCN, HCO⁺, and HNC in IRAS 20551−4250 is calculated to be

$$\begin{eqnarray}&&\mathrm{HCN}{\rm{:}}{\mathrm{HCO}}^{+}{\rm{:}}\mathrm{HNC}({P}_{\mathrm{IR}})(\mathrm{predicted})=1{\rm{:}}0.9{\rm{:}}27,\end{eqnarray} \tag{ 2 }$$

and the v₂ = 1 to v = 0 column density ratio among HCN, HCO⁺, and HNC is calculated to be

$$\begin{eqnarray}&&\mathrm{HCN}{\rm{:}}{\mathrm{HCO}}^{+}{\rm{:}}\mathrm{HNC}({N}_{{v}_{2}=1}/{N}_{v=0})(\mathrm{predicted})=1{\rm{:}}0.5{\rm{:}}9,\end{eqnarray} \tag{ 3 }$$

based on the Einstein A coefficients of HCN, HCO⁺, and HNC in Table 7, and the infrared 5–35 μm spectrum of IRAS 20551−4250 obtained with Spitzer IRS (Imanishi et al. 2016b). In essence, the number of infrared photons at ∼22 μm (corresponding to the vibrational-rotational transitions for HNC) is higher than at 12–14 μm (corresponding to the vibrational-rotational transitions of HCN and HCO⁺), such that it is naively expected that HNC is more vibrationally excited by infrared radiative pumping than HCN and HCO⁺.

The ratio of the column density at the v₂ = 1f J = 4 level (${N}_{{v}_{2}=1{\text{}}f,{\text{}}J=4}$) among HCN, HCO⁺, and HNC is derived with ∝$\tfrac{\mathrm{flux}}{{A}_{{ul}}}$, where the line flux from the upper (u) to lower (l) transition level is in units of (Jy km s⁻¹) and A_ul is the Einstein A coefficient for spontaneous emission from the upper (u) to lower (l) level (Goldsmith & Langer 1999; Izumi et al. 2013; Imanishi et al. 2016b). From observational fluxes in Table 4 and Einstein A coefficients in Table 7, we obtain the column density ratio at J = 4 at v₂ = 1

$$\begin{eqnarray}&&\mathrm{HCN}{\rm{:}}{\mathrm{HCO}}^{+}{\rm{:}}\mathrm{HNC}({N}_{{v}_{2}=1f,{\text{}}J=4})(\mathrm{observed})\\ &&\quad =\,1{\rm{:}}\lt 0.2{\rm{:}}0.7,\end{eqnarray} \tag{ 4 }$$

where we adopt that the actual flux of HNC v₂ = 1f J = 3–2 emission line is half of the Gaussian fitted result from Table 4 (see Section 3).

Assuming that the fraction of J = 4 level, relative to all J-levels at v = 0, does not differ significantly among HCN, HCO⁺, and HNC,⁵ and that v₂ = 1f J-transition emission is optically thin (Imanishi et al. 2016b), Equations (3) (prediction) and (4) (observation) can be reconciled if the column density at v = 0 (N_v = 0) (and thereby abundance) for HCN is higher than HCO⁺ and HNC by factors of >2.5 and ∼10, respectively. Thus, we obtain an abundance ratio of

$$\begin{eqnarray}&&\mathrm{HCN}{\rm{:}}{\mathrm{HCO}}^{+}{\rm{:}}\mathrm{HNC}(\mathrm{abundance})=1{\rm{:}}\lt 0.4{\rm{:}}0.1.\end{eqnarray} \tag{ 5 }$$

The HCN-to-HCO⁺ abundance ratio of >2.5 and HCN-to-HNC abundance ratio of ∼10 are comparable to those derived from J = 3–2 data (Imanishi et al. 2016b). The HCN-to-HNC abundance ratios derived in this manner at J = 4–3 and J = 3–2 (∼10) appear to be larger than that from isotopologue observations (∼2). This discrepancy could be partly explained by the possibly large ambiguity of the HN¹³C J = 3–2 flux, which is due to its marginal (∼3σ) detection (Section 4.2). Although a quantitatively accurate estimate of the actual HCN-to-HNC abundance ratio is subject to further investigation, two independent methods consistently suggest that HCN abundance is higher than HCO⁺ and HNC, by at least a few times.

To obtain a better understanding of the physical properties of these line-emitting molecular gases, we use the RADEX software (van der Tak et al. 2007), which treats the non-LTE analysis and predicts molecular line flux ratios at v = 0 by solving the excitation of rotational transitions at v = 0 by collision, based on the large velocity gradient (LVG) approximation. Free parameters are the number density of hydrogen molecule (${n}_{{\rm{H}}2}$), kinetic temperature of H₂ (${T}_{\mathrm{kin}}$), column density of a particular molecule divided by its line width (${N}_{\mathrm{mol}}/{\rm{\Delta }}v$), and background radiation temperature (${T}_{\mathrm{bak}}$), which is set as ∼3 K.

Infrared radiative pumping to v₂ = 1 should play a role for molecular rotational excitation at v = 0, in addition to collision, in IRAS 20551−4250 (Section 4.3). However, with the limited amount of available molecular rotational J-transition line data in our study, separating the contributions from collisional excitation and infrared radiative pumping with sufficient reliability is not easy. Even for collisional excitation alone, it is sometimes argued that (1) HCN and HCO⁺ probe substantially different phases of molecular gas, and (2) there exist multiple temperature components for HCN-, HCO⁺-, and HNC-emitting molecular gas in LIRGs (Greve et al. 2009). These collisional excitation models with increased numbers of parameters and infrared radiative pumping effects are virtually indistinguishable. Vollmer et al. (2017) estimated that, in LIRGs, the infrared radiative pumping has a measurable, but relatively small (20–30%), effect of enhancing the HCN J = 1–0 flux when compared to collisional excitation alone, although the effect for higher-J was not estimated and could be higher. This 20–30% level of flux alternation is comparable to the uncertainty of different J-transition line flux ratios, coming from ALMA's absolute flux calibration uncertainty in individual observations in each band (maximum ∼10%). Thus, we first derive rough values of the molecular gas parameters based on collisional excitation, and will later investigate how they change by including the infrared radiative pumping.

We vary the H₂ number density (${n}_{{\rm{H}}2}$) between 10³ and 10⁸ cm⁻³, because HCN, HCO⁺, and HNC are well-known to be dense (>10³ cm⁻³) molecular gas tracers, due to their high dipole moments and resultant high critical densities. The kinetic temperature range (${T}_{\mathrm{kin}}$) is set from 10 to 500 K, where the lowest kinetic temperature of 10 K corresponds to cool molecular gas far from energy sources, whereas the highest kinetic temperature of 500 K reflects compact (<10 pc), warm molecular gas in the close vicinity of the putative AGN. Because the molecular line widths of HCN, HCO⁺, and HNC at J = 1–0, 3–2, 4–3, and 8–7 are ∼300 km s⁻¹ at 10% of the peak intensity (Imanishi & Nakanishi 2013b; Imanishi et al. 2016b), we adopt Δv = 300 km s⁻¹ for the line width in our RADEX calculation. We include ${T}_{\mathrm{bak}}\sim 3$ K cosmic microwave background radiation in the calculations. The column density of individual molecular gas is not directly derived from observations. The XMM-Newton 2–10 keV X-ray observation of IRAS 20551−4250 shows that the hydrogen column density toward the putative X-ray-emitting obscured AGN is estimated to be ${N}_{{\rm{H}}}={7.9}_{-1.9}^{+6.9}$ × 10²³ cm⁻² (Franceschini et al. 2003). Because the abundance ratios of HCN, HCO⁺, and HNC, relative to H, in warm molecular gas in LIRGs are estimated to be ∼10⁻⁸ (Greve et al. 2009; T. Saito et al. 2017, in preparation), we first set the column density (${N}_{\mathrm{mol}}$) of HCN, HCO⁺, and HNC to be 1 × 10¹⁶ cm⁻².

We summarize our RADEX calculation results in Figure 12, where the ratios are derived from the comparison of observed fluxes in (Jy km s⁻¹) (Table 9). The parameter ranges that can reproduce the observed flux ratios of individual J-transitions are shown as thick lines. For ${n}_{{\rm{H}}2}={\rm{a}}\ \mathrm{few}\times {10}^{5}$ cm⁻³, ${T}_{\mathrm{kin}}=70\mbox{--}100$ K, Δv = 300 km s⁻¹, and HCN column density ${N}_{\mathrm{HCN}}=1\times {10}^{16}$ cm⁻², which can reproduce the observed HCN J = 3–2 to J = 1–0 flux ratio in Figure 12, our RADEX calculations show that the HCN J = 3–2 and J = 4–3 line optical depths are τ ∼ 3, roughly comparable to the observationally derived values from isotopologue line data (Section 4.2). Thus, the choice of ${N}_{\mathrm{HCN}}=1\times {10}^{16}$ cm⁻² is regarded as appropriate.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Because HCO⁺ and HNC abundance are estimated to be at least a few times lower than HCN in IRAS 20551−4250 (Section 4.3), we show the same RADEX calculation results in Appendix (Figure 19) for the HCO⁺ and HNC column density of 3 × 10¹⁵ cm⁻². The overall trends change little.

In interpreting these results, we must consider two possible caveats. First, the line-emitting volume can differ between the high-J and low-J transition lines. AGNs are typically located at galactic center regions where typical molecular gas density is expected to be higher than that in outer off-nuclear regions, due to nuclear molecular gas concentration. Additionally, in an AGN, the surrounding molecular gas and dust has a strong temperature gradient, such that the temperature at the inner part is higher than that at the outer part. Consequently, dense and warm molecular gas is concentrated in the inner part of the nuclear region, from whence high-J transition lines can be efficiently emitted. More diffuse, low-density cool gas is distributed in more spatially extended regions, which can strongly emit low-J transition lines, but high-J transition line emission is relatively limited. Thus, the probed effective molecular gas volume is very likely to be smaller at higher-J than lower-J transition lines, whereas we assume that both lines come from the same volume in the RADEX calculations. If we focus on the volume of high-J line-emitting molecular gas and compare flux ratios within that volume, then the actual high-J to low-J line flux ratios will be higher than the observed ratios. Thus, the derived temperature and density from the observed flux ratios and RADEX calculations are taken as lower limits for the actual values for the high-J line-emitting molecular gas.

Our second caveat is that, because vibrational transition is not included in our RADEX calculations, infrared radiative pumping is not properly taken into account. For IRAS 20551−4250, it has been quantitatively demonstrated by Imanishi et al. (2016b) and this study (Section 4.3) that infrared radiative pumping plays a role in rotational excitation at v = 0, through the decay from v₂ = 1 to v = 0, at least for HCN and HNC, and could realize higher rotational excitation at v = 0 than collisional excitation alone, for a given H₂ number density and kinetic temperature. Owing to infrared radiative pumping, lower values of H₂ number density and kinetic temperature may suffice to reproduce the observed flux ratios. Thus, the H₂ number density and kinetic temperature should be taken as upper limits. This second ambiguity works in opposition to the first, and it is not easy to disentangle the two in a quantitatively reliable manner.

Regarding the first issue, the derived number density and kinetic temperature using low-J lines largely reflect a large volume of spatially extended, low-density, low-temperature molecular gas, resulting in low values. If we compare the flux ratios at high-J transitions only, a small volume of dense and high-temperature molecular gas is preferentially probed, and the derived density and temperature can be high. In Figure 12, we find that the density and kinetic temperature that reproduce the observed flux ratios increase from "J = 3–2 to J = 1–0" to "J = 4–3 to J = 3–2" or "J = 8–7 to J = 4–3" or "J = 8–7 to J = 3–2" for HCN, HCO⁺, and HNC, i.e., the thick curved line moves from the lower-left to upper-right regions. This is a natural consequence—higher-J lines selectively probe dense and warm molecular gas with a small emitting volume.

In Figure 12, the J = 3–2 to J = 1–0 flux comparison is expected to reflect the largest molecular gas volume, including outer low-density and low-temperature gas. For the kinetic temperature of ${T}_{\mathrm{kin}}$ > 25 K, the H₂ number density of ${n}_{{\rm{H}}2}$=10⁵–10⁶ cm⁻³ (10⁴–10⁵ cm⁻³) seems to be required for HCN and HNC (HCO⁺). The lower required density for HCO⁺ is reasonable, given that HCO⁺ has a lower critical density than HCN and HNC under similar line opacity (Meijerink et al. 2007; Greve et al. 2009).

This possible difference among the probed molecular line-emitting volumes can be small if we compare the flux ratios between adjacent J-transitions, such as J = 4–3 and J = 3–2. The typical density and temperature of molecular gas that emits J = 4–3 and J = 3–2 lines are >10⁶ cm⁻³ and >50 K for HCN, HCO⁺, and HNC (Figure 12). We thus quantitatively confirm that these J = 4–3 and J = 3–2 emission lines selectively probe dense and warm molecular gas in the close vicinity of energy sources.

The derived density and kinetic temperature that reproduce the observed "J = 8–7 to J = 3–2" and "J = 8–7 to J = 4–3" flux ratios are lower than the "J = 4–3 to J = 3–2" flux ratios for HCN, HCO⁺, and HNC. This is opposite to the expected trend that the derived density becomes higher when we compare flux ratios at higher-J transition. A plausible explanation for this is that the J = 8–7 emitting volume is substantially smaller than those of J = 4–3 and J = 3–2. This can result in much smaller "J = 8–7 to J = 3–2" and "J = 8–7 to J = 4–3" flux ratios than expected from the same emitting volume assumption, and thereby smaller derived density and kinetic temperature.

The dense molecular mass, derived from the observed HCN J = 1–0 luminosity (Table 6), is estimated to be $\sim 2\times {10}^{9}\,{M}_{\odot }$, where we adopt the relationship between dense molecular mass and HCN J = 1–0 luminosity, ${M}_{\mathrm{dense}}$ = 10×HCN J = 1–0 luminosity [${M}_{\odot }$ (K km s⁻¹ pc²)⁻¹] (Gao & Solomon 2004). Conversely, the hydrogen column density toward the X-ray emitting AGN in IRAS 20551−4250 is ${N}_{{\rm{H}}}\sim 8\times {10}^{23}$ cm⁻² (Franceschini et al. 2003). We here assume that molecular gas distributes spherically and the AGN is located at the center. If all X-ray-obscuring material is in a dense molecular form that is probed with HCN J = 1–0, we obtain the following two relations:

$$\begin{eqnarray}&&r(\mathrm{pc})\times n({\mathrm{cm}}^{-3})\sim 8\times {10}^{23}({\mathrm{cm}}^{-2}),\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&4/3\times \pi \times {r}^{3}({\mathrm{pc}}^{3})\times n({\mathrm{cm}}^{-3})\times {m}_{{\rm{H}}2}\sim 2\times {10}^{9}\,{M}_{\odot },\end{eqnarray} \tag{ 7 }$$

where r (pc) is the radius of the sphere, n (cm⁻³) is volume averaged gas number density, and ${m}_{{\rm{H}}2}$ is the mass of molecular hydrogen (=3.3×10⁻²⁴ g). From Equations (6) and (7), we obtain r ∼ 200 (pc) and n ∼ 10³ (cm⁻³). However, if the kinetic temperature is >30 (K), the typical number density of HCN J = 1–0 emitting gas, derived from the observed HCN J = 3–2 to J = 1–0 flux ratio and RADEX calculation, is 10⁵–10⁶ cm⁻³ (Figure 12). These results can be reconciled if nuclear molecular clouds in IRAS 20551−4250 consist of dense (10⁵–10⁶ cm⁻³) clumps with a small volume filling factor of 0.1%–1%, as argued by Solomon et al. (1987). The slightly low estimated density of 10⁴–10⁵ cm⁻³ from HCO⁺ observations (Figure 12) could be reconciled if HCO⁺ probes more extended lower-density envelopes of each clump (Imanishi et al. 2016b), because HCO⁺ emission is observationally found to be more diffuse than HCN (Imanishi et al. 2007; Saito et al. 2015; T. Saito et al. 2017, in preparation). This kind of clumpy molecular structure can naturally explain similar profiles of various molecular emission lines with different line opacity (Solomon et al. 1987), as observed in IRAS 20551−4250 (Imanishi et al. 2016b).

Figure 13 plots the observed HCN-to-HNC and HCN-to-HCO⁺ flux ratios at J = 1–0, J = 3–2, J = 4–3, and J = 8–7. For J = 3–2 and J = 4–3, line-opacity-corrected intrinsic flux ratios are also shown. In the intrinsic flux ratios, IRAS 20551−4250 deviates even farther away from starburst-dominated regions, suggesting that the line opacity correction, particularly for high-HCN-abundance AGN-containing objects, makes our HCN-based molecular line flux ratio method even more powerful in distinguishing AGNs from starbursts.

Zoom In Zoom Out Reset image size

On the ordinate, the observed HCN-to-HCO⁺ flux ratios are comparable at J = 1–0, J = 3–2, and J = 4–3, but are slightly lower than that at J = 8–7. Because the critical density of HCN is higher than HCO⁺ at individual J-transitions, the higher HCN-to-HCO⁺ flux ratio at higher J-transition is difficult to explain with collisional excitation under the same line opacity. However, the estimated higher line opacity of HCN reduces its effective critical density with 1/τ. Furthermore, due to higher HCN abundance than HCO⁺, line-opacity-corrected intrinsic HCN-to-HCO⁺ flux ratios will be higher than those observed. The intrinsic HCN-to-HCO⁺ flux ratios at J = 3–2 and J = 4–3 are at least higher than the observed ratio at J = 8–7. Because the line opacity at J = 1–0 is estimated to be higher than J = 8–7 with RADEX under the parameter range of ${n}_{{\rm{H}}2}={10}^{5}$ cm⁻³, ${T}_{\mathrm{kin}}=50\mbox{--}100$ K, ${N}_{\mathrm{HCN}}=1\times {10}^{16}$ cm⁻², and Δv = 300 km s⁻¹, the intrinsic HCN-to-HCO⁺ flux ratio at J = 1–0 will move to the upper side more than at J = 8–7. Both of these factors could help to explain the result.

On the abscissa, the observed HCN-to-HNC flux ratio at J = 1–0 is higher than those at J = 3–2, J = 4–3, and J = 8–7. Figure 14 is the plot of the "high-J to low-J" flux ratios for HCN, HCO⁺, and HNC, where the combination of "J = 3–2 to J = 1–0," "J = 4–3 to J = 1–0," "J = 8–7 to J = 1–0," "J = 4–3 to J = 3–2," "J = 8–7 to J = 3–2," and "J = 8–7 to J = 4–3" are selected. The actual ratios are summarized in Table 9. We see a trend that the ratios are a factor of ∼2 higher for HNC than HCN and HCO⁺ for the first three combinations (including J = 1–0 as the low-J transition), but comparable for the latter three combinations (not including J = 1–0). HCN and HNC have comparable critical densities at individual J-transition under a similar line opacity (Greve et al. 2009). The derived higher HCN opacity (Sections 4.2 and 4.3) will make the effective HCN critical density lower than HNC, in which case the "high-J to low-J" flux ratios for HNC are expected to be lower than those for HCN, when collisional excitation is considered. This is opposite to the trends observed in Figure 14.

Zoom In Zoom Out Reset image size

The excitation temperatures of rotational J-transition at v = 0 are summarized in Table 8, where the differences of beam sizes are not considered, to be consistent with our previous calculations (Imanishi et al. 2016b). This assumes that the intrinsic emission size is the same among HCN, HCO⁺, and HNC, and is smaller than the beam sizes in all observations. The rotational excitation temperatures within v = 0 are slightly higher for HNC than for HCN (Table 8).

We suggest that the following two mechanisms could help reduce the discrepancy between observations and predictions. First, we find from our RADEX calculations that, at high density (${n}_{{\rm{H}}2}\gt 3\,\times \,{10}^{5}$ cm⁻³) and temperature (${T}_{\mathrm{kin}}=50\mbox{--}100$ K), a maser phenomenon occurs at J = 1–0 first for HCN and HCO⁺, and then later for HNC. For example, for the column density of 1 × 10¹⁶ cm⁻², Δv = 300 km s⁻¹, and ${n}_{{\rm{H}}2}$ = 1–3 × 10⁶ cm⁻³, the maser phenomenon is predicted to occur at ${T}_{\mathrm{kin}}\gt 30$ K for HCN J = 1–0, but only at ${T}_{\mathrm{kin}}\gt 70$ K for HNC J = 1–0. Further, with ${n}_{{\rm{H}}2}=3\,\times \,{10}^{5}$ cm⁻³, HCN J = 1–0 shows the maser phenomenon at ${T}_{\mathrm{kin}}=50\mbox{--}100$ K, but HNC J = 1–0 does not. The maser phenomenon could boost J = 1–0 emission. If the maser occurs only for HCN J = 1–0, but not for HNC J = 1–0, under some parameter range, the "high-J to J = 1–0" flux ratios for HCN could be smaller than HNC. HCO⁺ J = 1–0 also starts to show the maser phenomenon earlier than HNC from lower density and temperature. The higher values of "high-J to J = 1–0" flux ratios for HNC than HCO⁺ could also be explained by the maser phenomenon in some parameter range.

Our second scenario is infrared radiative pumping. In IRAS 20551−4250, it is calculated that the infrared radiative pumping rate is higher for HNC than for HCN and HCO⁺ (Section 4.3 and Imanishi et al. 2016b). This could boost high-J HNC emission more than HCN and HCO⁺, through decay back from v₂ = 1 to v = 0. Given that this infrared radiative pumping is estimated to play a role in IRAS 20551−4250, at least for HCN and HNC (Section 4.3), its effect must be evaluated in an appropriate manner.

Our ultimate goal is to quantitatively estimate (1) what fraction of molecular excitation at individual J-transition levels is due to collision and infrared radiative pumping, and (2) the increase of molecular line fluxes at individual J-transitions, brought by infrared radiative pumping, when compared to collisional excitation alone. However, this is difficult to achieve with currently available data, because of a large freedom of parameters for collisional excitation. First, the possibly different emitting volumes at different J-transitions (Section 4.4) are not constrained in a meaningful manner. Next, the typical temperature and density probed by different J-transition lines are different (Section 4.4), so different numbers of multiple temperature and density components need to be taken into account when we compare flux ratios among different J-transition lines. We thus investigate quantitatively how infrared radiative pumping reproduces the observed trend of higher values of "high-J to low-J" flux ratios for HNC than HCN shown in Figure 14, by fixing molecular gas physical parameters to canonical values. We use MOLPOP (Elitzur & Asensio Ramos 2006), which handles the exact solution of radiative transfer problem in multi-level atomic systems. The ∼3 K cosmic microwave background radiation is always included in MOLPOP. The LVG method is prepared as an option, and we adopt it to compare the output with our RADEX calculation results. We have quantitatively confirmed that, under LVG, both MOLPOP and RADEX outputs about molecular line flux ratios among J = 1–0, J = 3–2, J = 4–3, and J = 8–7 agree within 20–30% for various H₂ number density, H₂ kinetic temperature, and molecular column density divided by line width.

To calculate the role of infrared radiative pumping, we require Einstein A coefficient values for vibrational (v)-rotational (J) transitions between (v₂, J) = (1, J) and (0, J ± 1 or J). We find the required information for HCN and HNC, but not for HCO⁺, in the GEISA database (Jacquinet-Husson et al. 2016). We thus perform this calculation for HCN and HNC. Collision rates between J-transition levels at v₂ = 1 are not available, so we tentatively adopt the same values as v = 0. Because vibrational-rotational transitions are much faster than the transitions between J-levels within v₂ = 1 (Imanishi et al. 2016b), possible small ambiguities in the collision rates at v₂ = 1 do not significantly affect our final results. Collisional transitions between v₂ = 1 and v = 0 are also not taken into account, due to their very large energy gaps.

In MOLPOP, we first set ${n}_{{\rm{H}}2}={10}^{5}$ cm⁻³, ${T}_{\mathrm{kin}}=100$(K), Δv = 300 km s⁻¹, and the HCN and HNC column density of 10¹⁶ cm⁻², in order to calculate molecular line fluxes by collisional excitation. We confirm that the maser phenomenon does not happen at J = 1–0. We choose the small side of the estimated density range for HCN and HNC, ${n}_{{\rm{H}}2}={10}^{5}$ cm⁻³ for our calculations, because the maser phenomenon happens at higher ${n}_{{\rm{H}}2}$ values, which complicates interpretations. We then add infrared radiation at 5–30 μm in the power-law form of ${F}_{\nu }$ (Jy) = 0.0024 × λ(μm)^2.0 (the dashed straight line in Imanishi et al. (2016b)). In this form, the number of infrared photons at ∼22 μm is larger than at 12–14 μm. The integrated luminosity between 5–30 μm is estimated to be $\sim 2\times {10}^{11}\,{L}_{\odot }$ at a distance of 187 Mpc (z = 0.043). Plane-parallel geometry is adopted, and the infrared 5–30 μm radiation source is placed at 30 pc and 100 pc from the molecular gas. We then investigate if the high-J to low-J flux ratios of HNC are higher than HCN, particularly for "J = 3–2 to J = 1–0," "J = 4–3 to J = 1–0," and "J = 8–7 to J = 1–0," as observed in Figure 14. Figure 15(a) plots the values of "high-J to low-J flux ratio of HNC, divided by the same ratio of HCN" in the ordinate. Because the infrared radiative pumping effects are calculated to be higher for HNC than HCN, the high-J to low-J flux ratios for HNC are expected to increase more than those for HCN, by the inclusion of infrared radiation. Namely, the values in the ordinate are expected to be higher in the case of stronger infrared radiation than the case with no infrared radiation (i.e., collisional excitation only). We find that the inclusion of infrared radiative pumping better reproduces the observed trend of higher values of "high-J to low-J flux ratios" for HNC than HCN, for the combinations of "J = 3–2 to J = 1–0," "J = 4–3 to J = 1–0," and "J = 8–7 to J = 1–0" observed in Figure 14.

Zoom In Zoom Out Reset image size

In Figure 15(b), we display a similar plot, but with a factor of ∼3 smaller HNC column density than HCN, as indicated in Section 4.3. Because an additional parameter is introduced, the behavior becomes complicated, but the general trend of enhanced high-J to J = 1–0 flux ratios for HNC, as compared to HCN, is observed. We thus argue that the observed trend in IRAS 20551−4250 can be largely accounted for by infrared radiative pumping.

Although collisional excitation is usually adopted to interpret the observed molecular J-transition emission line fluxes in galaxies, it has been widely argued that the infrared radiative pumping mechanism can play an important role in AGN-containing galaxies, because AGNs can emit strong mid-infrared 5–30 μm continuum photons due to AGN-heated hot dust grains. In this AGN-containing ULIRG, IRAS 20551−4250, we have demonstrated that the infrared radiative pumping indeed plays a role for molecular rotational J-level excitation, at least for HCN and HNC, which suggests that this is the case in other AGN-containing ULIRGs as well. However, in general ULIRGs, the similar estimates of the infrared radiative pumping effects, through the v₂ = 1f emission line flux measurements, are prohibitively difficult because HCN J = 3–2 and J = 4–3 emission lines at the v₂ = 1 level are spectrally highly overlapped with the nearby, much brighter HCO⁺ J = 3–2 and J = 4–3 emission lines at the v = 0 level, respectively (Aalto et al. 2015b; Imanishi et al. 2016c). IRAS 20551−4250 is unique in that, thanks to its small observed molecular emission line widths, the flux estimates of the HCN v₂ = 1f J = 3–2 and J = 4–3 emission lines are possible, by clearly separating from the bright HCO⁺ (v = 0) emission lines. Further detailed investigations of IRAS 20551−4250 with additional emission lines at v = 0 and v₂ = 1 will help better constrain the quantitative role of infrared radiative pumping in this ULIRG as well as other AGN-containing ULIRGs.

Figure 16 shows the spectral energy distribution (SED) of IRAS 20551−4250 in the infrared and radio. The q-value, defined as the far-infrared-to-radio flux ratio, is q ∼ 2.7 (Table 1). This is not smaller than for starburst-dominated galaxies and many radio-quiet AGNs (q ∼ 2.3–2.4) (Condon et al. 1991; Barvainis et al. 1996; Crawford et al. 1996; Roy et al. 1998), suggesting that the infrared-to-radio SED of IRAS 20551−4250 is an example of a radio-quiet AGN and starburst composite.

Zoom In Zoom Out Reset image size

Continuum data points in ALMA bands 7 (275–373 GHz) and 9 (602–720 GHz) are well-fitted to extrapolation from the infrared dust thermal radiation at 2000–7500 GHz (40–150 μm). ALMA band 3 (84–116 GHz) continuum data points roughly agree with the extrapolation from lower-frequency synchrotron emission. The thermal free–free emission from starburst regions in IRAS 20551−4250 is also present (Imanishi & Nakanishi 2013b; Imanishi et al. 2016b), but mainly contributes to the frequency range at the minimum flux density (∼200 GHz), due to its relatively flat spectral shape. Overall, continuum flux measurements with our ALMA data in the (sub)millimeter wavelength range are consistent with those at other wavelengths taken with other observing facilities, further confirming that most activity (AGN and/or starburst) in IRAS 20551−4250 is spatially compact, and their emission is well-recovered with our ALMA interferometric data.