Circumnuclear Multiphase Gas in the Circinus Galaxy. II. The Molecular and Atomic Obscuring Structures Revealed with ALMA

Circinus is a spiral galaxy (SAb; de Vaucouleurs et al. 1991) with a high inclination angle on the galactic scale (∼65°; Freeman et al. 1977), and hosts the nearest (D = 4.2 Mpc, 1'' = 20.4 pc; Tully et al. 2009) type 2 Seyfert nucleus (Oliva et al. 1994; Moorwood et al. 1996). It shows a one-sided ionization cone and outflowing gas that extends up to kiloparsec scales in Hα, [O iii], and several other coronal emission lines (e.g., Marconi et al. 1994; Veilleux & Bland-Hawthorn 1997; Müller Sánchez et al. 2006). The genuine existence of an obscured type 1 nucleus was confirmed by the spectropolarimetric detection of a broad Hα line (≃3300 km s⁻¹; Oliva et al. 1998). Circinus is also undergoing modest star-forming activity, as measured using far-infrared luminosity (star formation rate (SFR) ∼a few ${M}_{\odot }$ yr⁻¹; Moorwood & Glass 1984; Elmouttie et al. 1998b) and several hydrogen recombination lines (e.g., Marconi et al. 1994; Wilson et al. 2000), while the ≲100 pc scale SFR is admittedly low (∼0.1 M_⊙ yr⁻¹; Esquej et al. 2014): thus, significant SN feedback might not be expected there.

As is typical in the case of a late-type galaxy, a huge amount of molecular gas (≳10^8–9 M_⊙) was detected in Circinus through observations of multi-J CO lines including isotopologues, as well as other dense gas tracers such as HCN and HCO⁺ (e.g., Aalto et al. 1991; Israel 1992; Curran et al. 1998; Elmouttie et al. 1998b; Curran et al. 1999, 2001, 2008; Hitschfeld et al. 2008; Zhang et al. 2014; Israel et al. 2015). While the intense H i emission is distributed in a ∼10 kpc radius region dominated by large-scale spiral arms (Jones et al. 1999), molecular gas is concentrated in the inner ∼1 kpc radius region (Curran et al. 1998, 2008; Elmouttie et al. 1998b). An expected morphology for that molecular gas structure, with the H₂ mass (${M}_{{{\rm{H}}}_{2}}$)⁹ of ∼4 × 10⁸ M_⊙ ($r\lt 440$ pc; Elmouttie et al. 1998b), has been a ring/disk-like one plus widespread outflows perpendicular to the morphological major axis of this galaxy (Curran et al. 1998, 1999). However, recent ∼2''–3'' resolution CO(1–0) mapping with ALMA imaged instead a more spiral arm-like gas distribution (Zschaechner et al. 2016), while revealing another molecular outflow at 35'' northwest of the nucleus.

With regard to the nuclear scale, X-ray spectra below 10 keV exhibit flat continuum and prominent 6.4 keV Fe Kα line emission, indicative of strong Compton scattering (Matt et al. 1996). Subsequent harder X-ray observations at >30 keV confirmed the existence of the Compton-thick nucleus (e.g., Guainazzi et al. 1999; Matt et al. 1999), with a line-of-sight obscuring column density of ${N}_{{\rm{H}}}=(6-10)\times {10}^{24}$ cm⁻², and an absorption-corrected 2–10 keV intrinsic luminosity of L_2–10keV = (2.3–5.1) × 10⁴² erg s⁻¹, respectively (Arévalo et al. 2014). Detections of H₂O mega-masers at both millimeter and submillimeter wavelengths at the heart of Circinus support the idea that a dense Keplerian disk is located there (Greenhill et al. 1997, 2003; Hagiwara et al. 2013), with the mass of the central supermassive black hole (SMBH) given as ${M}_{\mathrm{BH}}=(1.7\pm 0.3)\,\times {10}^{6}\,{M}_{\odot }$ (Greenhill et al. 2003). With this M_BH, Tristram et al. (2007) estimated the Eddington ratio of the AGN as ∼0.2.

Our prime reason for studying Circinus is that it clearly exhibits parsec-scale polar elongation at MIR continuum emission (Tristram et al. 2014), which makes it an ideal laboratory to test our multiphase dynamic torus model. We describe the ALMA Cycle 4 observations in Section 2. The continuum emission maps are shown in Section 3, while line emission distributions and their ratio are reported in Section 4. Details of the gas kinematics of both the CO(3–2) and the [C i](1–0)  lines are presented in Section 5. We compare the observed torus properties with the predictions of our multiphase dynamic torus model in Section 6. Finally, our conclusions are summarized in Section 7.

