EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 655 (November 2021) A&A, 655 (2021) A25 Full HTML
Free Access
Issue
A&A
Volume 655, November 2021
Article Number A25
Number of page(s) 13
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202039773
Published online 05 November 2021

© ESO 2021

1. Introduction

Understanding the relation between the feedback from accreting supermassive black holes (SMBHs) and the interstellar medium (ISM) of their host galaxies is still an open issue. Host bulge properties such as velocity dispersion, luminosity, and mass are tightly correlated with the mass of the SMBHs in the galactic centre, as shown by past seminal works (Gebhardt et al. 2000; Ferrarese & Ford 2005; Kormendy & Ho 2013; Shankar et al. 2016, 2017). The gas in the galactic bulge inflows towards the nucleus and gives rise to the growth of a central SMBH through the active luminous phase of the active galactic nucleus (AGN) (Fabian 2012; King & Pounds 2015). During their active phases, AGN can generate winds that interact with galaxy ISM, potentially altering both the star formation process and the nuclear gas accretion. When the black hole reaches a critical mass the AGN driven winds, the nuclear activity, and the SMBH growth are stopped giving rise to the SMBH-host bulge properties relations (Silk & Rees 1998; Fabian 1999; King 2003).

Outflows are ubiquitous in both luminous AGN and in local Seyfert galaxies, and occur on a wide range of physical scales, from highly ionised semi-relativistic winds and jets in the nuclear region at subparsec scales to galactic scale outflows seen in mildly ionised, molecular, and neutral gas (Morganti et al. 2016; Fiore et al. 2017; Fluetsch et al. 2019; Lutz et al. 2020; Veilleux et al. 2020, and references therein). In some cases molecular and ionised winds have similar velocities and are nearly co-spatial, suggesting a cooling sequence scenario where molecular gas forms from the cooling of the gas in the ionised wind (Richings & Faucher-Giguere 2017; Menci et al. 2019). Other AGN show ionised winds that are faster than the molecular winds, suggesting a different origin of the two phases (Veilleux et al. 2020, and references therein). The molecular phase is a crucial element of the feeding and feedback cycle of AGN because it constitutes the bulk of the total gas mass and it is the site of star formation activity. On galactic scales massive molecular winds are common in local Seyfert galaxies (e.g. Feruglio et al. 2010; Cicone et al. 2014; Dasyra et al. 2014; Morganti et al. 2015; García-Burillo et al. 2014, 2017, 2019); these winds likely suppress star formation (i.e. negative feedback) as they reduce the molecular gas reservoir by heating or expelling gas from the host-galaxy ISM. In late-type AGN-host galaxies, the gas kinematics appears complex at all scales, showing several components such as bars, rings, and (warped) discs, with high velocity dispersion regions (e.g. Shimizu et al. 2019; Feruglio et al. 2020; Fernández-Ontiveros et al. 2020; Alonso-Herrero et al. 2020; Aalto et al. 2020; Audibert et al. 2020). Accurate dynamical modelling of the molecular gas kinematics reveals kinematically decoupled nuclear structures, high velocity dispersion at nuclei, trailing spirals, and evidence of inflows and AGN-driven outflows. (e.g. Combes et al. 2019; Combes 2019, 2021). The outflow driving mechanism (wind shock, radiation pressure, or jet), their multiphase nature, and their relative weights and impact on the galaxy ISM are still open problems (Faucher-Giguère & Quataert 2012; Zubovas & King 2012; Richings & Faucher-Giguere 2017; Menci et al. 2019; Ishibashi et al. 2021). To date, far different outflow phases have been observed only for a handful of sources. Atomic, cold, and warm molecular outflows have been observed in radio galaxies (e.g. Morganti et al. 2007; Dasyra & Combes 2012; Dasyra et al. 2014; Tadhunter et al. 2014; Oosterloo et al. 2017). The nuclear semi-relativistic phase and the galaxy-scale molecular phase have been observed simultaneously in less than a dozen sources, with varied results: in some cases data suggest energy driven flows (Feruglio et al. 2015; Tombesi et al. 2015; Longinotti et al. 2018; Smith et al. 2019), in other cases data suggest momentum driven flows (e.g. García-Burillo et al. 2014; Feruglio et al. 2017; Fluetsch et al. 2019; Bischetti et al. 2019; Marasco et al. 2020). Fiore et al. (2017), using a compilation of local and high redshift winds, showed that there is a broad distribution of the momentum boost, suggesting that both energy- and momentum-conserving expansion may occur. Enlarging the sample of local AGN-host galaxies with outflows detected in different gas phases is important to understand the nature and driving mechanisms of galaxy-scale outflows.

In this paper we present a study of Markarian 509 (Mrk 509) based on ALMA observations. Mrk 509 is a Seyfert 1.5 galaxy located at D ∼ 142.9 Mpc (z = 0.034397, Huchra et al. 1993) implying a physical scale of 0.689 kpc arcsec−1, and hosts an AGN with L2 − 10 keV = 1044.16 erg s−1 (Shinozaki et al. 2006), a bolometric luminosity of Lbol = 1044.99 erg s−1 (Duras et al. 2020), and a black hole with mass Mbh = 108.04 M⊙ (Bentz et al. 2009). The AGN is hosted by a medium-size galaxy resembling a bulge in optical imaging (Ho & Kim 2014), but hosting both an ionised gas disc and a starburst ring (Phillips et al. 1983; Kriss et al. 2011; Fischer et al. 2015), where star formation is currently ongoing (star formation rate, SFR ∼ 5 M⊙ yr−1, Gruppioni et al. 2016).

Mrk 509 is a complex system showing evidence of both an ongoing minor merger and multiphase gas winds. The minor merger with a gas-rich dwarf is suggested by an ionised gas linear tail (Fischer et al. 2013, 2015), whereas Liu et al. (2015) found evidence of a [OIII] quasi-spherical wind powered by the AGN in the inner ∼2 kpc region, with velocity 290 km s−1 and mass outflow rate Ṁof = 5 M⊙ yr−1. The ionised wind appears spatially anti-correlated with the Hβ emission, suggesting suppression of star formation in this region. The gas kinematics suggests that the ionised wind is physically unrelated to the linear tail, and their apparent overlap is due to projection effects alone (Liu et al. 2015). On nuclear scales, highly ionised semi-relativistic gas winds are detected with velocity v ∼ 0.15−0.2c and outflow rates in the range 0.005−0.05 M⊙ yr−1 (Dadina et al. 2005; Cappi et al. 2009; Tombesi et al. 2011, 2012; Detmers et al. 2011). Simultaneous observations performed by the Chandra Telescope and Hubble Space Telescope with the Cosmic Origins Spectrograph (HST-COS) showed that the wind detected through UV absorption lines is located within ∼200 pc of the active nucleus, and that X-ray warm absorbers and UV absorbers are related (Ebrero et al. 2011; Kriss et al. 2019). Given the widespread evidence of winds from nuclear to kiloparsec scales and in different gas phases (Phillips et al. 1983; Kriss et al. 2011; Tombesi et al. 2011, 2012; Liu et al. 2015), Mrk 509 is an optimal candidate to investigate the wind driving mechanisms, stellar feedback versus AGN feedback, and their impact on the host-galaxy ISM. In this paper we exploit ALMA data of Mrk 509 to trace the distribution and kinematics of cold molecular gas and to look for evidences of molecular winds. The CO(2−1) and underlying continuum observations examined here are part of the ALMA survey of the IBISCO sample. The IBISCO sample is drawn from the International Gamma-Ray Astrophysics Laboratory (INTEGRAL) Imager on Board INTEGRAL Spacecraft (IBIS) catalogue of hard X-ray (20−100 keV) selected AGN (Malizia et al. 2012). The IBISCO sample is selected to have z < 0.05, AGN bolometric luminosity Lbol > 1043 erg s−1, and accurate black hole mass measurements. Being hard X-ray selected, the IBISCO sample is unbiased against nuclear obscuration, similarly to BASS, the BAT AGN Spectroscopic Survey (Koss et al. 2017, 2021). The aim of the ALMA survey of IBISCO is obtaining a reliable overview of the cold gas kinematics in the central ∼kpc of the AGN host galaxies with angular resolutions in the range 100−200 pc.

The paper is organised as follows. Section 2 presents the ALMA observational set-up and data reduction. Section 3 presents the observational results, in particular the continuum emission, the gas properties, and the CO kinematics. Section 4 describes the dynamical modelling, and Sect. 5 the properties of the molecular wind. In Sect. 6 we discuss our results, and in Sect. 7 we present our conclusions and future perspectives.

2. Observations and data reduction

Mrk 509 was observed with ALMA (program ID 2017.1.01439.S) during cycle 5 in September 2018 for a total integration time of 75 min in the frequency range 221.8−223.7 GHz (band 6), which covers the CO(2−1) line (observed redshifted frequency of 222.9 GHz) and the underlying 1.2 mm continuum. The obervations were conducted with the 12 m array in the C43-5 configuration with 44 antennas, minimum baseline of 14.6 m, and maximum baseline of 1397.9 m. This provided us with an angular resolution of about 0.4 arcsec and largest angular scale of ∼4.6 arcsec.

We used CASA 5.4.1 software (McMullin et al. 2007) to produce a set of calibrated visibilities and to generate clean data cubes. We calibrated the visibilities in the pipeline mode using the default calibrators provided by the observatory: bandpass and flux calibrators J2148+0657 and J2000–1748, J2025–0735 phase–amplitude calibrator. The absolute flux accuracy is better than 10%. To estimate the rest-frame 1.2 mm continuum emission we averaged the visibilities in the four spectral windows, excluding the spectral range covered by the CO(2−1) emission line. To estimate the continuum emission underlying the CO(2−1) line, we modelled the visibilities in the spectral window (spw) containing the emission line with the uvcontsub task, adopting a first-order polynomial. We verified that using a zero-order polynomial did not affect the results. We subtracted the continuum from the global visibilities and created the continuum-subtracted CO(2−1) visibility table. We imaged CO(2−1) and the continuum using the tclean task by adopting a natural weighting scheme, and produced a clean cube by using the hogbom cleaning algorithm with a detection threshold of three times the rms noise, a pixel size of 0.05 arcsec, and a spectral resolution of 10.6 km s−1. With this procedure we obtained a clean map of the continuum with a synthesised beam of 0.396 × 0.343 arcsec2 at PA = 51 deg, which corresponds to a spatial resolution of ∼285 × 244 pc2. The resulting rms noise in the clean continuum map is 0.014 mJy beam−1 in the aggregated bandwidth. From the continuum-subtracted CO(2−1) visibilities we produced a clean data cube with a beam size of 0.412 × 0.353 arcsec2 at a PA = 48 deg and a rms noise of 0.38 mJy beam−1 for a channel width of 10.6 km s−1.

3. Results

Figure 1, left panel, shows the 1.2 mm continuum map. The peak position of the continuum flux density map, obtained through 2D fitting in the image plane with CASA, is consistent with the AGN position reported by the NASA/IPAC Extragalactic Database (NED1). We measure a continuum flux density S1.2 mm = 2.17 ± 0.04 mJy, and a flux at the peak position Speak = 2.09 ± 0.02 mJy beam−1. According to our 2D fit the continuum emission is consistent with an unresolved source.

thumbnail Fig. 1.

Continuum emission, CO(2−1) integrated flux maps of Mrk 509 and the CO(2−1) spectrum extracted from the continuum-subtracted clean data cube. Left panel: 1.2 mm continuum map of Mrk 509. Regions with emission below 5σ have been blanked. The synthesised beam (0.396 × 0.343 arcsec2) is shown by the grey filled ellipse. Central panel: CO(2−1) integrated flux map, where a mask with threshold of 3σ has been applied. Black contours show the 1.2 mm continuum emission at (5, 30, 100) × σ. The synthesised beam is shown in the bottom left corner. Right panel: CO(2−1) spectrum (black crosses) extracted from the continuum-subtracted clean data cube from the regions where emission is above 3σ. The solid red line shows the multi-Gaussian fit; the magenta lines show the three Gaussian components.

Figure 1, central panel, shows that the CO(2−1) surface brightness (moment 0 map) is distributed over a region of approximately 8 arcsec (i.e. 5.5 kpc) and, in the inner ∼1.4 kpc region, shows an elongated shape along the north-west to south-east direction, crossing the AGN position. A CO-depleted region, which is a region with a CO(2−1) surface brightness lower than the neighbouring regions, is detected located at 1.4 kpc west of the AGN. The CO(2−1) spectrum (Fig. 1, right panel), extracted from the clean data cube in the region defined by a mask with 3σ threshold, shows a multiple peak profile reminiscent of disc-like kinematics. We fitted the emission line with a combination of three Gaussian profiles whose parameters are reported in Table 1. From the best fit we measured an integrated flux density of SCO = 36.7 ± 0.8 Jy km s−1 over a line width of 180 km s−1. We then derived the CO(2−1) luminosity L CO(2-1) ′ =(5.0±0.1)× 10 8 $ L^\prime_{\rm CO(2{-}1)} = (5.0 \pm 0.1)\times 10^8 $ K km s−1 pc2 using the relation of Solomon & Vanden Bout (2005). To estimate the proper conversion factor αCO, we used the stellar mass M* and the conversion function from Accurso et al. (2017). The M* is derived from broad-band spectral energy distribution (SED) decomposition including the emission of stars, dust heated by star formation and the AGN dusty torus (Gruppioni et al. 2016). We report in Appendix A details on the SED fitting. We derive a total stellar mass of M* = 1.2 × 1011 M⊙ and a αCO = 3.3 M⊙ (K km s−1 pc2)−1. Using these values, we derive a molecular gas mass MH2 = (1.65 ± 0.04)×109 M⊙ within a region of ∼5.2 kpc diameter. For the CO-depleted region we derived a 3σ upper limit molecular gas mass MH2 < 0.35 × 107 M⊙ assuming the same conversion factor.

Table 1.

Parameters of the Gaussian fit of the CO(2−1) line.

The mean-velocity map (moment 1 map, Fig. 2 top left panel) shows a gradient oriented north-east to south-west along a position angle (PA) of ∼250 deg, with range from −200 to 200 km s−1, which likely indicates an inclined rotating disc (PA in degrees is measured anti-clockwise from north). Figure 2 (top right panel) shows the CO(2−1) velocity dispersion (moment 2 map), σdisp, with values of 10−20 km s−1 in the outer regions. In the inner 1.4 × 1.4 kpc2 (corresponding to 2 × 2 arcsec2), σdisp has a butterfly-shaped morphology with a σdisp = 80 km s−1 towards the central beam, where beam-smearing effects can boost the σdisp (e.g. Davies et al. 2011). A more reliable estimate of the real σdisp in the nuclear region is 50−60 km s−1, as measured in an annulus of 0.3 < R < 0.6 arcsec. We detect a region with enhanced σdisp = 30−40 km s−1 located south of the CO-depleted region, suggesting that also at this location dynamic perturbations are likely present.

thumbnail Fig. 2.

Comparison between moment 1 and 2 maps generated from data, 3DBAROLO disc model and residuals. Top panels: CO(2−1) moment 1 (top left) and moment 2 (top right) maps. A threshold of 3σ has been applied. The black cross gives the AGN position. Central panels: moment 1 (centre left) and moment 2 (centre right) maps of the 3DBAROLO disc model. The black dashed lines give the orientation of the major and minor kinematic axes (PA = 248 and 338 deg, respectively). Bottom panels: moment 1 (bottom left) and moment 2 (bottom right) maps of the residuals obtained by subtracting the model velocity and velocity dispersion maps from the observed ones. A, B, C, and D indicate the main dynamical perturbations discussed in the text.

4. Molecular discs

Based on the moment maps we built a dynamical model of the system using the 3D-Based Analysis of Rotating Objects from Line Observations (3DBAROLO, Di Teodoro & Fraternali 2015).

Any deviation from rotating disc kinematics can then be identified by comparing the disc dynamical model with the observed data, through the residual data cube, which is obtained by subtracting the model from the observed data cube. We fitted a 3D tilted-ring model to the emission line data cube to provide a first-order dynamical model of the CO-emitting gaseous disc. In the first run we allowed four parameters to vary: rotation velocity, velocity dispersion, the disc inclination, and position angle. We fixed the kinematic centre to the 1.2 mm continuum position (Sect. 3). We then ran the model again fixing the inclination, so only PA, rotation velocity, and velocity dispersion are allowed to vary. The inclination from the minor/major axis ratio from UV imaging (GALEX archive) is 44 ± 10 deg, and we fixed it to 44 deg in the model. The first run produces amplitude residuals of about 5% in the moment 0 map that further increase by adding a radial velocity component. The second run, with fixed inclination, produces amplitude residuals of about 2%; therefore, we adopt the latter as best-fit disc model. The disc model has a position angle PA = 248 ± 5 deg, a vLOS that ranges from −200 to 200 km s−1, and a velocity dispersion that decreases from 70 km s−1 in the inner region to ∼10 km s−1 in the outer part of the disc (Fig. 2, central panels). Figure 3 shows the rotation velocity (left panel) of the best-fit disc model, and the virial dynamical mass (right panel) M dyn =r v rot 2 /2G $ M_{\rm dyn} = r{\it v}^2_{\rm rot}/2G $ as a function of the radius out to ∼2.5 kpc. The rotation velocity shows a maximum value of ∼300 km s−1 and the dynamical mass enclosed within the maximum radius is ∼2.1 × 1010 M⊙.