We observed Circinus with ALMA Band 7 and 8 during 2016–2017 using 42–47 antennas, as a Cycle 4 program (ID = 2016.1.01613.S, PI = T. Izumi). Table 1 summarizes the log of our observations. Observations were conducted in a single pointing with fields of view of 18'' (Band 7) and 13'' (Band 8), which fully covered the central ∼2'' of the CND (see Section 4). The expected maximum recoverable scales per pointing are ∼7'' (Band 7) and ∼5'' (Band 8). We set the phase tracking center to (${\alpha }_{{\rm{J}}2000.0}$, ${\delta }_{{\rm{J}}2000.0}$) = (14^h13^m09950, −65°20'210), which is one of the nuclear 22 GHz H₂O maser spots of Circinus (Greenhill et al. 2003). The angular separations between this tracking center and the phase calibrators are ∼3°. In both Band 7 and 8 observations, one of the four spectral windows (each with a width of 1.875 GHz) was used to fully cover the CO(3–2) or [C i](1–0) emission lines, both in the 2SB dual-polarization mode.

Data reduction, calibration, and analyses were performed with CASA version 4.7 (McMullin et al. 2007) in the usual way. The line and underlying continuum emission were reconstructed using the CLEAN task with Briggs weighting (robust = 0.5). The velocity spacings of the original data were 3.4 km s⁻¹ (Band 7) and 2.4 km s⁻¹ (Band 8) per channel, but several channels were binned to improve the signal-to-noise ratio, which resulted in a final and common velocity resolution of ∼10 km s⁻¹. Note that velocities are expressed in the optical convention in this work (local standard of rest [LSR] frame).

The achieved synthesized beams and rms sensitivities are listed in Table 2. The rms values in the line cubes were measured at channels free of line emission (i.e., thermal noise), whereas those of the continuum maps were measured in areas free of such emission. These continuum emission were subtracted in the uv plane before making the line cubes. Throughout this paper, the pixel scale of Band 7 maps is set to 003, whereas that of Band 8 maps is 01, and emission with <1.5σ are clipped to enhance the clarity of the images. The kinematic position angle (PA) of the major axis is defined to be on the receding half of the galaxy, taken anti-clockwise from the north direction on the sky. Given this definition, we added 180° to some morphological PAs reported in previous works to maintain consistency. The absolute flux uncertainty is ∼10% according to the ALMA Cycle 4 Proposer's Guide, but the displayed errors indicate only statistical ones unless otherwise mentioned. Some parts of our analyses were also performed with the MIRIAD package (Sault et al. 1995).

We retrieved the Hubble Space Telescope (HST) Planetary Camera images (PC-F656N, F606W, F814W) of Circinus from the Hubble Legacy Archive.¹⁰ The images were calibrated using the HST pipeline. Continuum emission was subtracted from the F656N image using the adjacent continuum following standard line extraction procedures to construct an image of the Hα emission line. The structures seen in the resultant image are consistent with those found in previous works (e.g., Marconi et al. 1994; Wilson et al. 2000). We also use the Ks-band image obtained with the NaCo/VLT (Mezcua et al. 2016) to compare the spatial distributions of the continuum emission.

Figure 2 shows the spatial distributions of the ALMA Band 7 (ν_rest = 351 GHz or λ_rest = 854 μm) and 8 (ν_rest = 485 GHz or λ_rest = 618 μm) continuum emission at the central 10'' × 10'' (∼200 pc×200 pc) and 2'' × 2'' (∼40 pc × 40 pc) regions of Circinus. The total fluxes are ∼66.5 mJy (Band 7) and ∼320 mJy (Band 8) within the 10'' box. Each of the emission peaks at α_ICRS = 14^h13^m09948, and δ_ICRS = −65°20'2105, which coincides with the AGN position identified by the very long baseline interferometry (VLBI) H₂O maser observations (α_J2000.0 = 14^h13^m09953, δ_J2000.0 = −65°20'21187; Greenhill et al. 2003) within positional uncertainties. Hereafter, we define our continuum peak position as the AGN location. The peak flux densities are 22.4 mJy beam⁻¹ (Band 7) and 87.8 mJy beam⁻¹ (Band 8), respectively.

Zoom In Zoom Out Reset image size

Although the achieved angular resolutions were admittedly different (Table 2), consistent spatial distributions can be identified between the Band 7 and 8 maps. The continuum-emitting regions appear to consist of (likely) three spiral arms which converge to the central bright CND or jointly constitute the CND: these structures are better echoed by the CO(3–2) and [C i](1–0) spatial distribution maps (see also Figure 6 as a reference). Note that we were unable to find a clear indication of the nuclear gaseous bar postulated by Maiolino et al. (2000) southeast, nor an expected counter bar northwest, of the nucleus. This may imply that an aligned configuration of dense molecular clouds in the spiral arms (Section 4), viewed with a high inclination, might have been mistaken as a bar.

The nature of the observed submillimeter continuum emission is of interest because the obscuring structure of this AGN would begin to be directly traced at the high resolutions obtained here. Although this galaxy possesses a radio jet (Harnett et al. 1990; Elmouttie et al. 1998a; Murphy et al. 2010), we first argue that contamination by such synchrotron emission is not significant at the ALMA Band 7 and 8 frequencies. For example, a Band 7 flux density at the central 14 region expected by extrapolating the λ = 6 cm emission (50 mJy) with a typical synchrotron spectral index¹¹ of α = −0.7 is only ∼2.5 mJy, which is ∼9 times smaller than the value measured at a much smaller beam (029 × 024) placed at the AGN position. Note that Elmouttie et al. (1998a) revealed a very flat radio spectral index at the nuclear region of Circinus at ${\nu }_{\mathrm{rest}}\lesssim 8\,\,\mathrm{GHz}$ primarily due to the optically thick nature at that frequency range. However, it is more reasonable to adopt the above α = −0.7 for higher frequencies, where synchrotron emission becomes optically thinner.

A Band 8 to 7 flux density ratio at the AGN position, after matching the Band 7 uv range and beam size to the Band 8 ones (uv range ≃25–470 kλ; the resultant Band 7 flux density = 26.3 mJy beam⁻¹), revealed a steep spectral index of 3.72, or dust emissivity index¹² of β = 1.72. This is consistent with a typical β index observed in local star-forming galaxies (∼1.8; Clemens et al. 2013). Thus, the nuclear submillimeter SED of Circinus is dominated by thermal dust emission. We found that the Band 8 continuum flux density at the AGN position is almost consistent with the value expected by the ${T}_{\mathrm{dust}}\sim 300$ K blackbody model in Tristram et al. (2014), which was introduced to describe the observed IR SED at the (circum)nuclear region. Meanwhile, the Band 7 flux density at the same position falls by a factor of a few below the model prediction, which seems to reflect the fact that the dust continuum emission tends to follow the modified blackbody spectrum more, rather than the blackbody model, at a longer wavelength.