thumbnail Fig. 3.

Rotational velocity and dynamical mass obtained from the best-fit disc model of the main disc. Left panel: rotational velocity vrot vs. radius of the best-fit disc model of the main disc. Right panel: dynamical mass vs. radius within the inner ∼2.5 kpc radius, derived from the relation M dyn =r v rot 2 /2G $ M_{\rm dyn} = r{\it v}^2_{\rm rot}/2G $.

We computed the residual velocity and velocity dispersion maps by subtracting the disc model from the data (Fig. 2, bottom panels). Residual maps indicate that the model nicely describes the observed large-scale velocity gradient, from r ∼ 1 arcsec to the outer boundary. Dynamical perturbations are detected in the residual velocity map (Fig. 2, bottom left panel) at positions A, C (redshifted) and B, D (blueshifted). In the velocity dispersion residual map (Fig. 2, bottom right panel) residuals are present only at positions A and B. The position-velocity (PV) plots shown in Fig. 4 capture the main kinematics components seen in the moment maps and in the residual maps (Fig. 2). The slices are taken through the continuum peak position, using a slit 0.45 arcsec wide (i.e. approximately the FWHM size of the synthetic beam major axis). The PV diagram along the kinematic major axis (PA = 248 deg) shows kinematic perturbations at offset 0.5 arcsec with VLOS ∼ 250 km s−1 (A), and at offset 2 arcsec with VLOS ∼ 150 km s−1 (B), in addition to the typical disc rotation pattern. The PV diagram along the kinematic minor axis (PA = 338 deg) shows both blue-shifted and red-shifted non-rotational motions with velocities that reach 300 km s−1 and −200 km s−1 in the inner ±0.5 arcsec region (C, D). The perturbations C and D suggest the presence of a warped disc or bar. The nuclear region is barely resolved at the current angular resolution. Here we attempt to model it with a nuclear disc, as detailed in Appendix B. The modelling with a disc of inclination fixed to 66 deg and PA fixed to 330 deg (i.e. along the two brightest emission regions in the residual intensity map, C and D; black contours in Fig. 5, left panel) is able to account for and nicely model the residual emission (Fig. 5, right panel, and Fig. 6, top panel), suggesting the presence of a nuclear warped disc. We note that observations with higher angular resolution would be required to better assess the kinematic properties in this region. The modelling of the nuclear region does not impact the analysis of the molecular wind presented in the next section.

thumbnail Fig. 4.

Position-velocity diagrams along the kinematic major and minor axis. Left panel: PV plot along PA = 248 deg (the kinematic major axis). Right panel: PV plot along PA = 338 deg (the kinematic minor axis). The slit width is set to 0.45 arcsec (approximately the FWHM size of the synthetic beam major axis). The white contours refer to data, and are drawn at (1, 2, 4, 8, 16, 32) × σ; the black dots give the VLOS value from the 3DBAROLO model.

thumbnail Fig. 5.

Moment 1 map and PV diagram from the residual cube. Left panel: velocity map generated from the residual cube obtained subtracting the disc model from the data cube (a mask with a 2σ threshold on the residual intensity map was applied). The black contours show the CO(2−1) emission in residual intensity map at (3, 4, 5, 6, 7) × σ. The black dashed lines represent the slice adopted for the PV diagram (right panel) and comprise the two main kinematic perturbations in the residual intensity map, along PA = 330 deg. The black cross gives the AGN position. Right panel: PV diagram along PA = 330 deg (kinematic major axis). The red contours and yellow points represent the disc model of this inner region; the blue contours represent the data (contours are drawn at (2, 4, 6) × σ). A, C, and D indicate the main dynamical perturbations.

5. Molecular wind

Regions A and B show residuals both in the velocity and in the velocity dispersion maps. Figure 6, top panel, shows the intensity map of the residuals. The spectra extracted from A and B are shown in Fig. 6, bottom panels. The gas in region A has a peak velocity of 220 km s−1, and a red tail reaching 300 km s−1, whereas the gas in B is characterised by a peak velocity of about 170 km s−1. We note that the disc rotational velocity at A and B is 70 and 210 km s−1 (dashed lines), respectively, thus showing that gas in these regions is not participating in the disc rotation; instead, it can be associated with a wind. The wind mass outflow rate is computed using the relation from Fiore et al. (2017), assuming a spherical sector:

M ˙ of = 3 M of v of R of · $$ \begin{aligned} \dot{M}_{\rm of} = 3 \frac{M_{\rm of}\,{ v}_{\rm of}}{R_{\rm of}}\cdot \end{aligned} $$(1)

thumbnail Fig. 6.

Intensity map and spectra of regions A and B. Top panel: CO(2−1) residual intensity map after the best-fit disc model has been subtracted from the data cube. The red and blue circles indicate regions A and B respectively. The black contours show the emission at (2, 3, 4, 5) × σ. Bottom panels: spectra extracted from region A (red histogram) and B (blue histogram). The dashed black lines indicate the value of the disc rotational velocity at A and B based on our best-fit dynamical model (vdisc = 75 km s−1 around region A, and vdisc = 210 km s−1 around region B).

Here vof = v98 is the wind velocity, the velocity enclosing 98% of the cumulative velocity distribution of the outflowing gas (Bischetti et al. 2019), Mof is the wind molecular mass, and Rof is the distance from the active nucleus reached by the wind. Table 2 gives v98, Mof, and Rof for regions A and B. Because we do not know a priori whether the gas in the wind is optically thin or thick, we derived the molecular mass range in A and B adopting a conversion factor αCO = 0.3 M⊙ (K km s−1 pc2)−1 for the optically thin case (Morganti et al. 2015; Dasyra et al. 2016) and αCO = 0.8 M⊙ (K km s−1 pc2)−1 (Bolatto et al. 2013; Morganti et al. 2015) for the optically thick case. These assumptions imply wind masses in the ranges MH2 = (2.1 − 5.6)×106 M⊙ and MH2 = (2.0 − 5.4)×106 M⊙, thus a molecular outflow rate in the range Ṁof = 5.5 − 14.7 M⊙ yr−1 and Ṁof = 0.86 − 2.30 M⊙ yr−1 for region A and B, respectively. The lowest boundaries are pertinent to optically thin gas, and thus represent lower limits to the molecular mass and outflow rate. These estimates are derived assuming AGN-driven winds because Rof is the distance from the AGN.

Table 2.

Parameters of the molecular wind components.

We evaluated the escape velocity at A and B using the relation

v esc = G M dyn R of , $$ \begin{aligned} { v}_{\rm esc} = \sqrt{\frac{GM_{\rm dyn}}{R_{\rm of}}}, \end{aligned} $$(2)

where we used the dynamical mass derived from the best-fit disc model (Fig. 3, right panel). We find that vesc ∼ 140 km s−1 at A, which is significantly less than the wind velocity, suggesting that the outflowing gas can leave the innermost 300 pc region. The escape velocity at B is ∼227 km s−1, similar to wind B velocity. Computing the time needed for the gas to travel from the centre to the wind position and assuming a constant velocity equal to v98, we obtain 1.16 Myr to reach A, and 6.85 Myr to reach B. Therefore, these two outflowing gas components might represent two consecutive episodes of winds launched at different times, or the same event with a decelerating wind (e.g. Fischer et al. 2013; Veilleux et al. 2020). Current data cannot discriminate between these two scenarios.

6. Discussion

Mrk 509 is classified as a spheroid based on optical imaging (Ho & Kim 2014), and it hosts an ionised gas ring and disc, where star formation is currently active (Kriss et al. 2011; Fischer et al. 2015). With its SFR of 5.1 ± 0.5 M⊙ yr−1 and the stellar mass of M* = (1.2 ± 0.1)×1011 M⊙ (Gruppioni et al. 2016, Appendix A), it lies on the main sequence of star-forming galaxies. Our observations show that Mrk 509 hosts a molecular gas reservoir of MH2 = (1.7 ± 0.1)×109 M⊙ within 5.2 kpc. The stellar mass within the same region was derived by scaling the total stellar mass assuming an exponential profile and the size of the galaxy in K band (33 kpc, NED). It turns out that M*(R < 2.6 kpc)∼3 × 1010 M⊙. The molecular gas fraction μ = MH2/M* in the inner ∼5.2 kpc is then μ ∼ 5%, consistent with that derived for local star-forming galaxies with the same stellar mass (e.g. Saintonge et al. 2011, 2017). We find that the CO(2−1) moment 1 map shows a velocity gradient consistent with a inclined rotating disc with inclination consistent with that derived from UV imaging of the host galaxy. The derived dynamical mass is Mdyn = (2.0 ± 1.1) × 1010 M⊙ within a ∼2.6 kpc radius. The uncertainty is dominated by the uncertainty on the inclination. The dynamical mass is thus consistent with the M*(R < 2.6 kpc). Our dynamical modelling also suggests a possible warp in the disc towards the centre, as reported by Combes (2021) in other nearby Seyfert galaxies. This nuclear region, however, remains marginally resolved in current observations, and we cannot exclude more complex structures, such as bars or rings.

The molecular disc kinematics and inclination are consistent with those reported for the ionised gas phase by Fischer et al. (2015) and Liu et al. (2015), who mapped the ionised disc on scales of up to 4 kpc (see Fig. 6 in Liu et al. 2015). We compare our CO(2−1) moment 0 map with the HST image in FQ508N filter, which traces ionised gas emitting in the [OIII]λ5007 Å, to investigate the relation between cold molecular and ionised gas distributions. Figure 7 left panel shows that the molecular gas (red contours) is detected at the position of the starburst ring. The [OIII] image shows a linear tail extending out from the starburst ring towards the north-west. The ‘jut’ of tail that seems to bend towards the nucleus is also visible. These two structures have been interpreted as being due to a minor merger with a gas-rich dwarf galaxy (Fischer et al. 2015). In our ALMA data we do not detect CO-emitting gas at the position and velocity of the linear tail. The 3σ upper limit on the H2 mass in the tail is < 107 M⊙, based on the rms of our observation, adopting a Milky Way conversion factor and FWHM = 50 km s−1 for the CO line (see Aalto et al. 2001). The limit appears consistent with a minor merger scenario (Knierman et al. 2013; van de Voort et al. 2018). The tidal feature and the presence of a molecular disc with an ongoing starburst are both in agreement with the scenario where galaxy interactions and mergers produce gas destabilisation, and feed star formation and AGN activity (e.g. Menci et al. 2014, and references therein). Assuming that most star formation occurs within the starburst ring, we derive a star formation efficiency SFE = SFR/MH2 = 6.2 × 10−9 yr−1, consistent with that found for local early-type galaxies with recent minor mergers (Saintonge et al. 2012; Davis et al. 2015).

thumbnail Fig. 7.

HST [OIII]λ5007 Å image and Toomre Q-parameter map of Mrk 509. Right panel: HST image in the FQ508N filter, probing rest-frame [OIII]λ5007 Å emission (greyscale). The red contours represent CO(2−1) emission (Fig. 1centre panel, contours are drawn at (1, 2, 3, 4, 5, 7, 10, 15, 20, 25, 30, 37) × σ). The blue contours represent CO(2−1) emission in region A and B (Fig. 6 top panel, contours are drawn at (3, 4, 5) × σ). Image from the Mikulski Archive for Space Telescopes (MAST), HST Proposal 12212. Right panel: Toomre Q-parameter map for the cold molecular gas component. The black contours represent the moment 0 CO emission at (2, 5, 15, 30, 45) × σ. The white cross gives the continuum peak position. The white dashed ellipses indicate projected annuli with radius R = 0.3, 0.8, 1.4, 2.2 arcsec (see Table 3).

We quantify the disc stability against gravitational fragmentation by deriving a map of the Toomre Q-parameter (Toomre 1964) for the cold molecular gas component, Qgas. The Toomre parameter is related to the local gas velocity dispersion σdisp, the circular velocity vrot, and the gas surface density Σgas at any given radius R, through the relation

Q gas = a v rot σ disp π R G Σ gas , $$ \begin{aligned} Q_{\rm gas} = \frac{a { v}_{\rm rot} \sigma _{\rm disp}}{\pi RG \Sigma _{\rm gas}}, \end{aligned} $$(3)

where a is a constant that varies from 1 for a Keplerian rotation curve and 1.4 for a flat rotation curve to 2 for a solid-body rotation curve (see Genzel et al. 2011, 2014 and references therein). To obtain a map of Qgas we use vrot, derived through the 3DBAROLO model of the whole disc, using the assumed inclination and best-fit PA. We adopt a = 1 up to R = 1.1 arcsec, where the rotation curve is steep, and a = 1.4 at larger radii. The velocity dispersion and Σgas are calculated pixel by pixel from the moment 2 and 0 maps, respectively; the radius R is deprojected. The map of Qgas is shown in Fig. 7 right panel, while in Table 3 we report σdisp, Σgas, and Qgas in five annuli between R = 0.3 and 0.8 arcsec, 0.8 and 1.4 arcsec, 1.4 and 2.2 arcsec (which includes a ring-like feature of enhanced σdisp visible in Fig. 2), R > 2.2 arcsec, and in the whole disc excluding the central region of 0.3 arcsec radius, to avoid the beam-smearing effect. We find values that range from Qgas ≈ 0.5 up to 10, and we note that Qgas exceeds 1 across most of the molecular disc. A similar Qgas range has been reported in nearby Seyfert galaxies (e.g. García-Burillo et al. 2003; Sani et al. 2012; Romeo & Fathi 2016), in M 51 (Hitschfeld et al. 2009), and in local spirals (Romeo et al. 2020). Where the value of the Q-parameter is higher than the critical value Qcrit the disc is unstable against gravitational collapse. The commonly adopted critical value of Qgas is in the range 1−3 (Genzel et al. 2014; Leroy et al. 2008), whereas if the total Q tot − 1 = Q gas − 1 + Q star − 1 $ Q_{\mathrm{tot}}^{-1} = Q_{\mathrm{gas}}^{-1} + Q_{\mathrm{star}}^{-1} $ is considered, including both the gas and stellar contribution, the critical value is ∼1 (Aumer et al. 2010). The Qgas we find in Mrk 509 implies that the disc is unstable in its outer regions across the starburst ring, whereas it is stable against fragmentation at the nucleus and in a lopsided ring-like structure located at R ∼ 2 arcsec from the AGN. To draw more robust conclusions on the role that disc gravitational instability plays in Mrk 509, a multi component analysis of Qgas also taking into account the stellar component is needed, as shown by Romeo & Falstad (2013), Romeo & Mogotsi (2017). We find significant kinematic perturbations at different locations across the disc where the molecular gas shows deviations from the disc rotation pattern. The first (region A) is located at a distance of ∼300 pc from the nucleus, the second (region B) further out at a distance of about 1.4 kpc. In the proximity of region B we detect a CO-depleted region (Fig. 6). This is close (in projection) to the location where the linear tail related to the merger event observed with HST crosses the disc. The CO-depleted region may be related to the merger event locally decreasing the CO in the disc; however, the optical image does not show any lack of ionised gas there. In Sect. 5 we show that regions A and B both show non-circular gas motions, which are likely due to outflowing gas. In region A integral field spectroscopy (IFS) spectroscopy shows a [OIII] wind with velocity similar to the molecular value (Liu et al. 2015), supporting the notion that this region hosts a multiphase wind. The ionised and molecular winds are co-spatial, within the accuracy allowed by the angular resolution of the data, suggesting a cooling sequence (Richings & Faucher-Giguere 2017; Menci et al. 2019). Region B has a small projected distance from the tidal tail seen in [OIII]. The ionised gas velocity in the tidal tail is −100 km s−1, according to Liu et al. (2015), while it is 200 km s−1 for the molecular gas in B, and 250 km s−1 in A. Therefore, it is unlikely that the molecular gas in A and in B is related to the tidal tail.

Table 3.

Average Qgas and velocity dispersion in 4 annulii.