It is remarkable that the extended spiral arm-like or filamentary structure seen around the CND in our submillimeter continuum maps traces the region that is radially farther from the center than the nuclear bright V-shaped structure seen in the Ks/F814 flux ratio map (Figure 3: the latter seems to collimate the ionized cone seen in the Hα line emission distribution) that enhances the distribution of warm dust emission (Mezcua et al. 2016). This indicates the existence of temperature-dependent dusty structures around this AGN. This relative geometry, viewed with a high inclination angle (∼65°–75°, Freeman et al. 1977; Elmouttie et al. 1998b; Tristram et al. 2014), is also consistent with the prediction of our multiphase dynamic torus model where the warm dust is embedded in outflows at the inner part of the disk, and the cold dust is located in the extended disk (Schartmann et al. 2014; Wada et al. 2016, see also Figure 1). Within the framework of our model, the nuclear warm dust seen in the Ks/F814 ratio map is thus considered to be distributed in a geometrically thick volume, whereas the cold dust is in a thinner structure. This warm dust carries only ∼1/3 of the obscuring matter required to make a type 1 Seyfert nucleus appear as an obscured type 2 nucleus such as Circinus (Mezcua et al. 2016). Thus, an additional absorber is required to completely obscure the nucleus: we suggest that the extended cold material in the CND studied here contributes to that obscuration (see also Section 4).¹³

Zoom In Zoom Out Reset image size

To further investigate the detailed spatial distribution of the submillimeter continuum emission, we applied the UVMULTIFIT procedure (Martí-Vidal et al. 2014) to the Band 7 continuum visibility data, in which we have sampled long uv distance visibilities denser than the case of the Band 8 data. This model fitting method to the direct interferometric observable (visibility) is preferable to avoid systematic uncertainty in the image plane fitting due to, e.g., a nonlinear deconvolution algorithm, particularly when we fit beam-unresolved components.

The achieved visibility data is shown in Figure 4(a) as a function of uv distance (UVD), after averaging the visibilities over the observation time and frequency (representative frequency = 350.2 GHz). The almost constant amplitude at the UVD ≳ 600 kλ manifests the existence of a compact component, while the short UVD components reflect spatially extended structures. Here, we performed a two-dimensional double Gaussian fit to this visibility data: the use of double components is motivated by the existence of a parsec-scale MIR polar elongation (Tristram et al. 2007, 2014), in addition to the extended CND seen in Figure 2. The results of the fit are summarized in Table 3 and the corresponding model is displayed in Figure 4(a) as well. As expected in the above, our best fit indeed supports the existence of a compact component at the AGN position (Component-1 in Table 3), in addition to the apparently extended structure (Component-2). The latter would correspond to a bright part of the CND (see also Figure 2). The fact that Component-2 is brighter than Component-1 is consistent with the prediction of our multiphase dynamic torus model that cold dust visible at longer FIR to submillimeter wavelengths is predominantly located in an extended disk, while more centrally concentrated components including warm dusty outflows are prominent at shorter NIR to MIR wavelengths (Schartmann et al. 2014). Note that, however, we also found that the residual amplitude remains high at the short UVD (≲400 kλ; Figure 4(b)). This is due to the fact that there are some other spatially extended circumnuclear dusty structures in Circinus (e.g., spiral arms; Figure 2). As we would like to focus on the very nuclear structures in this subsection, fitting all components is beyond our scope.

Zoom In Zoom Out Reset image size

The spatial distributions of several (circum)nuclear major components that we are interested in are illustrated in Figure 5. The axis ratio of Component-2 suggests that this structure is moderately inclined if it is a circular disk (≃55°). This inclination angle is slightly smaller than the galactic-scale one (≃65°) possibly because the spatially extended/elongated emission along the northeast–southwest direction, particularly at the southwestern side of the nucleus, is faint (see also Figure 2(a)). This will decrease the inclination angle with our simple method. Such an extended structure is clearly recognized in, e.g., the CO(3–2) integrated intensity map (Section 4), which indeed yields a better agreeing inclination angle with the above-mentioned galactic-scale one. Meanwhile, Component-1 appears elongated along the polar direction of the H₂O maser disk. This immediately brings to mind the parsec-scale MIR polar elongation (Tristram et al. 2014). The derived source size (1.6 ± 0.3 pc ×1.2 ± 0.1 pc) and PA ($298^\circ \pm 21^\circ$) of Component-1 are roughly consistent with those of the MIR polar elongation (∼1.9 pc and ∼287°, respectively). Therefore, we suggest that Component-1 traces the same physical structure as the MIR polar elongation. Although we need higher-resolution direct imaging to robustly reveal this compact structure, our finding here will support the actual existence of such elongated structures at the hearts of AGNs (e.g., López-Gonzaga et al. 2016), which created a new challenge to the classic torus paradigm.

Zoom In Zoom Out Reset image size

Figure 6(a) shows the CO(3–2) velocity-integrated intensity map of Circinus, over the central 25'' (=500 pc) region. The structures seen here are the counterparts to the postulated 300 pc-radius ring or disk in Elmouttie et al. (1998b). After convolving the map to an 18'' aperture and correcting for the primary beam attenuation, we found that a roughly ∼50% of the total CO(3–2) line flux has been resolved out compared to the APEX single-dish data (Zhang et al. 2014). As for the dust continuum emission, we identify two main structures in this circumnuclear region through closer inspection of the map: these are illustrated in Figure 6(b). A zoomed-in view (central 10'' = 200 pc) is shown in Figure 7(a). We describe each structure below.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