It is likely that the molecular gas in A and B is involved in an AGN-driven wind, whose total molecular outflow rate is in the range 6.4−17.0 M⊙ yr−1, for the optically thin and thick cases respectively. The contribution of region B to the global outflow budget is a factor of ∼10 smaller than region A; therefore, including or excluding region B does not change the main result of the work. We computed the wind momentum load, which is the ratio of the outflow momentum flux, Ṗof = Ṁof · v98, to the radiative momentum flux, Ṗrad = Lbol/c. We adopt a bolometric luminosity Lbol = 1044.99 erg s−1 derived by Duras et al. (2020) from SED fitting and consistent with that derived by Gruppioni et al. (2016). The wind momentum load values obtained are reported in Fig. 8 for the different wind components: the three ultra-fast ouflow (UFO) components from Tombesi et al. (2010, 2011, 2012), the ionised gas winds from Liu et al. (2015), and the molecular winds from this work. The molecular wind A and the ionised gas wind show similar velocities and arise from approximately the same region, and also show similar mass outflow rates, suggesting that they might be part of a multiphase wind. The molecular wind A and the ionised wind are consistent with a momentum-conserving driving (black dashed line), and are well below the expectation for energy-conserving winds. Wind B has a much lower mass outflow rate. From current data we cannot determine whether A and B components of the molecular wind are due to consecutive episodes or to the same event. Assuming that the molecular components and the ionised wind detected by Liu et al. (2015) are both part of the same multiphase wind, we derive a total mass outflow rate, molecular plus ionised, of ∼22 M⊙ yr−1. For comparison, the accretion rate onto the SMBH is Ṁacc = Lbol/ϵc2 ∼ 0.176 M⊙ yr−1, assuming an efficiency ϵ = 10%. The total wind momentum load, Ṗof/Ṗrad ∼ 1.4 ± 0.6 (magenta point) is below the energy-conserving wind expectation, and compatible with a momentum-conserving wind powered by the nuclear semi-relativistic wind probed by the UFOs. This momentum boost can also be attained by radiation pressure on dust (Ishibashi et al. 2021). The wind mass-loading factor is η = Ṁof/SFR = 4.4, in agreement with the value found for molecular winds in other Seyfert galaxies (e.g. Fiore et al. 2017).

thumbnail Fig. 8.

Wind momentum load, Ṗof/Ṗrad, as a function of the wind velocity for Mrk 509 (Tombesi et al. 2011, 2012; Liu et al. 2015, and this work), and a compilation of local AGN (PDS456: Nardini et al. 2015; Bischetti et al. 2019, Mrk 231: Feruglio et al. 2015, and IRASF11119+3257: Tombesi et al. 2015). The red symbols are UFOs; the green star is the ionised [OIII] wind for Mrk 509; the blue symbols are molecular winds, and the magenta star is the total wind (ionised plus molecular phases) for Mrk 509. The black dashed line indicates the expectation for momentum-conserving wind driven by the AGN in Mrk 509. The red shaded area represents the prediction for an energy-conserving wind with Ṗmol/Ṗrad = vUFO/v.

In the following we discuss the possibility that the outflows are driven by energy injected by supernova (SN) explosions occurring in the starburst ring. To derive the mass outflow rate at B in this case, we adopt a radius of R = 0.4 kpc, which is half the projected size of the starburst ring. This implies Ṁ = 3 − 8 M⊙ yr−1 by adopting Eq. (1), and Ṗ = 0.4 − 1.0 × 1041 erg s−1. Assuming that the total SFR (5 M⊙ yr−1) is homogeneously distributed across the ring, we consider a portion of the starburst ring equal to one-tenth of the ring in area, and we estimate a kinetic power of a SN-driven wind Ėkin,SN = 3.5 × 1041 erg s−1, assuming a SN rate of 0.02(SFR/M⊙ yr−1) (e.g. Veilleux et al. 2005). Outflow B is then energetically consistent with a SN-driven wind originating from the starburst ring (Venturi et al. 2018; Leaman et al. 2019). Regarding wind A, the SFR at its location cannot be resolved in current data, but it is most likely a small fraction of the total SFR which occurs mainly in the starburst ring. Therefore outflow A is most likely driven by the AGN.

Figure 9 shows the outflow rate, Ṁof, versus the AGN bolometric luminosity for Mrk 509 and a compilation of AGN adapted from Fiore et al. (2017). Mrk 509 shows an ionised wind and X-ray winds that are broadly consistent with the correlations found by Fiore et al. (2017), whereas its molecular outflow rate is significantly below the best-fit correlation of molecular winds. The ionised outflow rate is similar to the molecular rate, unlike the average AGN with similar Lbol. We note, however, that the correlation is derived from a biased sample which collects mainly high luminosity objects and known massive winds, and therefore is probably biased towards positive detection, boosting the outflow rates. This highlights the importance of performing statistical studies on unbiased AGN samples to derive meaningful scaling relations. This study will be carried out by using the whole IBISCO sample in a separate publication.

thumbnail Fig. 9.

Wind mass outflow rate vs. AGN bolometric luminosity for Mrk 509 and a compilation of AGN from (Fiore et al. 2017). For Mrk 509 we plot the Ṁof of the total molecular wind (this work, light blue circle), the ionised wind (green circle Liu et al. 2015), and the UFO (Tombesi et al. 2010, 2011, 2012, red circle). The size of the circles represents the uncertainty in Ṁof. The small blue symbols represent molecular gas winds; the small green symbols represent ionised gas winds; and the small red symbols represent highly ionised gas winds detected in X-rays. The dashed blue, green, and red lines are the best-fit correlations of the molecular, ionised, and X-ray absorbers samples, respectively. Adapted from Fiore et al. (2017).

7. Conclusions

We presented an analysis of the CO(2−1) line and 1.2 mm continuum of Mrk 509, a Seyfert 1.5 galaxy drawn from the IBISCO sample of hard-X-ray selected local AGN. Mrk 509 is optically classified as a spheroid and hosts ionised winds in different gas phases, from warm ionised gas wind traced by [OIII] to highly ionised UFOs seen in X-rays. We report the following findings:

  1. The galaxy hosts a molecular gas reservoir of MH2 = (1.7 ± 0.1)×109 M⊙ within the inner 5.2 kpc, and partially overlapping with the starburst ring seen in HST optical imaging. In this region we estimate a molecular gas fraction μ ∼ 5%, consistent with that of local star-forming galaxies with the same stellar mass (e.g. Saintonge et al. 2011, 2017). The signatures of a minor merger and the presence of a molecular disc with an ongoing starburst are both in agreement with the scenario where galaxy interactions and mergers produce gas destabilisation, and feed both star formation and AGN activity. The star formation efficiency across the molecular disc and starburst ring, SFE = 6.2 × 10−9 yr−1, is consistent with that found for local early-type galaxies with recent minor mergers.

  2. We quantify the disc stability by estimating the spatially resolved Toomre Q-parameter for the cold molecular gas component, Qgas. We find that Qgas varies across the molecular disc in the range ∼0.5 − 10. The disc is unstable across the starburst ring and stable against fragmentation at the nucleus, in a lopsided ring-like structure located at R ∼ 2 arcsec from the AGN and at the location of the molecular winds.

  3. The main velocity gradient detected in CO(2−1) is well modelled by a disc with Mdyn = (2.0 ± 1.1) × 1010 M⊙ within ∼5.2 kpc and inclination 44 ± 10 deg with respect to the line of sight. The CO kinematics in the nuclear region, r ∼ 700 pc, is barely resolved at the current angular resolution and may be affected by a warped nuclear disc. Higher angular resolution observations may help to better constrain the gas distribution and kinematics in the nuclear region.

  4. We find significant perturbations of the molecular gas kinematics at two different locations in the disc, where the molecular gas shows deviations from the disc rotation, that we interpret as molecular winds. Wind A has a velocity of v98 = 250 km s−1 and is located at a distance of ∼300 pc from the AGN, in the same region where an ionised gas wind is detected. Wind B is found at a distance of about 1.4 kpc, overlapping on the line of sight with the starburst ring, and at a small projected distance from the tidal tail. Its velocity v98 = 200 km s−1 suggests that B is not related to the tidal tail. The total molecular outflow rate for A and B is in the range 6.4 − 17.0 M⊙ yr−1, for the optically thin and thick cases, respectively, with outflow B constituting only a small fraction of the total outflow rate (about 10%). The wind kinetic energy is consistent with a momentum-conserving wind driven by the AGN. The spatial overlap of the ionised and molecular winds in A, and their similar velocity, suggest a cooling sequence within a multi-phase AGN-driven wind, whereas outflow B is consistent with a SN-driven wind originating in the starburst ring.

Acknowledgments

We thank the referee for his/her careful reading and insightful suggestions that helped improving the paper. This paper makes use of the ALMA data from project: ADS/JAO.ALMA#2017.1.01439.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. We acknowledge financial support from PRIN MIUR contract 2017PH3WAT, and PRIN MAIN STREAM INAF “Black hole winds and the baryon cycle”.

References

  1. Aalto, S., Hüttemeister, S., & Polatidis, A. G. 2001, A&A, 372, L29 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Aalto, S., Falstad, N., Muller, S., et al. 2020, A&A, 640, A104 [EDP Sciences] [Google Scholar]
  3. Accurso, G., Saintonge, A., Catinella, B., et al. 2017, MNRAS, 470, 4750 [NASA ADS] [Google Scholar]
  4. Alonso-Herrero, A., Pereira-Santaella, M., Rigopoulou, D., et al. 2020, A&A, 639, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Audibert, A., Combes, F., García-Burillo, S., & Dasyra, K. 2020, Proc. IAU, 15, 307 [CrossRef] [Google Scholar]
  6. Aumer, M., Burkert, A., Johansson, P. H., & Genzel, R. 2010, ApJ, 719, 1230 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bentz, M. C., Peterson, B. M., Netzer, H., Pogge, R. W., & Vestergaard, M. 2009, ApJ, 697, 160 [Google Scholar]
  8. Bischetti, M., Piconcelli, E., Feruglio, C., et al. 2019, A&A, 628, A118 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207 [CrossRef] [Google Scholar]
  10. Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
  11. Cappi, M., Tombesi, F., Bianchi, S., et al. 2009, A&A, 504, 401 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Cicone, C., Maiolino, R., Sturm, E., et al. 2014, A&A, 562, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Combes, F. 2019, ArXiv e-prints [arXiv:1911.00707] [Google Scholar]
  14. Combes, F. 2021, IAU Proc., 359, 312 [Google Scholar]
  15. Combes, F., García-Burillo, S., Audibert, A., et al. 2019, A&A, 623, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. da Cunha, E., Charlot, S., & Elbaz, D. 2008, MNRAS, 388, 1595 [Google Scholar]
  17. Dadina, M., Cappi, M., Malaguti, G., Ponti, G., & de Rosa, A. 2005, A&A, 442, 461 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Dasyra, K. M., & Combes, F. 2012, A&A, 541, L7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Dasyra, K. M., Combes, F., Novak, G. S., et al. 2014, A&A, 565, A46 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Dasyra, K. M., Combes, F., Oosterloo, T., et al. 2016, A&A, 595, L7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Davies, R., Förster Schreiber, N. M., Cresci, G., et al. 2011, ApJ, 741, 69 [NASA ADS] [CrossRef] [Google Scholar]
  22. Davis, T. A., Rowlands, K., Allison, J. R., et al. 2015, MNRAS, 449, 3503 [NASA ADS] [CrossRef] [Google Scholar]
  23. Detmers, R. G., Kaastra, J. S., Steenbrugge, K. C., et al. 2011, A&A, 534, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Di Teodoro, E. M., & Fraternali, F. 2015, MNRAS, 451, 3021 [Google Scholar]
  25. Duras, F., Bongiorno, A., Ricci, F., et al. 2020, A&A, 636, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Ebrero, J., Kriss, G. A., Kaastra, J. S., et al. 2011, A&A, 534, A40 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Fabian, A. C. 1999, MNRAS, 308, L39 [NASA ADS] [CrossRef] [Google Scholar]
  28. Fabian, A. C. 2012, ARA&A, 50, 455 [Google Scholar]
  29. Faucher-Giguère, C.-A., & Quataert, E. 2012, MNRAS, 425, 605 [Google Scholar]
  30. Feltre, A., Hatziminaoglou, E., Fritz, J., & Franceschini, A. 2012, MNRAS, 426, 120 [Google Scholar]
  31. Fernández-Ontiveros, J. A., Dasyra, K. M., Hatziminaoglou, E., et al. 2020, A&A, 633, A127 [Google Scholar]
  32. Ferrarese, L., & Ford, H. 2005, Space Sci. Rev., 116, 523 [Google Scholar]
  33. Feruglio, C., Maiolino, R., Piconcelli, E., et al. 2010, A&A, 518, L155 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Feruglio, C., Fiore, F., Carniani, S., et al. 2015, A&A, 583, A99 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. Feruglio, C., Ferrara, A., Bischetti, M., et al. 2017, A&A, 608, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Feruglio, C., Fabbiano, G., Bischetti, M., et al. 2020, ApJ, 890, 29 [Google Scholar]
  37. Fiore, F., Feruglio, C., Shankar, F., et al. 2017, A&A, 601, A143 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Fischer, T. C., Crenshaw, D. M., Kraemer, S. B., & Schmitt, H. R. 2013, ApJS, 209, 1 [Google Scholar]
  39. Fischer, T. C., Crenshaw, D. M., Kraemer, S. B., et al. 2015, ApJ, 799, 234 [NASA ADS] [CrossRef] [Google Scholar]
  40. Fluetsch, A., Maiolino, R., Carniani, S., et al. 2019, MNRAS, 483, 4586 [NASA ADS] [Google Scholar]
  41. Fritz, J., Franceschini, A., & Hatziminaoglou, E. 2006, MNRAS, 366, 767 [Google Scholar]
  42. García-Burillo, S., Combes, F., Hunt, L. K., et al. 2003, A&A, 407, 485 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. García-Burillo, S., Combes, F., Usero, A., et al. 2014, A&A, 567, A125 [Google Scholar]
  44. García-Burillo, S., Viti, S., Combes, F., et al. 2017, A&A, 608, A56 [Google Scholar]
  45. García-Burillo, S., Combes, F., Ramos Almeida, C., et al. 2019, A&A, 632, A61 [Google Scholar]
  46. Gebhardt, K., Bender, R., Bower, G., et al. 2000, ApJ, 539, L13 [Google Scholar]
  47. Genzel, R., Newman, S., Jones, T., et al. 2011, ApJ, 733, 101 [Google Scholar]
  48. Genzel, R., Förster Schreiber, N. M., Lang, P., et al. 2014, ApJ, 785, 75 [NASA ADS] [CrossRef] [Google Scholar]
  49. Gruppioni, C., Berta, S., Spinoglio, L., et al. 2016, MNRAS, 458, 4297 [Google Scholar]
  50. Hitschfeld, M., Kramer, C., Schuster, K. F., Garcia-Burillo, S., & Stutzki, J. 2009, A&A, 495, 795 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Ho, L. C., & Kim, M. 2014, ApJ, 789, 17 [NASA ADS] [CrossRef] [Google Scholar]
  52. Huchra, J., Latham, D. W., da Costa, L. N., Pellegrini, P. S., & Willmer, C. N. A. 1993, AJ, 105, 1637 [NASA ADS] [CrossRef] [Google Scholar]
  53. Ishibashi, W., Fabian, A. C., & Arakawa, N. 2021, MNRAS, 502, 3638 [NASA ADS] [CrossRef] [Google Scholar]
  54. King, A. 2003, ApJ, 596, L27 [Google Scholar]
  55. King, A., & Pounds, K. 2015, ARA&A, 53, 115 [NASA ADS] [CrossRef] [Google Scholar]
  56. Knierman, K. A., Scowen, P., Veach, T., et al. 2013, ApJ, 774, 125 [NASA ADS] [CrossRef] [Google Scholar]
  57. Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511 [Google Scholar]
  58. Koss, M., Trakhtenbrot, B., Ricci, C., et al. 2017, ApJ, 850, 74 [Google Scholar]
  59. Koss, M. J., Strittmatter, B., Lamperti, I., et al. 2021, ApJS, 252, 29 [CrossRef] [Google Scholar]
  60. Kriss, G. A., Arav, N., Kaastra, J. S., et al. 2011, A&A, 534, A41 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  61. Kriss, G. A., Arav, N., Edmonds, D., et al. 2019, A&A, 623, A82 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Leaman, R., Fragkoudi, F., Querejeta, M., et al. 2019, MNRAS, 488, 3904 [NASA ADS] [CrossRef] [Google Scholar]
  63. Leroy, A. K., Walter, F., Brinks, E., et al. 2008, AJ, 136, 2782 [Google Scholar]
  64. Liu, G., Arav, N., & Rupke, D. S. N. 2015, ApJS, 221, 9 [NASA ADS] [CrossRef] [Google Scholar]
  65. Longinotti, A. L., Vega, O., Krongold, Y., et al. 2018, ApJ, 867, L11 [CrossRef] [Google Scholar]
  66. Lutz, D., Sturm, E., Janssen, A., et al. 2020, A&A, 633, A134 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Malizia, A., Bassani, L., Bazzano, A., et al. 2012, MNRAS, 426, 1750 [Google Scholar]
  68. Marasco, A., Cresci, G., Nardini, E., et al. 2020, A&A, 644, A15 [EDP Sciences] [Google Scholar]
  69. McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in Astronomical Data Analysis Software and Systems XVI, eds. R. A. Shaw, F. Hill, & D. J. Bell, ASP Conf. Ser., 376, 127 [Google Scholar]
  70. Menci, N., Gatti, M., Fiore, F., & Lamastra, A. 2014, A&A, 569, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  71. Menci, N., Fiore, F., Feruglio, C., et al. 2019, ApJ, 877, 74 [NASA ADS] [CrossRef] [Google Scholar]
  72. Morganti, R., Holt, J., Saripalli, L., Oosterloo, T. A., & Tadhunter, C. N. 2007, A&A, 476, 735 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Morganti, R., Oosterloo, T., Oonk, J. B. R., Frieswijk, W., & Tadhunter, C. 2015, A&A, 580, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Morganti, R., Veilleux, S., Oosterloo, T., Teng, S. H., & Rupke, D. 2016, A&A, 593, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Nardini, E., Reeves, J. N., Gofford, J., et al. 2015, Science, 347, 860 [Google Scholar]
  76. Oosterloo, T., Raymond Oonk, J. B., Morganti, R., et al. 2017, A&A, 608, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  77. Phillips, M. M., Baldwin, J. A., Atwood, B., & Carswell, R. F. 1983, ApJ, 274, 558 [NASA ADS] [CrossRef] [Google Scholar]
  78. Richings, A. J., & Faucher-Giguere, C. A. 2017, https://doi.org/10.5281/zenodo.833792 [Google Scholar]
  79. Romeo, A. B., & Falstad, N. 2013, MNRAS, 433, 1389 [NASA ADS] [CrossRef] [Google Scholar]
  80. Romeo, A. B., & Fathi, K. 2016, MNRAS, 460, 2360 [NASA ADS] [CrossRef] [Google Scholar]
  81. Romeo, A. B., & Mogotsi, K. M. 2017, MNRAS, 469, 286 [NASA ADS] [CrossRef] [Google Scholar]
  82. Romeo, A. B., Agertz, O., & Renaud, F. 2020, MNRAS, 499, 5656 [NASA ADS] [CrossRef] [Google Scholar]
  83. Saintonge, A., Kauffmann, G., Kramer, C., et al. 2011, MNRAS, 415, 32 [NASA ADS] [CrossRef] [Google Scholar]
  84. Saintonge, A., Tacconi, L. J., Fabello, S., et al. 2012, ApJ, 758, 73 [Google Scholar]
  85. Saintonge, A., Catinella, B., Tacconi, L. J., et al. 2017, ApJS, 233, 22 [Google Scholar]
  86. Sani, E., Davies, R. I., Sternberg, A., et al. 2012, MNRAS, 424, 1963 [NASA ADS] [CrossRef] [Google Scholar]
  87. Shankar, F., Bernardi, M., Sheth, R. K., et al. 2016, MNRAS, 460, 3119 [NASA ADS] [CrossRef] [Google Scholar]
  88. Shankar, F., Bernardi, M., & Sheth, R. K. 2017, MNRAS, 466, 4029 [NASA ADS] [Google Scholar]
  89. Shimizu, T. T., Davies, R. I., Lutz, D., et al. 2019, MNRAS, 490, 5860 [Google Scholar]
  90. Shinozaki, K., Miyaji, T., Ishisaki, Y., Ueda, Y., & Ogasaka, Y. 2006, AJ, 131, 2843 [NASA ADS] [CrossRef] [Google Scholar]
  91. Silk, J., & Rees, M. J. 1998, A&A, 331, L1 [NASA ADS] [Google Scholar]
  92. Smith, R. N., Tombesi, F., Veilleux, S., Lohfink, A. M., & Luminari, A. 2019, ApJ, 887, 69 [NASA ADS] [CrossRef] [Google Scholar]
  93. Solomon, P. M., & Vanden Bout, P. A. 2005, ARA&A, 43, 677 [NASA ADS] [CrossRef] [Google Scholar]
  94. Tadhunter, C., Morganti, R., Rose, M., Oonk, J. B. R., & Oosterloo, T. 2014, Nature, 511, 440 [Google Scholar]
  95. Tombesi, F., Cappi, M., Reeves, J. N., et al. 2010, A&A, 521, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  96. Tombesi, F., Cappi, M., Reeves, J. N., et al. 2011, ApJ, 742, 44 [Google Scholar]
  97. Tombesi, F., Cappi, M., Reeves, J. N., & Braito, V. 2012, MNRAS, 422, L1 [Google Scholar]
  98. Tombesi, F., Meléndez, M., Veilleux, S., et al. 2015, Nature, 519, 436 [NASA ADS] [CrossRef] [Google Scholar]
  99. Toomre, A. 1964, ApJ, 139, 1217 [Google Scholar]
  100. van de Voort, F., Davis, T. A., Matsushita, S., et al. 2018, MNRAS, 476, 122 [NASA ADS] [CrossRef] [Google Scholar]
  101. Veilleux, S., Cecil, G., & Bland-Hawthorn, J. 2005, ARA&A, 43, 769 [NASA ADS] [CrossRef] [Google Scholar]
  102. Veilleux, S., Maiolino, R., Bolatto, A. D., & Aalto, S. 2020, A&ARv, 28, 2 [NASA ADS] [CrossRef] [Google Scholar]
  103. Venturi, G., Nardini, E., Marconi, A., et al. 2018, A&A, 619, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  104. Zubovas, K., & King, A. 2012, ApJ, 745, L34 [Google Scholar]