4.1.1. Spiral Arms

4.1.2. The CND

The [C i](1–0)  line emission is successfully detected both at the spiral arms and the CND (Figure 8(a)), with a clear peak at the AGN position (154 Jy beam⁻¹ km s⁻¹). The [C i](1–0)  distribution is consistent with that of the simultaneously obtained Band 8 continuum emission (Figure 8(b)). This [C i](1–0)  map shows the global distribution of cold and low-density (critical density ${n}_{\mathrm{crit}}=5\times {10}^{2}$ cm⁻³, Kaufman et al. 1999) gas in Circinus. It is difficult to estimate the missing flux for this [C i](1–0)  observation, as the field of view of ALMA Band 8 (128) is even smaller than the mapped area (18'') with APEX reported in Zhang et al. (2014). However, we found that our observation recovered 45% of that single-dish flux after correcting for the primary beam attenuation.

Zoom In Zoom Out Reset image size

The overall resemblance of the [C i](1–0)  distribution and the CO(3–2) distribution (Figure 9) supports the previous argument that [C i](1–0)  essentially traces the same area as low-J CO lines, which has previously been reported for Galactic molecular clouds (e.g., Ikeda et al. 1999; Plume et al. 2000; Oka et al. 2005) and extragalactic objects (Krips et al. 2016; Miyamoto et al. 2018), and has been found in numerical models (e.g., Glover et al. 2015). However, we also found slight spatial inconsistencies in the detailed gas distributions, both at the spiral arms and the CND (Figure 9). Therefore, the [C i](1–0)  also traces different gas volume from what the CO(3–2) does.

Zoom In Zoom Out Reset image size

With regard to the spiral arms, we suggest that the spatial offsets (≲10 pc) of the [C i](1–0)  peaks and the CO(3–2) peaks are due primarily to excitation conditions, as CO(3–2) requires ∼100× denser gas than [C i](1–0) to be excited. Miyamoto et al. (2018) also reported spatial offsets between [C i](1–0)-bright knots and CO(3–2)-bright knots in the circumnuclear starburst ring of NGC 613. Meanwhile, our particular interest in this paper is centered on the [C i](1–0) to CO(3–2) spatial offset inside the CND. The [C i](1–0)  peaks exactly at the AGN location, while the CO(3–2) peak is ∼1'' northeast of that position (Figure 9); this is discussed in more detail in Section 4.3.

We also performed a first-order estimate of ${M}_{{{\rm{H}}}_{2}}$ from the [C i](1–0)  line luminosity (${L}_{\mathrm{CI}(1-0)}^{{\prime} }$) at the CND, defined by the beam-deconvolved source size (362 × 166) from the CO(3–2) analyses (Section 4.1). The obtained value there is 1.42 × 10⁶ K km s⁻¹ pc². We first estimate the atomic carbon mass (${M}_{\mathrm{CI}}$) by following Ikeda et al. (2002) and Weißet al. (2005) as,

$$\begin{eqnarray}&&{M}_{\mathrm{CI}}=5.71\times {10}^{-4}Q({T}_{\mathrm{ex}})\displaystyle \frac{1}{3}{e}^{23.6/{T}_{\mathrm{ex}}}{L}_{\mathrm{CI}(1-0)}^{{\prime} }\,[{M}_{\odot }],\end{eqnarray} \tag{ 1 }$$

where $Q({T}_{\mathrm{ex}})=1+3{e}^{-{T}_{1}/{T}_{\mathrm{ex}}}+5{e}^{-{T}_{2}/{T}_{\mathrm{ex}}}$ is the C i partition function, ${T}_{1}=23.6$ K and ${T}_{2}=62.5$ K are the level energies above the ground state, and ${T}_{\mathrm{ex}}$ is the excitation temperature. This equation assumes optically thin line emission under the LTE condition.

Even under the LTE condition, we would need to correct for line opacity, as there is no guarantee that the [C i](1–0)  line emission in Circinus is optically thin. However, we expect that its line opacity is moderate (${\tau }_{\mathrm{CI}(1-0)}\lesssim 1$) given its peak line intensity (∼1144 mJy beam⁻¹ or ∼12.3 K) at the AGN position and its likely high ${T}_{\mathrm{ex}}$ such as ≳30 K (typical value observed in high-redshift quasars, Walter et al. 2011). If ${\tau }_{\mathrm{CI}(1-0)}\lesssim 1$, then the correction factor for flux attenuation (${\tau }_{\mathrm{CI}(1-0)}$/$1-\exp (-{\tau }_{\mathrm{CI}(1-0)})$) is ≲1.6, i.e., the Equation (1) gives a good estimate of ${M}_{\mathrm{CI}}$. Then we achieve ${M}_{\mathrm{CI}}\simeq 1770\,{M}_{\odot }$ within the CND by assuming that ${T}_{\mathrm{ex}}=30$ K; this mass does not change significantly as long as ${T}_{\mathrm{ex}}\gtrsim 30$ K. If we adopt a relative C i abundance with respect to H₂ as 8 × 10⁻⁵, which is also the typical value observed in high-redshift quasars (Walter et al. 2011), we obtain ${M}_{{{\rm{H}}}_{2}}$ of ≃3.7 × 10⁶ M_⊙ in the CND. This is consistent with the CO(3–2)-based ${M}_{{{\rm{H}}}_{2}}$ derived in Section 4.1.2.