Appendix A: Stellar mass

The rest-frame SED is fitted with a three-component model using the approach described in Gruppioni et al. (2016). The adopted code combines three components simultaneously: the stellar component, the dusty star formation part, and the AGN–torus emission. The adopted libraries are the Bruzual & Charlot (2003) stellar library, the da Cunha et al. (2008) IR dust-emission library, and the library of AGN tori by Fritz et al. (2006), updated by Feltre et al. (2012). The latter includes both the emission of the dusty torus heated by the central AGN engine and the emission of the accretion disc. The considered AGN models are based on a continuous distribution of dust across the torus. The torus emission is shown as a green dashed line. The decomposition code provides a probability distribution function (PDF) in each photometric band and for each fit component, allowing an estimate of the uncertainty related to each decomposed contribution. In Fig. A.1 we show the observed SED of Mrk509 with the galaxy and AGN components. We defer to Gruppioni et al. (2016) for further details about the SED decomposition method. This SED decomposition delivers the star formation rate SFR = 5.1 ± 0.5 M⊙ yr−1 (Gruppioni et al. 2016), and the stellar mass M* = (1.2 ± 0.1)×1011 M⊙ (this work).

thumbnail Fig. A.1.

Mrk509 observed SED decomposed into stellar, AGN, and SF components using the technique described by Gruppioni et al. (2016). The black filled circles with error bars are data, the blue dotted lines show the unabsorbed stellar component, the red dashed lines show the combination of extinguished stars and dust IR emission, while the long-dashed green lines show the dusty torus emission. The pale blue dot-dashed lines show the dust re-emission, while the black solid lines are the sum of all the components (total emission).

Appendix B: Inner disc model

Based on the residual maps of the disc model (Section 3), we perform a dynamical modelling of the inner ∼ 1.5 kpc region, where the model of the large-scale disc produces significant blue- and red-shifted residuals (Fig. 2). This nuclear region is barely resolved with the current angular resolution, but we attempt to model it with a nuclear disc using the technique explained in Section 3. Given the limited signal-to-noise ratio and angular resolution we fix the inclination to 66 deg, derived from the minor/major axis ratio of the inner emitting region measured from the zeroth moment of the emission line data cube. We fix the PA value to 330 deg, across the two brightest emission regions in the residual intensity map (Fig. 5). We allow two parameters to vary: the rotational velocity vrot and the velocity dispersion σgas, with first guess values set to 80 and 50 km/s, respectively (Fig. 5, right panel). The moment 1 and 2 maps of this model are reported in the central panels of Fig. B.1. The right panel of Fig. 5 shows the position velocity diagrams cut along the kinematic major axis (330 deg), where the greyscale and blue contours refer to the data and the red contours to the model.