Here, we discuss the [C i](1–0)/CO(3–2) integrated intensity ratio ($\equiv {R}_{\mathrm{CI}/\mathrm{CO}};$ brightness temperature scale) at both the CND and the spiral arms by comparing them with previous measurements in external galaxies and Galactic star-forming regions. We first adjust the uv range and the beam size of the CO(3–2) data to those of the [C i](1–0)  data (25–470 kλ, 071 × 066 with PA = 96°). The resultant CO(3–2) cube has a 1σ sensitivity of 0.54 mJy beam⁻¹ at a velocity resolution of 10 km s⁻¹.

Figure 10 shows the R_CI/CO in the central region of Circinus. The typical value of R_CI/CO at the spiral arms is ∼0.3, although there are two spots with high ratios; we expect that the local gas densities there are too low to efficiently excite the CO(3–2) line. On the other hand, it is notable that the R_CI/CO at the AGN position is ∼3 times higher (∼0.9) than the spiral arms. Given the much higher ${n}_{\mathrm{crit}}$ of CO(3–2) than that of [C i](1–0), this ratio would be hard to explain by simple gas excitation, as it would be naively expected that the gas density profile is centrally peaked. The same situation is also demonstrated in the line spectra (Figure 11): while the [C i](1–0)  and the CO(3–2) show totally comparable line flux densities when measured at the central 10'' box (i.e., CND + spiral arms), the [C i](1–0)  emission clearly becomes brighter than the CO(3–2) emission at the AGN position. The relevant R_CI/CO measured at several spatial scales are listed in Table 4; the high R_CI/CO at the AGN position is uncommon in star-forming galaxies, which could be explained by an efficient CO dissociation due to hard X-ray irradiation from the AGN (XDR, Maloney et al. 1996; Meijerink & Spaans 2005).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Indeed, based on the one-zone XDR model of Maloney et al. (1996), we can estimate the X-ray energy deposition rate per particle (H_X/n) with

$$\begin{eqnarray}&&{H}_{X}\sim 7\times {10}^{-22}{L}_{44}{r}_{2}^{-2}{N}_{22}^{-1}\,\mathrm{erg}\,{{\rm{s}}}^{-1},\end{eqnarray} \tag{ 2 }$$

where L₄₄ is the 1–100 keV X-ray luminosity in units of 10⁴⁴ erg s⁻¹, r₂ is the distance from the AGN to the point of interest in units of 100 pc, and N₂₂ is the X-ray attenuating column density in units of 10²² cm⁻², respectively. Here, we adopted L₄₄ = 0.13 (with a photon index =1.31 and cut-off energy = 160 keV) from Arévalo et al. (2014), as well as ${n}_{{{\rm{H}}}_{2}}={10}^{5}$ cm⁻³ as a typical value in the molecular part of CNDs of Seyfert galaxies (e.g., Izumi et al. 2013; Viti et al. 2014). Then, we estimated H_X/n at r = 7 pc from the center (i.e., the half-major axis of the [C i](1–0)  synthesized beam) as $\mathrm{log}({H}_{X}/n)=-27.3$. According to the models in Maloney et al. (1996), this is almost exactly the rate at which the fractional abundance of C i (${X}_{\mathrm{CI}}$) becomes equal to that of CO (${X}_{\mathrm{CO}}$). At the regions where $\mathrm{log}({H}_{X}/n)\gt -27.3$, we can expect ${X}_{\mathrm{CI}}/{X}_{\mathrm{CO}}\gt 1$. Thus, we suggest that the molecular dissociation due to X-ray irradiation is significant at the close vicinity of the Circinus AGN, which could cause the high R_CI/CO there.

However, the actual three-dimensional (3D) geometrical structure of the XDR requires careful consideration: it does not necessarily form a simple one-zone structure. Indeed, single Gaussian fits to the observed CO(3–2) and [C i](1–0)  spectra measured at the AGN position (Figure 11(b) and Table 5) revealed that the [C i](1–0)  spectrum has a wider full width at half maximum (FWHM) than the CO(3–2) spectrum, as well as deviated components from the single Gaussian profile at about ±70 km s⁻¹ to the systemic velocity (446 km s⁻¹). Both of these imply that the [C i](1–0)  emission also traces a different gas volume than the CO(3–2) emission. This is consistent with the multiphase dynamic torus model (Section 5 and 6), where the complex interplay of gas chemistry and dynamics indeed occurs. Furthermore, there could be several other possibilities that can cause high R_CI/CO, including severe self-absorption in the CO(3–2) line and very peculiar gas excitation conditions. Observations of multiple transition lines and multiple species, including isotopologues, are required to better understand the origin of this enhancement.

Finally, we compare the R_CI/CO of Circinus measured with three apertures with those of other galaxies (compiled from Zhang et al. 2014) in Figure 12. The comparison data consist of high-z quasar-host galaxies and submillimeter galaxies (Walter et al. 2011), as well as nearby galaxies, including NGC 6946 and M83 (Israel & Baas 2001), IC 342 (Israel & Baas 2003), Henize 2-10 and NGC 253 (Bayet et al. 2004), and M51 (Israel et al. 2006). As stated in Zhang et al. (2014), the R_CI/CO of the nearby star-forming galaxies are ∼0.1 to 0.2, while two AGN-host galaxies, M51 and Circinus (18'' aperture), show slightly higher values of ∼0.3, which already implies a sort of AGN influence on the line ratio. The averaged R_CI/CO of the high-z objects is 0.32 ± 0.13 (Walter et al. 2011). As these are essentially starburst galaxies (e.g., Casey et al. 2014), one plausible explanation for the higher R_CI/CO compared to nearby galaxies is that their mean interstellar radiation field over the entire galaxy scale is as high as those of the central ∼1 kpc of the nearby AGNs, which would also lead to efficient CO dissociation. Compared to those nearby and high-z samples, the increasing trend of R_CI/CO toward the center of Circinus, which would imply XDR chemistry, is remarkable.

Zoom In Zoom Out Reset image size

First, we show the global gas dynamics (intensity-weighted line-of-sight velocity and velocity dispersion) probed by both the CO(3–2) and the [C i](1–0)  emission lines in Figure 13. These maps were made with the CASA task IMMOMENTS with 10σ clipping to avoid noisy pixels. The patterns traced by the two lines are reasonably consistent with each other, although we would need to beam-deconvolve them to achieve the actual gas dynamics (i.e., intrinsic rotation and dispersion; Section 5.2). The gas motion is clearly dominated by the galactic rotation with an overall northeast–southwest orientation (Figures 13(a), (d)), which is consistent with previous studies at larger spatial scales (e.g., Curran et al. 2008; Zschaechner et al. 2016). On the other hand, it is also evident in the spiral arms that streaming motions are superposed on the rotation pattern (Figures 13(b), (e)). These streaming motions, as well as spatial structures in each velocity channel (see the Appendix), may help the reader to recognize the likely three spiral arms that we have postulated and illustrated in Figure 6(b).

Zoom In Zoom Out Reset image size

To extract basic beam-deconvolved dynamical information, particularly rotational velocity (V_rot) and dispersion (σ_disp), we fitted concentric titled rings to the observed velocity structures with the ^3D Barolo code (Di Teodoro & Fraternali 2015). The main parameters were dynamical center, systemic velocity (V_sys), V_rot, σ_disp, galactic inclination (i), and PA, all of which can be varied in each ring. Here, we fixed the dynamical center to the AGN position and the V_sys to 446 km s⁻¹ based on the single Gaussian fits to the observed nuclear spectra (Table 5); V_rot, σ_disp, i, and PA were thus the major free parameters and were evaluated by minimizing $|$ model—observed data $|$. Given the significantly different spatial resolution between the CO(3–2) cube and the [C i](1–0)  cube, we first present the results using the full-resolution CO(3–2) cube in Section 5.2.1 to describe the details of the gas dynamics as much as possible, and then present those with the uv- and beam-matched CO(3–2) and [C i](1–0)  cubes in Section 5.2.2.

5.2.1. The CO(3–2) Full-resolution Data

The full-resolution CO(3–2) cube (θ = 029 × 024, dV = 10 km s⁻¹) was used. We modeled 50 concentric rings with a separation of Δr = 015, which is roughly half the size of the major axis of the synthesized beam. The modeled mean velocity field (moment 1) and the residual image after subtracting this model from the observed data are displayed in Figures 14(a) and (b), respectively. The residuals at the spiral arms are ∼±40 km s⁻¹ at arms 1 and 2, which are attributed to streaming motions. As the southeastern part is the near side of this galaxy (Wilson et al. 2000), these motions are likely inflows toward the vicinity of this AGN. Meanwhile, the residual is only a few km s⁻¹ at the AGN position, which is consistent with our torus model tuned for Circinus (Wada et al. 2016), where the inflow velocity through the dense mid-plane of the CND is as slow as ≲10 km s⁻¹ (see also Figure 2 of Wada 2012).

Zoom In Zoom Out Reset image size

Figures 15(a) and (b) present the modeled position–velocity diagrams (PVDs) along the global kinematic major and minor axes, respectively, overlaid on the observed PVDs. The overall structures are well reproduced by a combination of gas rotation and dispersion, although streaming motions are evident as well along the minor axis (offset ∼2'' and ∼4''), which positionally correspond to the spiral arms.

Zoom In Zoom Out Reset image size

Figures 16(a), (b), and (c) show the radial profiles of the decomposed V_rot and σ_disp, i, and PA of our model, respectively. We found that variation in the i and the PA is small (within ∼15°) around $i\sim 65^\circ$ and PA ∼216 These values are consistent with those on the galactic scale (e.g., Freeman et al. 1977; Elmouttie et al. 1998b). Nevertheless, we found a slight increasing trend in the i toward the AGN position. The inferred i ≳ 70° at r ≲ 10 pc would be consistent with the $i\gt 75^\circ$ estimated for the 1 pc scale nuclear MIR disk (Tristram et al. 2014), as well as with the $i\sim 75^\circ$ required to reproduce the nuclear IR SED of this AGN based on our torus model (Wada et al. 2016). This increasing i will eventually reach ∼90° as Circinus hosts the 22 GHz H₂O maser disk at the center (Greenhill et al. 2003).

Zoom In Zoom Out Reset image size

The radial profile of the σ_disp/V_rot ratio is shown in Figure 16(d). This ratio can be considered as a scale height (H) to scale radius (R) ratio (H/R, i.e., the aspect ratio of a disk), under the hydrostatic equilibrium condition. From this panel, one can see that the σ_disp/V_rot at r ≳ 15 pc is ∼0.25, while the ratio at the central r ≲ 10 pc increases to ∼0.4 (except for the innermost r = 3 pc ring), i.e., the dense molecular disk becomes geometrically thicker. Within the qualitative framework of our model, this moderate thickness of the dense gas disk at the central ≲10 pc can be interpreted such that the dense gas and dust are puffed up due to the turbulence induced by the AGN-driven failed winds (Wada 2012; Wada et al. 2016, see Figure 1).