thumbnail Fig. B.1.

Same as Fig. 2 but for the inner disc model. Top panels: CO(2-1) moment 1 (top left) and moment 2 maps (top right) computed from the data cube, after subtraction of the disc model (Section 3). A threshold of 3σ has been applied. The black cross gives the AGN position. Central panels: Moment 1 (central left) and moment 2 (central right) maps of the disc model of the inner r = 700 pc region. Bottom panels: Moment 1 (bottom left) and moment 2 (bottom right) maps of the residuals obtained by subtracting the model velocity and velocity dispersion maps (central panels) from the observed ones (top panels).

All Tables

Table 1.

Parameters of the Gaussian fit of the CO(2−1) line.

In the text
Table 2.

Parameters of the molecular wind components.

In the text
Table 3.

Average Qgas and velocity dispersion in 4 annulii.

In the text

All Figures

thumbnail Fig. 1.

Continuum emission, CO(2−1) integrated flux maps of Mrk 509 and the CO(2−1) spectrum extracted from the continuum-subtracted clean data cube. Left panel: 1.2 mm continuum map of Mrk 509. Regions with emission below 5σ have been blanked. The synthesised beam (0.396 × 0.343 arcsec2) is shown by the grey filled ellipse. Central panel: CO(2−1) integrated flux map, where a mask with threshold of 3σ has been applied. Black contours show the 1.2 mm continuum emission at (5, 30, 100) × σ. The synthesised beam is shown in the bottom left corner. Right panel: CO(2−1) spectrum (black crosses) extracted from the continuum-subtracted clean data cube from the regions where emission is above 3σ. The solid red line shows the multi-Gaussian fit; the magenta lines show the three Gaussian components.
In the text
thumbnail Fig. 2.

Comparison between moment 1 and 2 maps generated from data, 3DBAROLO disc model and residuals. Top panels: CO(2−1) moment 1 (top left) and moment 2 (top right) maps. A threshold of 3σ has been applied. The black cross gives the AGN position. Central panels: moment 1 (centre left) and moment 2 (centre right) maps of the 3DBAROLO disc model. The black dashed lines give the orientation of the major and minor kinematic axes (PA = 248 and 338 deg, respectively). Bottom panels: moment 1 (bottom left) and moment 2 (bottom right) maps of the residuals obtained by subtracting the model velocity and velocity dispersion maps from the observed ones. A, B, C, and D indicate the main dynamical perturbations discussed in the text.
In the text
thumbnail Fig. 3.

Rotational velocity and dynamical mass obtained from the best-fit disc model of the main disc. Left panel: rotational velocity vrot vs. radius of the best-fit disc model of the main disc. Right panel: dynamical mass vs. radius within the inner ∼2.5 kpc radius, derived from the relation M dyn =r v rot 2 /2G $ M_{\rm dyn} = r{\it v}^2_{\rm rot}/2G $.
In the text
thumbnail Fig. 4.

Position-velocity diagrams along the kinematic major and minor axis. Left panel: PV plot along PA = 248 deg (the kinematic major axis). Right panel: PV plot along PA = 338 deg (the kinematic minor axis). The slit width is set to 0.45 arcsec (approximately the FWHM size of the synthetic beam major axis). The white contours refer to data, and are drawn at (1, 2, 4, 8, 16, 32) × σ; the black dots give the VLOS value from the 3DBAROLO model.
In the text
thumbnail Fig. 5.

Moment 1 map and PV diagram from the residual cube. Left panel: velocity map generated from the residual cube obtained subtracting the disc model from the data cube (a mask with a 2σ threshold on the residual intensity map was applied). The black contours show the CO(2−1) emission in residual intensity map at (3, 4, 5, 6, 7) × σ. The black dashed lines represent the slice adopted for the PV diagram (right panel) and comprise the two main kinematic perturbations in the residual intensity map, along PA = 330 deg. The black cross gives the AGN position. Right panel: PV diagram along PA = 330 deg (kinematic major axis). The red contours and yellow points represent the disc model of this inner region; the blue contours represent the data (contours are drawn at (2, 4, 6) × σ). A, C, and D indicate the main dynamical perturbations.
In the text
thumbnail Fig. 6.

Intensity map and spectra of regions A and B. Top panel: CO(2−1) residual intensity map after the best-fit disc model has been subtracted from the data cube. The red and blue circles indicate regions A and B respectively. The black contours show the emission at (2, 3, 4, 5) × σ. Bottom panels: spectra extracted from region A (red histogram) and B (blue histogram). The dashed black lines indicate the value of the disc rotational velocity at A and B based on our best-fit dynamical model (vdisc = 75 km s−1 around region A, and vdisc = 210 km s−1 around region B).
In the text
thumbnail Fig. 7.

HST [OIII]λ5007 Å image and Toomre Q-parameter map of Mrk 509. Right panel: HST image in the FQ508N filter, probing rest-frame [OIII]λ5007 Å emission (greyscale). The red contours represent CO(2−1) emission (Fig. 1centre panel, contours are drawn at (1, 2, 3, 4, 5, 7, 10, 15, 20, 25, 30, 37) × σ). The blue contours represent CO(2−1) emission in region A and B (Fig. 6 top panel, contours are drawn at (3, 4, 5) × σ). Image from the Mikulski Archive for Space Telescopes (MAST), HST Proposal 12212. Right panel: Toomre Q-parameter map for the cold molecular gas component. The black contours represent the moment 0 CO emission at (2, 5, 15, 30, 45) × σ. The white cross gives the continuum peak position. The white dashed ellipses indicate projected annuli with radius R = 0.3, 0.8, 1.4, 2.2 arcsec (see Table 3).
In the text
thumbnail Fig. 8.

Wind momentum load, Ṗof/Ṗrad, as a function of the wind velocity for Mrk 509 (Tombesi et al. 2011, 2012; Liu et al. 2015, and this work), and a compilation of local AGN (PDS456: Nardini et al. 2015; Bischetti et al. 2019, Mrk 231: Feruglio et al. 2015, and IRASF11119+3257: Tombesi et al. 2015). The red symbols are UFOs; the green star is the ionised [OIII] wind for Mrk 509; the blue symbols are molecular winds, and the magenta star is the total wind (ionised plus molecular phases) for Mrk 509. The black dashed line indicates the expectation for momentum-conserving wind driven by the AGN in Mrk 509. The red shaded area represents the prediction for an energy-conserving wind with Ṗmol/Ṗrad = vUFO/v.
In the text
thumbnail Fig. 9.

Wind mass outflow rate vs. AGN bolometric luminosity for Mrk 509 and a compilation of AGN from (Fiore et al. 2017). For Mrk 509 we plot the Ṁof of the total molecular wind (this work, light blue circle), the ionised wind (green circle Liu et al. 2015), and the UFO (Tombesi et al. 2010, 2011, 2012, red circle). The size of the circles represents the uncertainty in Ṁof. The small blue symbols represent molecular gas winds; the small green symbols represent ionised gas winds; and the small red symbols represent highly ionised gas winds detected in X-rays. The dashed blue, green, and red lines are the best-fit correlations of the molecular, ionised, and X-ray absorbers samples, respectively. Adapted from Fiore et al. (2017).
In the text
thumbnail Fig. A.1.

Mrk509 observed SED decomposed into stellar, AGN, and SF components using the technique described by Gruppioni et al. (2016). The black filled circles with error bars are data, the blue dotted lines show the unabsorbed stellar component, the red dashed lines show the combination of extinguished stars and dust IR emission, while the long-dashed green lines show the dusty torus emission. The pale blue dot-dashed lines show the dust re-emission, while the black solid lines are the sum of all the components (total emission).
In the text
thumbnail Fig. B.1.

Same as Fig. 2 but for the inner disc model. Top panels: CO(2-1) moment 1 (top left) and moment 2 maps (top right) computed from the data cube, after subtraction of the disc model (Section 3). A threshold of 3σ has been applied. The black cross gives the AGN position. Central panels: Moment 1 (central left) and moment 2 (central right) maps of the disc model of the inner r = 700 pc region. Bottom panels: Moment 1 (bottom left) and moment 2 (bottom right) maps of the residuals obtained by subtracting the model velocity and velocity dispersion maps (central panels) from the observed ones (top panels).
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.