It is particularly noteworthy that the innermost r = 3 pc ring shows a very low ${\sigma }_{\mathrm{disp}}/{V}_{\mathrm{rot}}\sim 0.1$, which is not observed at the other radii. The implied very thin disk geometry, as well as the high i, would suggest that the geometrically thin Keplerian disk has finally started to be captured (i.e., the gas motion is governed by the gravity of the central SMBH) with this cold molecular emission line. Indeed, in the zoomed-in view of the radial profile of V_rot (Figure 16(e)), we found a turnover trend of V_rot at r = 3 pc. Although the V_rot at that radius remains comparable to those at larger radii within uncertainties, it is close to the extension of the nuclear Keplerian rotation curve around the central SMBH (${M}_{\mathrm{BH}}=1.7\times {10}^{6}\,{M}_{\odot }$), traced by the VLBI 22 GHz H₂O maser observations (Greenhill et al. 2003). Although it is difficult at this moment to achieve a firmer conclusion given the large errors associated with the inner rings, we can observe a clear difference in V_rot between e.g., r = 1 pc and r = 3 pc if the CO gas dynamics genuinely follow the Keplerian pattern. This can be tested with higher-resolution ALMA observations.

5.2.2. The uv- and Beam-matched Data

Next, we model the gas dynamics using the [C i](1–0)  cube and the uv- and beam-matched CO(3–2) cube (1σ = 0.54 mJy beam⁻¹ at dV = 10 km s⁻¹). With this treatment, we expect to reduce the potential systematic effects on the resultant dynamical properties stemming from the unmatched uv ranges and beam sizes of the input data. Given the larger angular resolution (071 × 066), we modeled only 20 concentric rings with Δr = 025 and repeated the procedure in Section 5.2.1.

Figure 17(a) displays the resultant radial profiles of V_rot and σ_disp traced by the CO(3–2) and the [C i](1–0)  emission lines. We found that both the V_rot and the σ_disp of the CO(3–2) derived here are consistent with those from the full-resolution CO(3–2) cube in Section 5.2.1. Thus, the ${\sigma }_{\mathrm{disp}}/{V}_{\mathrm{rot}}$ profile of the CO(3–2) shown in Figure 17(b) is essentially the same as that in Figure 16(d). Furthermore, the V_rot and the σ_disp of the [C i](1–0)  are consistent with those of the CO(3–2) derived here within uncertainties at most radii, except for the σ_disp at the innermost one (025 = 5 pc). Considering our model, the similarity in the ${\sigma }_{\mathrm{disp}}/{V}_{\mathrm{rot}}$ ratios at $r\sim 10-20$ pc between the two tracers (Figure 17(b)) in turn suggests that the turbulence induced by the failed winds, which can act both on the atomic and molecular gas, plays a major role in determining the disk geometry there, i.e., the contribution from the outflows preferentially seen at ionized/atomic gas (Wada et al. 2016) would not be very significant. The commonly geometrically thin nature at r ≳ 20 pc, where the radiation-driven fountain will not work well for the case of Circinus (Wada et al. 2016), is also reasonable because SN feedback (an efficient mechanism to puff up the disk) should be weak given the low SFR (∼0.1 M_⊙ yr⁻¹; Esquej et al. 2014).

Zoom In Zoom Out Reset image size

Meanwhile, the discrepancy of the ${\sigma }_{\mathrm{disp}}/{V}_{\mathrm{rot}}$ ratios at r ≲ 10 pc (∼0.15 in the CO(3–2) and ∼0.4 in the [C i](1–0)) is more evident, which indicates the existence of multiphase dynamical structures around the AGN. Note that it is difficult to distinguish coherent outflow motion that can widen the line profile from isotropic turbulent motion due to failed winds, based simply on the observed σ_disp that mixes these motions. However, as the outflows are preferentially observed at ionized and/or atomic gas at the very nuclear region in our model (Wada et al. 2016), we expect that such outflows are the prime source of the geometrically thicker nature of the C i disk than the CO disk at this very central region. Indeed, we found that adding two more components (i.e., triple Gaussians) can fit the observed nuclear [C i](1–0)  spectrum more faithfully, as shown in Figure 18 (see also Table 6). The central Gaussian component has almost the same centroid velocity and FWHM as the CO(3–2) spectrum (Table 5), suggesting that this traces the same region as the CO(3–2) line (the rotating disk mid-plane in our model). We also found that the blueshifted and redshifted components appear almost symmetrically with respect to the central component, with a velocity offset of ∼75 km s⁻¹ (blue) and ∼60 km s⁻¹ (red), respectively; this is suggestive of coherent atomic gas motion around the AGN. Therefore, these results support the view that the diffuse atomic gas is distributed in a geometrically thicker volume than the dense molecular gas around the AGN due to atomic outflows, as expected in our multiphase dynamic torus model. As this argument is rather qualitative, however, we perform a more quantitative comparison between the observed [C i](1–0)  and CO(3–2) line properties and our model predictions in the next section.

Zoom In Zoom Out Reset image size

Table 6.  Results of the Three-component Gaussian Fit to the Observed [C i](1–0)  Spectrum

	Peak	Centroid	FWHM
	(mJy beam−1)	(km s−1)	(km s−1)
Central	1107.9 ± 1.6	453.7 ± 5.6	97.7 ± 0.6
Blueshifted	466.4 ± 5.6	379.1 ± 0.2	56.7 ± 0.5
Redshifted	365.9 ± 4.4	512.3 ± 0.1	28.3 ± 0.4

Note. These fits were performed for the spectra at the AGN position, measured with the single synthesized beam (071 × 066). The results of the single Gaussian fit are summarized in Table 5.

Download table as:  ASCIITypeset image

To interpret the physical/dynamical nature of the circumnuclear obscuring structures from the observed line cubes, we used a snapshot of our multiphase dynamic torus model tuned for Circinus (Wada et al. 2016), where non-equilibrium XDR chemistry (Meijerink & Spaans 2005) was solved, to predict the emission line properties of multiple species.

The simulated CND has a circular gas distribution for simplicity with $r\simeq 16$ pc, with the total ${M}_{{{\rm{H}}}_{2}}$ of $2\times {10}^{6}\,{M}_{\odot }$ (128³ grids; 0.25 pc resolution), which is roughly consistent with the observed ${M}_{{{\rm{H}}}_{2}}$ in Circinus measured with a slightly larger aperture (∼3–4 × 10⁶ M_⊙ in the $r\sim 35$ pc disk with α_CO = 0.8 M_⊙ (K km s⁻¹ pc²)⁻¹; Section 4). The central AGN has the same M_BH and Eddington ratio as those observed (Section 1). The ultraviolet radiation from the central thin accretion disk is angle-dependent, whereas we assumed spherically symmetric X-ray radiation from the AGN corona. Self-gravity of the gas is ignored, as it is not essential for the dynamics in the fountain scheme (see also Wada 2012). After the XDR chemistry (abundance) was solved, post-processed non-LTE radiative transfer calculations were performed. Once the radiation field and optical depth are determined in each grid, we can observe line emission from an arbitrary direction (for further details see Wada et al. 2018). For example, Wada et al. (2018) presented the resultant properties of multiple ¹²CO transition lines from J = 0 to 15. Their procedure was fully applied to the case of the [C i](1–0)  in this work.

Figure 19 shows the edge-on distributions of the simulated integrated intensities of CO(3–2) and [C i](1–0)  emission lines. This reinforces our argument that the diffuse atomic gas traced by the [C i](1–0)  is more extended along the z-direction of the disk than the dense molecular gas traced by the CO(3–2) at this spatial scale (r ≲ 10 pc), which is due to the atomic outflows in our model. On the other hand, one can recognize that the most intense part of the [C i](1–0)  distribution essentially traces the same region as the CO(3–2) distribution. Note that the hydrogen (H₂ + H i) volume density in the disk mid-plane (∼10⁵ cm⁻³) is ≳100 times higher than that of the spatially/vertically extended region (≲10³ cm⁻³, Wada et al. 2016). Thus, the line-emitting region can be divided into the following two parts:

• 1.  
  $| {\rm{\Delta }}z| \lesssim 3$ pc: atomic gas and molecular gas coexist. Turbulence due to the failed winds determines the geometrical thickness of this region. The observed CO(3–2) line profile, as well as the central Gaussian component of the [C i](1–0)  line profile (Figure 18), reflect this component.
• 2.  
  $| {\rm{\Delta }}z| \gtrsim 3$ pc: preferentially seen in the [C i](1–0)  emission, not in the CO(3–2) emission. The AGN-driven outflows, which are selectively seen in the ionized and atomic gas phase, are responsible for producing this geometrically thick structure in the [C i](1–0)  distribution.

Zoom In Zoom Out Reset image size

Figure 20 shows simulated channel maps of (a) CO(3–2) intensity and (b) [C i](1–0)  intensity at the systemic velocity of this galaxy with i = 75° (the best-fit value to reproduce the IR SED of Circinus with our model; Wada et al. 2016). For both the CO(3–2) and the [C i](1–0), the bulk of the emission comes from the region outlined by the CO(3–2) distribution. This region corresponds to the geometrically thin mid-plane of the CND, if viewed from the edge-on angle (Figure 19); therefore, we call this the inclined thin disk. Owing to the high inclination adopted here, we are still able to see the spatially extended (or vertically elongated along the z-direction due to outflows) component beside the inclined thin disk in the [C i](1–0)  map (Figure 20(b)). As the outflows would have coherent velocities, we expect that blueshifted and redshifted components with such velocities can be superposed on a simple Gaussian profile in the simulated [C i](1–0)  line profile, which is confirmed in the analysis below.

Zoom In Zoom Out Reset image size

Figure 21(a) therefore exhibits how the [C i](1–0)  and CO(3–2) line profiles appear when viewed with i = 75° based on our model (see also Table 7 for the results of the Gaussian fits). The simulated CO(3–2) profile that basically reflects the rotating mid-plane component of the inclined thin disk can be explained by a single Gaussian, as expected, whereas deviation from a single Gaussian is evident in the case of the [C i](1–0) profile. The [C i](1–0)  outflows are then clearly manifested in the three components Gaussian fit to the spectrum (Figure 21(b)). The blueshifted and redshifted features emerge at around the coherent outflow velocity of ∼50–60 km s⁻¹ in our model. Hence, we here suggest that (i) this qualitative difference between the CO(3–2) and the [C i](1–0)  line profiles, i.e., the existence of the additional two components to the single Gaussian profile in the case of the [C i](1–0), as well as (ii) the centroid velocities of these outflow components, are consistent with the observed profiles (Figure 11(b) and 18). Note that the amplitudes of the model components are smaller than, and the relative strengths among the three [C i](1–0)  components are somewhat different from, the observed ones (see Tables 6 and 7). These discrepancies can be reconciled by adjusting some input parameters of the model (e.g., AGN luminosity, M_H2), which are generally associated with considerable uncertainties. Lastly, we found that the outflow velocities of the [C i](1–0)  implied by our observations (∼60 and ∼75 km s⁻¹; Figure 18) will not significantly exceed the escape velocity (v_esc) from the vicinity of this AGN.¹⁵ Thus, a significant fraction of the gas carried by those winds will eventually fall back to the disk as failed winds. Therefore, given these consistencies, we conclude that our results support the view that a radiation-driven fountain including outflows indeed functions in the central region of Circinus to form a geometrically thick structure, which would explain the physical origin of its torus.

Zoom In Zoom Out Reset image size