© ESO 2022

1. Introduction

The study of the processes that power star formation (SF) in galaxies is paramount to understanding how galaxies form and evolve. If we assume that the gas scale-height is constant, the power law relating the gas volume density and the star-formation rate (SFR) volume density, originally proposed by Schmidt (1959), finds its equivalent in terms of the corresponding surface densities of the SFR (Σ_SFR) and the gas (Σ_gas) in the expression

$$\begin{array}{c}\hfill {\mathrm{\Sigma }}_{\mathrm{SFR}}=A{\mathrm{\Sigma }}_{\mathrm{gas}}^N.\end{array}$$(1)

In this equation A is a normalisation constant and N is the power-law index. Under the hypothesis that the relevant timescale for SF is the local free-fall time for the gas, the theoretically predicted value for N is 1.5. This relation is known as the Kennicutt–Schmidt (KS) law (Schmidt 1959; Kennicutt 1998a). In addition to the empirical power-law relation between Σ_gas and Σ_SFR, another key parameter in SF studies is the star-formation efficiency (SFE), defined as

$$\begin{array}{c}\hfill \mathrm{SFE}={\mathrm{\Sigma }}_{\mathrm{SFR}}/{\mathrm{\Sigma }}_{\mathrm{gas}},\end{array}$$(2)

which represents the inverse of the depletion time (T_dep) of the gas that is consumed by SF.

Different proxies for Σ_SFR and Σ_gas have been chosen for observations carried out during the last few decades on different galaxy populations. First, from spatially unresolved galaxy-scale global measurements, which used CO and HI as neutral gas tracers in different galaxy samples, observers found a single law with a range of indexes ≃1.2 − 1.7 (Kennicutt 1998a; Yao et al. 2003; Bouché et al. 2007; Daddi et al. 2010; Genzel et al. 2010; Liu et al. 2015; Kennicutt & De Los Reyes 2021). However, there is mounting evidence that global KS laws show signs of bimodality or multi-modality, which tend to separate the branches of normal galaxies from that of more extreme merger systems (Daddi et al. 2010; Genzel et al. 2010; Liu et al. 2015; Kennicutt & De Los Reyes 2021).

Furthermore, high-resolution imaging of neutral gas in galaxies has allowed the analysis of the KS relation at kiloparsec and sub-kiloparsec scales in a growing number of galaxies (e.g. Kennicutt et al. 2007; Bigiel et al. 2008; Leroy et al. 2008, 2017; Blanc et al. 2009; Casasola et al. 2015). In particular, Bigiel et al. (2008) used CO(2–1) as a tracer of molecular gas in a sample of 18 nearby star forming galaxies (SFGs) and found a single linear KS relation (N ≃ 1), which holds for gas surface densities Σ_gas ≥ 10 M_⊙ pc⁻². This dividing line identifies the transition from atomic to molecular gas. Leroy et al. (2008) derived a radial dependence in N with a decreasing SFE at larger radii within individual galaxies. Moreover, Casasola et al. (2015) derived KS relations in the nuclear regions of four low luminosity active galactic nuclei (AGN) using interferometric CO images on spatial scales between 20 pc and 200 pc. The KS relations were found to be sub-linear, but also super-linear, with a wide range of slopes, ∼[0.5 − 1.3]. Leroy et al. (2017) studied the local dynamical state of molecular gas in M 51 and found that the gas with stronger self-gravity forms stars at a higher rate. The variability in resolved KS relations also suggests a higher SFE in lower mass, low metallicity galaxies (e.g. Schruba et al. 2011; Leroy et al. 2013) and in late Hubble types (e.g. Colombo et al. 2018; Ellison et al. 2021).

The high-resolution CO observations of M 33 published by Onodera et al. (2010) show that at spatial scales similar to those of giant molecular clouds (GMCs; ∼50 − 100 pc), the correlation between Σ_SFR and Σ_gas is lost. In qualitative agreement with this picture, Schruba et al. (2010) observed a breakdown of the SF relation at scales of ≤300 pc in M 33. Moreover, the recent work of Williams et al. (2018) found significant correlations in M 33 down to scales of 100 pc, while the measured Schmidt index shows a marked dependence on the spatial scale. The breakdown of the KS relation observed below a ‘critical’ spatial scale can be attributed to the need of averaging over sufficiently large scales in order to have a statistical sampling of star forming sites at different evolutionary stages (e.g. Kruijssen & Longmore 2014).

The most recent studies of the SF relation reported above, which use low-J CO lines (sometimes in combination with HI data) as tracers of the bulk of neutral gas, cast doubts on the existence of a ‘universal’ or ‘uni-modal’ KS relation at all spatial scales. Dense molecular gas, namely gas with volume densities typically exceeding 10^4 − 5 cm⁻³, is believed to condense into GMCs and be therefore more directly related to recent and massive SF (Lada et al. 2010; André et al. 2010). In particular, observations of HCN(1–0) and HCO⁺(1–0) lines, which have associated critical densities of n_crit[HCN(1–0)] ∼ 1.7 × 10⁵ cm⁻³ and n_crit[HCO⁺(1–0)] ∼ 2.9 × 10⁴ cm⁻³ (Shirley 2015), are well suited for fairly tracing the dense molecular gas mass (M_dense) in galaxies. Dense gas probes have been used to study ‘galaxy-scale’ KS relations (e.g. Gao & Solomon 2004a,b; Solomon & Vanden Bout 2005; Graciá-Carpio et al. 2006, 2008; García-Burillo et al. 2012; Liu et al. 2015). As expected, the SF relations derived for the dense gas show a less scattered linear correlation (i.e. with N ≃ 1) compared to the global KS laws obtained from low-J CO lines. However, the residual but nevertheless significant scatter present in the SFR–M_dense plane has been interpreted as indicative of different average physical properties of the dense gas in normal SFGs and mergers (Graciá-Carpio et al. 2008; García-Burillo et al. 2012).

High-resolution (kiloparsec and sub-kiloparsec) single-galaxy studies of the SF relations of the dense molecular gas have started to resolve the degeneracy in the SFR–M_dense parameter space by showing how scaling laws change for different dynamical environments within a galaxy, including the Milky Way (Longmore et al. 2013; Kruijssen et al. 2014; Murphy et al. 2015; Usero et al. 2015; Bigiel et al. 2015, 2016; Chen et al. 2017; Viaene et al. 2018; Querejeta et al. 2019; Jiménez-Donaire et al. 2019; Bešlić et al. 2021). In particular, Usero et al. (2015) observed HCN(1–0) and CO(1–0) lines at several positions in the disks of 29 SFGs and found that SFE_dense, derived from the IR/HCN ratio, is ∼6 − 8 times lower near galaxy centres than in the outer regions of the disks. Furthermore, Querejeta et al. (2019) found that SFE_dense values measured on ∼100 pc scales from radio continuum-to-HCN ratios vary by more than 1 dex among the different dynamical environments of the disk of M 51. More recently, Bešlić et al. (2021) found significant differences in the ∼100 pc-scale SFE_dense values between the galaxy centre, bar, and bar-end regions of the nearby barred galaxy NGC 3627. These results contradict models that rely on a universal gas density threshold for SF (e.g. Gao & Solomon 2004b; Wu et al. 2005; Lada et al. 2010, 2012; Evans et al. 2014) and suggest instead that the dynamical environment of the dense molecular gas in GMCs can determine its efficiency at forming stars. This is supported by models of turbulent SF (e.g. Krumholz & McKee 2005; Krumholz & Thompson 2007; Hennebelle & Falgarone 2012; Meidt et al. 2013, 2018, 2020; Federrath 2015; Meidt 2016).

In this paper we study the spatially resolved SF relations of the dense molecular gas in the starburst (SB) ring of the nearby (D ≃ 14 Mpc; Bland-Hawthorn et al. 1997) Seyfert 2 barred galaxy NGC 1068, a target considered as an archetype of the composite SB and AGN classification. Previous interferometer images have resolved the large-scale distribution of molecular gas in the disk of the galaxy (Helfer & Blitz 1995; Schinnerer et al. 2000; Krips et al. 2011; Tsai et al. 2012; García-Burillo et al. 2014, 2017, 2019; Takano et al. 2014; Viti et al. 2014; Scourfield et al. 2020). Molecular line and dust continuum emissions are detected from a r ∼ 200 pc circumnuclear disk (CND), from the ∼2.6 kpc diameter stellar bar region, and from a ring where molecular gas is accumulating and feeding a SB episode. The SB ring is formed by a tightly wound two-arm spiral structure that starts from the ends of the stellar bar and unfolds in the disk over ∼180° in azimuth, forming a pseudo-ring at r ∼ 18″(1.3 kpc).

We used new images of the distribution of dense molecular gas (n_H2 ≥ 10^4 − 5 cm⁻³) obtained by the Atacama Large Millimeter Array (ALMA) in the 1−0 transitions of HCN and HCO⁺ with a native resolution of $\sim 0\stackrel{″}{.}8$ (56 pc). This spatial resolution is comparable to the typical size of GMCs. We also use high-resolution ($\sim 1\stackrel{″}{.}5=100$ pc) CO (1–0) images of the galaxy obtained by the IRAM array (Schinnerer et al. 2000), as well as available CO(3–2) and dust continuum images obtained by ALMA at a spatial resolution of $\sim 0\stackrel{″}{.}6$ (40 pc) (García-Burillo et al. 2014). The ensemble of these observations allows us to probe a wide range of molecular gas densities (n_H2 ≃ 10^3 − 5 cm⁻³) in the SB ring. To probe SF, we use Paα line images obtained by the Hubble Space Telescope (HST). We analyse how SF relations change depending on the adopted spatial resolution, on the choice of molecular gas tracer, and on the particular dynamical environment throughout the SB ring.

The paper is organised as follows: Sect. 2 presents the new ALMA observations and accompanying ancillary data. We describe in Sect. 3 the conversion factors adopted. Section 4 describes the molecular gas and Paα images used in this work. We study the different KS relations derived in NGC 1068 in Sect. 5. Section 6 explores a different prescription of SF relations and analyses the environmental dependence of SFE_dense as a function of a set of physical parameters in the different regions of the SB ring. We describe a scenario for the SF in the SB ring in Sect. 7. The main conclusions of this work are summarised in Sect. 8.

2. Observations

In this section we present an overview of the different datasets used to probe the distribution of molecular gas (Sect. 2.1) and the recent SF (Sect. 2.2) required to derive the SF relations in the SB ring of NGC 1068.

2.1. Molecular gas tracers

2.1.1. New ALMA data

We used ALMA to map the emission of HCN(1–0) and HCO⁺(1–0) in the central r ∼ 2.5 kpc of the NGC 1068 disk. Observations were executed during Cycle 2 in one track in August 2015 (project-ID: # 2013.1.00055.S, PI: S. García-Burillo). We used band 3 receivers and a single pointing with a field of view (FOV) of ∼70″ (∼5 kpc), covering the CND and the SB ring of the galaxy. Observations made use of 34 antennas of the array with projected baselines ranging from 12 m to 1430 m. The phase tracking centre was set to α₂₀₀₀ = 02^h42^m40.771^s, δ₂₀₀₀ = −00° 00′47.84″, which is the centre of the galaxy according to SIMBAD taken from the Two Micron All Sky Survey (2MASS; Skrutskie et al. 2006). The tracking centre is offset by ≤1″ relative to the AGN position: α₂₀₀₀ = 02^h42^m40.71^s, δ₂₀₀₀ = −00° 00′47.94″ (Gallimore et al. 1996, 2016, 2004; García-Burillo et al. 2014, 2016; Imanishi et al. 2016). The galaxy has a systemic velocity of v_sys(HEL) ∼ 1130 km s⁻¹ (García-Burillo et al. 2014, 2019).

Four spectral windows were placed, two in the lower sideband and two in the upper sideband. All the sub-bands have a spectral bandwidth of 1.875 GHz. The setup allowed us to simultaneously observe HCN(J = 1 − 0) (88.632 GHz at rest) and HCO⁺(J = 1 − 0) (89.189 GHz at rest) in the higher frequency lower sideband, as well as H¹³CN(J = 1 − 0) (86.340 GHz at rest) and H¹³CO⁺(J = 1 − 0) (86.754 GHz at rest) in the lower frequency lower sideband. The two spectral windows in the upper sideband were centred around the CS(J = 2 − 1) (97.981 GHz at rest) line and the continuum emission around 100 GHz, respectively. The CS(2–1) map was published by Scourfield et al. (2020). The 86.6 GHz-continuum map of the galaxy was published by García-Burillo et al. (2017).

The data were calibrated using the ALMA reduction package CASA¹. The calibrated uv-tables were exported to GILDAS²-readable format (Guilloteau & Lucas 2000) in order to perform the mapping and cleaning steps as detailed below. We estimated that the absolute flux accuracy is about 5%, which is in line with the goal of standard ALMA observations at these frequencies. The synthesised beam obtained using natural weighting is $1^{″}\times 0\stackrel{″}{.}6$ (70 pc × 42 pc) at a position angle PA = 69°. The line data cube was binned to a frequency resolution of 2.92 MHz (∼10 km s⁻¹). We estimated a 1σ sensitivity of 0.4 mJy beam⁻¹ per channel of 10 km s⁻¹ using line-free emission areas in the data. The conversion factor between Jy beam⁻¹ and K is 247.6 K Jy⁻¹ per beam in both emission lines. The spectral line maps were obtained after subtraction of the continuum emission performed in the (u, v) plane using the GILDAS tasks uv-average and uv-subtract.

We obtained the zeroth, first and second moment maps from the line data cubes using the GILDAS task moments, adopting a velocity window |v − v_sys| = 250 km s⁻¹, which is enough to cover the span of velocities due to rotation in the disk of the galaxy (García-Burillo et al. 2014). The total uncertainty on the velocity-integrated emission maps, ΔI, was derived from

$$\begin{array}{c}\hfill \mathrm{\Delta }I=\sqrt{\mathrm{\Delta }{I}_{\mathrm{noise}}^2+\mathrm{\Delta }{I}_{\mathrm{calib}}^2}\end{array}$$(3)

where

$$\begin{array}{c}\hfill \mathrm{\Delta }{I}_{\mathrm{noise}}=\sigma \times \delta V\times \sqrt{N_{\mathrm{window}}}.\end{array}$$(4)

Here, ΔI_noise is the velocity-integrated intensity error, which results from propagating the error of individual channels, σ = 0.4 mJy beam⁻¹, of width δV = 10 km s⁻¹, to the number of channels considered in the integration window, N_window = 50. Furthermore, ΔI_calib is the uncertainty due to the absolute flux calibration error, which is about ≃5%×I for ALMA band 3 observations (see ALMA Technical Handbook³). All in all, for a typical flux integrated value ≥5σ characteristic of the regions in the SB ring of NGC 1068 studied in this work, the total uncertainty on I as well as on all the related parameters, namely line luminosities and gas masses, amounts at most to ±0.09 dex in logarithmic units⁴.

The largest angular scale (LAS) of our observations is ∼4″ (∼280 pc). Since our observations do not contain short-spacing correction, the flux can start to be filtered out on scales larger than the LAS. The HCN emission, and very likely also the HCO⁺ emission, are expected to arise from a highly clumpy medium consisting of an ensemble of dense cloud cores. This particular hierarchy of the dense molecular gas probed by HCN and HCO⁺ helped by the velocity structure observed in the molecular disk of NGC 1068 (e.g. see García-Burillo et al. 2014) allows us to foresee that the amount of flux filtered on the spatial scales that are the most relevant for this paper is kept low in both lines.

To validate this hypothesis we estimated the maximum percentage of missing flux in the HCN image of NGC 1068 by comparing the fluxes measured by ALMA and by the IRAM 30-m telescope (Usero et al. 2021, priv. comm.) at different locations of the disk. With this aim we derived the spatially integrated fluxes using the single-dish aperture sizes ∼30″ ∼ 2 kpc (see Appendix A for details). The result of this comparison indicates that the maximum percentage of missing flux in the HCN ALMA map is about 25% on scales of 2 kpc. As the spatial scales relevant for this paper are much smaller (40−700 pc), we can conclude that 25% is a conservative upper limit on the missing flux for HCN (and very likely also for HCO⁺).

2.1.2. Ancillary data

We used the CO(3–2) line and 349 GHz (859 μm) continuum emission images of the galaxy obtained by García-Burillo et al. (2014) with ALMA during the Cycle 0 of the array in band 7 (project-ID: # 2011.0.00083.S, PI: S. García-Burillo). The angular resolution of these data is $\sim 0\stackrel{″}{.}6\times 0\stackrel{″}{.}5$ (42 pc × 35 pc) at a position angle of ∼60°. We refer to García-Burillo et al. (2014) for a detailed description of the data reduction steps. As the observations of García-Burillo et al. (2014) do not contain short-spacing correction, we expect that a non-negligible amount of flux may start to be filtered out on scales beyond the reported LAS ≥6″ (420 pc) for the continuum and CO(3–2) emission images. Based on a comparison between the fluxes measured by ALMA and different single-dish telescopes using a set of apertures, García-Burillo et al. (2014) estimated that their interferometer images may be filtering up to 20−30% and 65% of the total flux on spatial scales of about 1 kpc for the CO(3–2) and continuum emission, respectively. However, the clumpy distribution of the gas and also (in the case of the CO line) the velocity structure of the emission are expected to favour the recovery of most of the flux in the line and continuum maps on smaller apertures (≤6″ ∼ 420 pc) centred on the brightest emission spots of the SB ring.

We also used the CO(1–0) line map of the galaxy obtained by the IRAM array on the Plateau de Bure Interferometer (PdBI), published by Schinnerer et al. (2000). The CO(1–0) line allows us to study the bulk of the molecular gas reservoir in the disk of NGC 1068. The angular resolution of these data is $\sim 1\stackrel{″}{.}8\times 1\stackrel{″}{.}0$ (126 pc × 70 pc) at a position angle of ∼24°. Schinnerer et al. (2000) estimated that the CO PdBI map misses about 25 − 30% of the total flux on scales ∼55″ (3.8 kpc), based on the comparison between the CO flux measured by the PdBI and the flux derived from the CO maps obtained through the combination of the Berkeley-Illinois-Maryland-Association (BIMA) array and the 12 m Kitt Peak single-dish data published by Helfer & Blitz (1995). We derived a new upper limit on the missing flux in the CO(1–0) map based on a comparison between the fluxes measured by PdBI and the IRAM 30-m telescope in Appendix A (∼40−45% on scales ∼24″ ∼ 1.7 kpc). As in this paper the relevant spatial scales used in our analysis are smaller, we can therefore expect that the missing flux factors reported above can be taken as strict upper limits.

2.2. Star-formation tracer

We used the emission of the Paα hydrogen recombination line at 1.875 μm to image the distribution of recent SF in the disk of NGC 1068. We used the HST/NICMOS (NIC3) narrow-band (F187N, F190N) images of the galaxy retrieved from the Hubble Legacy Archive (HLA)⁵ to derive the continuum-subtracted Paα map, and followed the calibration and continuum subtraction steps detailed in Sect. 2.2 of García-Burillo et al. (2014). The pixel size of the HLA images is $0\stackrel{″}{.}1\times 0\stackrel{″}{.}1$. The angular resolution (full width at half maximum) of the Paα image is $0\stackrel{″}{.}26\times 0\stackrel{″}{.}26$ (∼18 pc × 18 pc), as determined from the estimated size of the point spread function in the observations. Uncertainties on flux calibration are at the 15−20% level (Böker et al. 1999; Alonso-Herrero et al. 2006), which implies associated uncertainties ≤ ± 0.09 − 0.10 dex for the regions of the SB ring of NGC 1068 examined in this work. The Paα line traces ionised gas produced by associations of massive (≳30 M_⊙) and young (≲8 − 10 Myr) stars. The main advantage of the near-infrared (NIR) recombination line compared to its optical counterpart, namely Hα, resides in the significantly lower extinction by dust of Paα (Kennicutt 1998b; Calzetti et al. 2007). We may therefore neglect any dust-extinction correction when we derive the SFR from the Paα fluxes. The validity of this hypothesis is examined in Appendix B, where we compare the Paα fluxes measured by HST over a number of hot spots of the SB ring with those measured in Hα using the ground-based image of the galaxy published by Díaz et al. (2000). We estimate an overall low extinction correction at 1.875 μm for the SB ring knots: A_Paα shows a median value ∼0.03 mag, compatible with optically thin emission (see Appendix B for details). If we allow for a ∼20% uncertainty in the flux scales due to absolute calibration errors, we conclude that the Paα map of the SB ring does not require any significant correction for dust extinction.

3. Conversion to physical parameters

3.1. Molecular gas masses

To derive the distribution of ‘local’ dense molecular gas mass (M_dense) from the HCN(1–0) velocity-integrated luminosities ($L_{\text{dense}}^{\prime }$) we assumed the standard conversion factor commonly applied for ‘global’ scales, α_HCN = 10 M_⊙ (K km s⁻¹ pc²)⁻¹ for HCN, following Gao & Solomon (2004a), for the different local scales used in this work. This factor includes a correction for helium. A canonical ‘constant’ conversion factor is usually adopted in the literature for both lines, which are considered as reliable tracers of the dense molecular gas phase above densities n ≈ 10⁴ cm⁻³ (see, however, García-Burillo et al. 2012; Evans et al. 2020).

This fixed conversion factor assumes that the HCN(1–0) emission is originated from gravitationally-bound ‘cores’ or clumps with volume-averaged density n_H2 ∼ 3 × 10⁴ cm⁻³ and a brightness temperature T_b ∼ 35 K. However, under the hypothesis that the emission of the HCN line is mostly optically thick, if the volume-averaged density of the emitting clumps is lower than 3 × 10⁴ cm⁻³ or if the T_b is larger than 35 K, the α_HCN could be smaller than the factor suggested by Gao & Solomon (2004a). In this context, Wu et al. (2005) estimated the value of α_HCN in Galactic star-forming cores, finding a slightly lower conversion factor at smaller scales: α_HCN = 7 ± 2 M_⊙ (K km s⁻¹ pc²)⁻¹. This value differs only by 30% from the global conversion factor used by Gao & Solomon (2004a), which is the one adopted in this work. We nevertheless note that adopting a lower value of the conversion factor for HCN would result in slightly lower molecular gas surface density values, particularly in the regions of the SB ring that show comparatively higher SF activity and SFE values (see discussion in Sect. 6).

We therefore derived M_dense as

$$\begin{array}{c}\hfill M_{\mathrm{dense}}[M_{\odot }]={\alpha }_{\mathrm{HCN}}{L}_{\mathrm{dense}}^{\prime }[\mathrm{K}\mathrm{km}{\mathrm{s}}^{-1}{\mathrm{pc}}^2].\end{array}$$(5)

The line luminosity, L′, is defined following Solomon et al. (1997) as

$$\begin{array}{c}\hfill L_{\mathrm{dense}}^{\prime }[\mathrm{K}\mathrm{km}{\mathrm{s}}^{-1}{\mathrm{pc}}^2]=3.25\times {10}^7\times S_{\mathrm{dense}}\mathrm{\Delta }V{\nu }_{\mathrm{obs}}^{-2}{D}_{\mathrm{L}}^2{(1+z)}^{-3},\end{array}$$(6)

where the velocity-integrated fluxes S_denseΔ V are in Jy km s⁻¹ particularised for each line, the observed frequency ν_obs is in GHz and the luminosity distance D_L is in Mpc units.

We obtained face-on values of the dense molecular gas surface densities (Σ_dense in M_⊙ pc⁻² units) from

$$\begin{array}{c}\hfill {\mathrm{\Sigma }}_{\mathrm{dense}}=\frac{M_{\mathrm{dense}}}{A_{\mathrm{aperture}}}\times cos(i),\end{array}$$(7)

where i  =  40° is the inclination of the disk of NGC 1068 (Bland-Hawthorn et al. 1997; Brinks et al. 1997; García-Burillo et al. 2014) and A_aperture is the area of the aperture used in pc².

Similarly, we transformed the measured CO(1–0) luminosities into molecular gas masses and surface densities using the conversion prescription of Bolatto et al. (2013), which assumes a standard Galactic conversion factor α_CO = 4.4 M_⊙ (K km s⁻¹ pc²)⁻¹, which already includes a correction for helium.

We used the GILDAS task gauss-smooth to convolve the initial resolution versions of the molecular line data cubes with the appropriate Gaussian kernels adapted to generate all the image versions for the common set of spatial resolutions used in this work, which range from ∼40 − 56 pc (∼100 pc for CO(1–0)) up to ∼700 pc (see Sect. 5.1).

3.2. Star-formation rates

We adopted the prescription proposed in Kennicutt & Evans (2012) regarding the conversion factor used to calculate the local SFR map from the Paα line luminosities. In particular, we assumed a Kroupa initial mass function (Kroupa 2001) as well as an intrinsic ratio for Hα/Paα ∼ 7.81 (Hummer & Storey 1987), which applies for the case B recombination at T_e = 5000 K and n_e = 10³ cm⁻³. These conditions are found in starbursting galaxies (Roy et al. 2008; Rieke et al. 2009).

We obtained the SFR from the expression

$$\begin{array}{c}\hfill \mathrm{SFR}(M_{\odot }{\mathrm{yr}}^{-1})=4.15\times {10}^{-41}L(\mathrm{Pa}\alpha ,{\mathrm{erg}\mathrm{s}}^{-1})\end{array}$$(8)

We derived the corresponding SFR surface densities (Σ_SFR) in units of M_⊙ yr⁻¹ pc⁻² from

$$\begin{array}{c}\hfill {\mathrm{\Sigma }}_{\mathrm{SFR}}=\frac{\mathrm{SFR}}{A_{\mathrm{aperture}}}\times cos(i).\end{array}$$(9)

We estimated the SFR for the range of spatial resolutions analysed in this work following the same procedure described in Sect. 3.1, which uses the GILDAS task gauss_smooth to convolve the initial resolution images with the appropriate Gaussian kernels.

4. Dense molecular gas and SF maps

4.1. The HCN and HCO⁺ maps

Figure 1 shows the HCN(1–0) and HCO⁺(1–0) velocity-integrated intensity maps of NGC 1068 obtained by ALMA in the central r ∼ 2.5 kpc of the disk. Overall, the distribution of dense molecular gas is similar to that shown by other molecular gas tracers as seen in previous interferometer images of the galaxy (Schinnerer et al. 2000; Krips et al. 2011; Tsai et al. 2012; García-Burillo et al. 2014, 2019; Viti et al. 2014; Scourfield et al. 2020). In particular, the bulk of the HCN(1–0) and HCO⁺(1–0) emission stems from three main regions.

Fig. 1. HCN(1–0) and HCO+(1–0) integrated intensity maps of NGC 1068. Left panel: velocity-integrated intensity map of HCN(1–0) of NGC 1068 obtained with ALMA. The map is shown in a logarithmic colour scale with contour levels of 3σ, 5σ, 7σ, 12σ, 24σ, 50σ, 120σ, and 200σ, where 1σ = 0.028 Jy beam−1 km s−1. We highlight the location of several representative regions of the emission, namely the CND, the bar, and the SB ring. Right panel: same as the left panel but for the HCO+(1–0) line map, with the same list of contour levels but here truncated at 120σ. The position of the AGN is identified by the dark grey star marker. The black bar in the top-right corner of each panel shows the spatial scale (in pc). The magenta filled ellipses at the bottom-left corners of the panels represent the ALMA beam size (1$\stackrel{″}{.}$0 × 0$\stackrel{″}{.}$6 at PA = 69°).

The first is the CND. Described as an asymmetric elliptical ring of 6″ × 4″–size (r ∼ 200 pc), the CND shows two emission knots located ∼1″ east and ∼1.5″ west of the AGN. The CND ring is off-centred relative to the AGN locus. The two emission knots are bridged by weaker emission north and south of the AGN. The morphology of the HCN and HCO⁺ maps of the CND is to a large extent similar to that of the ALMA CO and CS maps (García-Burillo et al. 2014, 2019; Scourfield et al. 2020).

The second is the bar. There is HCN and HCO⁺ emission in the region occupied by the ∼2.6 kpc diameter stellar bar, which is oriented along PA = 46° ±2° (Scoville et al. 1988; Schinnerer et al. 2000). As shown by other molecular gas tracers, the emission in this region from both lines (especially for HCN), tends to accumulate along the leading edges of the bar. We also detect significant emission in the ‘bow-shock arc’ feature identified in the CO(3–2) and continuum dust emission maps of García-Burillo et al. (2014) on the north-east side of the disk at r ∼ 4″ − 7″ (300 pc–500 pc).

The third is the SB ring. Most of the dense molecular gas in the disk concentrates in a ring of r ∼ 18″(1.3 kpc) formed by two tightly wound spiral arms, which unfold over ∼180° in azimuth in the disk from the ends of the stellar bar. The SB ring concentrates also most of the massive star forming complexes in the disk identified in the Paα image (see Sect 4.2). The emission of HCN and HCO⁺ is unevenly distributed azimuthally over the SB ring: in both lines the emission is strongest around two regions, located at r ∼ 18″ and PA ∼ 15° −75° ( ± 180°), where the ring is connected to the stellar bar ends. The emission of HCN and HCO⁺ in the SB ring is clumpy and it appears to be organised as coming from molecular cloud associations of ≥50 pc size. The SB ring is connected at larger radii to two emission lanes located at the edge of the HCN and HCO⁺ maps shown in Fig. 1 along PA ∼ 50° ±180°.

Figure 2 shows the HCN-to-HCO⁺ line brightness ratio (R_dense) in the disk of the galaxy. The R_dense map was obtained assuming a common 3σ threshold on the integrated intensities of both lines. The R_dense ratio changes significantly across the different regions of the disk identified above. In particular, R_dense ≃ 1.5 − 3.5 in the CND and the bow-shock arc region, for which we estimated a mean value ⟨R_dense⟩≃2.2. On the other hand, R_dense ≃ 0.5 − 1.5 in the SB ring and the corresponding mean value ⟨R_dense⟩≃1.1.

Fig. 2. HCN(1–0)/HCO+(1–0) brightness temperature ratio map (Rdense) derived at the common spatial resolution of ≃56 pc. The brown circle of 8″ in radius (≃560 pc) locates the inner region of the galaxy disk where the molecular outflow signature has been identified in the kinematics of molecular gas (García-Burillo et al. 2014, 2019). This region, purposely excluded from our analysis of SF relations, shows significantly higher line ratios (Rdense ≃ 1.5 − 3.5) compared to the SB ring (Rdense ≃ 0.5 − 1.5). Other symbols are the same as in Fig. 1.

The high R_dense values measured in the CND and the bow-shock arc are related to the molecular outflow signature identified in the kinematics of molecular gas in these regions (García-Burillo et al. 2014, 2019). The outflow is thought to be driven by the interaction of the AGN wind and the radio jet with the molecular gas in the disk in a mostly coplanar geometry. Besides leaving a distinct kinematic signature in the CND and the bow-shock arc, the outflow has left its imprint on the excitation and the chemistry of molecular gas, which is under the influence of large-scale shocks and a strong UV irradiation in these regions (Viti et al. 2014; García-Burillo et al. 2017).

4.2. The Paα map

Figure 3 shows the Paα image of NGC 1068 obtained at the initial resolution of the HST/NICMOS camera: $0\stackrel{″}{.}26\times 0\stackrel{″}{.}26$ (∼18 pc × 18 pc). The central region of the image reveals strong emission stemming from an asymmetric bipolar nebula of ionised gas. Both the morphology and the kinematics of the gas in this structure have been modelled in terms of an AGN-driven wind (e.g. Crenshaw & Kraemer 2000; Cecil et al. 2002; Das et al. 2006; Müller-Sánchez et al. 2011; Barbosa et al. 2014; Miyauchi & Kishimoto 2020). The AGN wind occupies a hollow bicone that extends up to a radius r ∼ 8″ (550 pc) on its northern side. The bicone feature is oriented along PA ∼ 30°, and is characterised by a wide opening angle (∼80°). In our subsequent analysis of the SF relations we therefore screen the central r ≤ 8″ region of the galaxy, where Paα cannot be considered as a reliable tracer of SF.

Fig. 3. HST/NICMOS continuum-subtracted Paα image of NGC 1068. The map is shown in colour scale with contour levels of 3σ, 9σ, 15σ, 32σ, 50σ, 120σ, and 700σ, where 1σ = 0.212 × 10−16 erg s−1 cm−2 pixel−1. The colour scale range is shown in units of 10−16 erg s−1 cm−2 pixel−1. The (truncated) square region identifies the outer edge of the HST/NICMOS FOV. The angular resolution of the image is 0$\stackrel{″}{.}$26 × 0$\stackrel{″}{.}$26. Other symbols are the same as in Fig. 2.

Outside the bright AGN bicone structure, most of the Paα emission comes from the SB ring. Similar to the distribution of the dense molecular gas, SF traced by Paα is not uniformly distributed throughout the SB ring. As for HCN and HCO⁺, the brightest SF complexes are located at the north-east section, and most particularly, at the south-west section of the SB ring (at PA ∼ 195° −255°). In either case these are the two regions where the ring is connected to the ends of the stellar bar. Rico-Villas et al. (2021) studied the 147 GHz free-free emission associated with SF and identified a similar concentration of massive super star clusters in the bar-ring interface region.

Figure 4 shows the overlay of the HST/NICMOS Paα emission image on the HCN(1–0), HCO⁺(1–0), and CO(1–0) maps. For a proper comparison we degraded the Paα image to the common spatial resolution of HCN and HCO⁺ (∼56 pc), and to that of CO (∼100 pc). A visual inspection of this figure illustrates that the Paα maxima do not always coincide with the strongest emission peaks in HCN, HCO⁺ or CO throughout the SB ring.

Fig. 4. Overlay of the HST/NICMOS Paα image on the molecular gas tracer maps. Left panel: overlay of the HST/NICMOS Paα image (contours) on the ALMA HCN(1–0) map (colour scale). Contours have a logarithmic spacing from 3σ to 1536σ, in steps of 0.3 dex, where 1σ = 0.180−16 erg s−1 cm−2 pixel−1. The colour scale is shown in units of Jy beam−1 km s−1. Middle panel: same as the left panel but showing the overlay of the HST/NICMOS Paα contours on the ALMA HCO+(1–0) map (colour scale). Right panel: same as the left panel but showing the overlay of the HST/NICMOS Paα contours on the PdBI CO(1–0) map of Schinnerer et al. (2000) (colour scale). Contour spacing is the same as in the left and right panels, but with 1σ = 0.162 × 10−16 erg s−1 cm−2 pixel−1. The magenta filled ellipses at the bottom-left corners in all panels represent the beam sizes of the molecular line maps: 1$\stackrel{″}{.}$0 × 0$\stackrel{″}{.}$6 at PA = 69° (left and middle panels) and 1$\stackrel{″}{.}$84 × 1$\stackrel{″}{.}$09 at PA = 24° (right panel). The HST/NICMOS Paα images in each panel have been degraded to the corresponding spatial resolutions of the molecular gas tracers used in this comparison. Other symbols are the same as in Fig. 2.

5. Star-formation relations in NGC 1068

5.1. Kennicutt–Schmidt laws

We use in this section the images described in Sects. 2.1 and 2.2 to obtain different versions of the pixel-wise spatially resolved KS relations in the SB ring of NGC 1068 for a set of seven spatial resolutions, ranging from ∼40 − 56 pc (∼100 pc for CO(1–0)) up to ∼700 pc. To derive Σ_SFR and Σ_gas for each spatial scale we degraded our datasets to the selected resolutions using the GILDAS task gauss-smooth. In our analysis we only consider pixels with a signal ≥3σ in Σ_SFR and Σ_gas. Furthermore, to minimise the redundancy in the scatter plots we explore the Σ_SFR − Σ_gas plane using a grid with Nyquist sampling adapted for each spatial resolution.

Table 1 lists the power-law indexes, as well as the Pearson correlation (ρ_ps) and Spearman rank (ρ_sp) parameters obtained for the KS laws (in logarithmic space) for the different spatial resolutions and tracers used in this work. We carried out the fits of log(Σ_SFR) versus log(Σ_gas) using the orthogonal distance regression (ODR) method. We considered that a correlation is noteworthy and statistically significant when its two-sided p-value < 1% and both ρ_ps and ρ_sp are ≥0.4. Validated correlations are highlighted in boldface in Table 1. Furthermore, in order to counterbalance the dependence of the estimated p-values on the size of the sample, we derived these using a common number of randomly selected points for the different spatial scales.

As an illustration of the wide variance of SF relations resulting from this analysis, we show in Fig. 5 the different versions of the KS law representing Σ_SFR as a function of Σ_dense (here derived from HCN) for all the spatial scales explored in this work, namely from the ‘initial resolution’ (56 pc) up to 700 pc. This figure shows that the correlation becomes looser with higher resolution, and it is hardly visible in the plot with a resolution of 56 pc. The KS relation for the dense gas at the initial resolution is highly scattered: Σ_SFR and Σ_dense span 1.5 dex each and show scarce evidence of correlation. We fit a super-linear KS relation with a power-law slope N = 3.74 ± 0.22. The correlation parameters at the initial resolution have low values (ρ_sp = 0.38 and ρ_ps = 0.42) and their corresponding two-sided p-values > 1%. At scales of 100 pc we obtain a power-law slope N = 2.17 ± 0.19 with correlation parameters of ρ_sp = 0.35 and ρ_ps = 0.43. In contrast, at 400 pc and 700 pc, the scatter in the KS relation is significantly reduced: Σ_SFR and Σ_dense span 1 dex each and the derived best-fit KS relation yields a power-law index N = 1.01 ± 0.22 and N = 1.08 ± 0.12, respectively, with correlation parameters ρ_sp = 0.43, ρ_ps = 0.45 and ρ_sp = 0.85, ρ_ps = 0.82, with associated p-values < 1%. We represent the KS plots for all the tracers in Figs. C.1–C.4.

Fig. 5. SFR surface density (ΣSFR) as a function of the dense molecular gas surface density (Σdense) derived from HCN(1–0) for the different resolutions explored in this work: 56 pc, 100 pc, 200 pc, 300 pc, 400 pc, 500 pc, 600 pc and 700 pc. The red lines show the ODR fits to the data. The dashed lines identify the correlations that are not found to be statistically significant. We only show on both axes the data points above 3σ that are considered in the fits. We indicate the Spearman’s rank correlation coefficients (ρsp) and the power-law indexes (N) of the best fits. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.

Figure 6 shows how ρ_ps and ρ_sp change as a function of the spatial resolution used for the different molecular gas and dust tracers used in this comparison. The grey-shaded region in Fig. 6 identifies the range of spatial scales where the correlation is judged not to be statistically significant. For spatial scales ≥300 − 400 pc the correlation improves monotonically as a function of the aperture size for all gas tracers excluding CO(1–0). The ρ_ps and ρ_sp parameters for CO(1–0) show values ≤0.3 for the entire range of spatial resolutions explored. We nevertheless find that for any spatial resolution the correlation parameters derived from the high-density tracers (CO(3–2), HCN(1–0) and HCO⁺(1–0)) are about a factor of two to three larger than that derived from CO(1–0). Dust continuum emission shows a behaviour similar to that of the rest of high-density tracers with the important particularity that the correlation is significant already at 100 pc scales. This result confirms that continuum emission is mostly sensitive to the column densities of the dust that is being directly heated by recent SF activity in the SB ring.

Fig. 6. Pearson and Spearman correlation coefficients (ρps: left panel; ρsp: right panel) of the KS laws derived for the different gas tracers used in this work as a function of the spatial resolution. Error bars represent the 67% confidence interval around the mean values. For spatial scales larger than the range identified by the grey-coloured region (≤350 pc) all the gas tracers (leaving aside CO(1–0)) show statistically significant correlations in their KS laws, defined by p-values < 1% and ρps, ρsp > 0.4.

The reported breakdown of the KS relations observed in the SB ring below a critical spatial scale of ∼300 − 400 pc is in qualitative agreement with the findings of Onodera et al. (2010), Schruba et al. (2010) in M 33 and Kreckel et al. (2018) in NGC 628. The existence of a ‘critical scale’ for KS laws can be explained by the diverse evolutionary states of the GMC population, which can be singled out in high spatial resolution observations. The exact value of this critical scale may change from galaxy to galaxy. It may also depend on the criterion adopted to choose the data points used to generate the scatter plots: either a ‘blind’ pixel-wise Nyquist sampling or a ‘biased’ selection of apertures centred either around SF or gas emission peaks (e.g. see discussion in Williams et al. 2018). Different GMC states can reflect an ordered time sequence of SF determined by the large-scale dynamics in galaxy disks or a more stochastic pattern due to the local dispersal of molecular gas by stellar feedback.

Figure 7 shows how the power-law index N obtained from the fit to the KS relation changes as a function of the spatial resolution for the different molecular gas and dust tracers. As in Fig. 6, the grey-shaded region in Fig. 7 identifies the range of spatial scales where the correlation is not statistically significant. Overall, we find a strong scale dependence of N. The power-law index for HCN and HCO⁺ shows a fairly systematic decrease with the spatial resolution from N = 3.5 ± 0.3 (at the initial resolution) to N = 1.0 ± 0.1 (at 700 pc). The value of N stays around ∼1.0 ± 0.1 for the whole range of spatial scales where the correlation is significant (≥300 − 400 pc). For CO(3–2) N shows values marginally below unity in the 300 pc–500 pc range. Similarly, the value of N for the dust continuum indicates a sub-linear relation (N < 1) within the range 100 pc–600 pc. The power law becomes nevertheless linear at 700 pc as for most of the high-density tracers. The reported slightly different behaviour of CO(3–2) and dust continuum relative to HCN or HCO⁺ can be attributed to the fact that, although all these tracers are sensitive to the presence of dense molecular gas, CO(3–2) and dust continuum are also to a large extent mostly sensitive to the presence of comparatively hotter molecular gas, characterised by high kinetic temperatures⁶.

Fig. 7. Power-law indexes (N) of the KS relations obtained from the ODR fit to the data for the different gas tracers used in this work as a function of the spatial resolution adopted. Symbols and colour codes are the same as in Fig. 6.

The linear behaviour of the power law observed for all the dense gas tracers in NGC 1068 is in agreement with the results obtained in other galaxies (e.g. Gao & Solomon 2004a,b; Graciá-Carpio et al. 2008; Wu et al. 2010; García-Burillo et al. 2012; Usero et al. 2015; Liu et al. 2015; Chen et al. 2017; Williams et al. 2018; Querejeta et al. 2019). The weak correlation found between Paα and CO(1–0) at all scales indicates that the distribution of the general molecular gas traced by CO(1–0) is not strongly correlated with the current location of recent SF in the SB ring. While we do not expect to have missed a high fraction of the CO(1–0) flux in the PdBI map on scales ≤200 − 300 pc (see Sect. 2.1.2 and Appendix A), the likely increasing fraction of flux filtered on scales close to the ∼1 kpc limit could be an explanation for the poor correlation shown by CO(1–0) in Fig. 6.

5.2. The Kennicutt–Schmidt laws of NGC 1068 in context

We compare in this section the KS relations derived from HCN(1–0) in the SB ring of NGC 1068 with those obtained by previous works in different populations of galaxies. Table 2 lists the references used in this comparison. These works comprise galaxy-scale studies of SFGs and (ultra-)luminous infrared galaxies ((U)LIRGs) (e.g. see the compilation by García-Burillo et al. 2012, and references therein), as well as spatially resolved studies of individual galaxies (Usero et al. 2015; Murphy et al. 2015; Bigiel et al. 2016; Chen et al. 2017; Viaene et al. 2018; Querejeta et al. 2019). These references used the equations of Sect. 3 to convert HCN(1–0) intensities into (deprojected) surface densities of dense gas and adopted the same HCN conversion factor used in this paper. However, the SF tracers chosen in these works are different, as detailed in Table 2.

Figure 8 compares the different versions of the KS law derived from HCN(1–0) in NGC 1068 at three spatial resolutions (56 pc, 100 pc, and 700 pc) with those obtained in the references listed in Table 2. NGC 1068 data lie within the linear power-law branch occupied by the rest of the galaxies shown in Fig. 8. Galactic dense cores, at sub-parsec scale, align along the same relationship (Wu et al. 2005, 2010; Rosolowsky et al. 2011; Stephens et al. 2016; Shimajiri et al. 2017). Taken at face value, this suggests that the SFE (SFR per unit dense molecular mass) is nearly constant on average. As expected, the internal scatter in the distribution NGC 1068 data points is lower as we move to larger apertures. Furthermore, NGC 1068 data points shift towards lower surface densities within the power law as we move to larger apertures.

Fig. 8. Comparison of the KS laws derived from HCN(1–0) in NGC 1068 with those obtained in different populations of galaxies for a range of spatial scales: global measurements (around a few kiloparsec) of the disks of SFGs and (U)LIRGs (García-Burillo et al. 2012), ≃0.3−1 kiloparsec-scale regions of the disks of nearby SFGs (Murphy et al. 2015; Usero et al. 2015; Bigiel et al. 2016), and 0.1−0.2 kpc size regions in the disk of M 51 (Chen et al. 2017; Viaene et al. 2018; Querejeta et al. 2019). For the sake of a better comparison with the reference sample, the scatter plots for NGC 1068 (magenta circles and isodensity contours) are displayed for three spatial resolutions: 56 pc (left panel), 100 pc (middle panel), and 700 pc (right panel). The scatter plot for M 51 at 100 pc resolution is displayed by the green circles and isodensity contours. The open squares stand for the median values of the distributions in NGC 1068 and M 51.

Compared at a common scale of 100 pc, the SB ring of NGC 1068 appears as a more extreme environment relative to the SF regions of M 51 studied by Querejeta et al. (2019). As illustrated in Fig. 8, the median values of the distributions of Σ_SFR and Σ_dense are about a factor of three to five higher in NGC 1068: Σ_SFR[NGC 1068] ≃ 1.5 M_⊙ Myr⁻¹ pc⁻² ≃ 3 × Σ_SFR[M 51] and Σ_dense[NGC 1068] ≃ 174 M_⊙ pc⁻² ≃ 5 × Σ_dense[M 51]. Furthermore, when examined at scales of 700 pc, NGC 1068 occupies in the KS plot a position intermediate between that of normal galaxies and (U)LIRGs, a result that hints at the relatively extreme conditions in the SB ring.

5.3. Environmental dependence of the star-formation efficiency of the dense gas

The observed linear relation between Σ_SFR and Σ_dense shown in Fig. 8 has been considered as evidence of the validity of density-threshold models of SF. For these models the SFE of dense molecular gas, defined as SFE_dense = Σ_SFR/Σ_dense or its inverse, which represents the depletion time of the dense gas ($T_{\mathrm{dep}}^{\mathrm{dense}}={\mathrm{SFE}}_{\mathrm{dense}}^{-1}$), are about constant in different populations of galaxies and also for different dynamical environments within galaxies (Gao & Solomon 2004b; Wu et al. 2005; Lada et al. 2010, 2012; Evans et al. 2014). In this section we use the NGC 1068 data to explore the existence of an environmental dependence of SFE_dense inside the SB ring.

Figure 9 (left panel) overlays the HST/NICMOS Paα image on the HCN(1–0) map of NGC 1068 obtained at the initial resolution of 56 pc. We identify the positions used to extract the fluxes of Paα and HCN(1–0) over the SB ring from the Nyquist-sampled grid⁷. SFE_dense values span almost 1.5 dex and show a highly scattered distribution as a function of $L_{\text{HCN}}^{\prime }$ around an ‘apparently’ constant mean value of about 0.01 Myr⁻¹, equivalent to a $T_{\mathrm{dep}}^{\mathrm{dense}}\sim 100$ Myr (see right panel of Fig. 9). This result is similar to the findings of Gallagher et al. (2018) and Querejeta et al. (2019) obtained from their high spatial resolution images of a sample of galaxies. Notwithstanding that we might attribute part of the scatter in the SFE_dense − $L_{\text{HCN}}^{\prime }$ plot to possible small-scale variations in the α_HCN conversion factor, any plausible range for these potential variations would nevertheless fall short of accounting for the bulk of the 1.5 dex span shown in Fig. 9.

Fig. 9. Distribution of the HCN and Paα peaks and the apertures showing the highest SFEdense in NGC 1068. Left panel: overlay of the HST/NICMOS Paα image (contours) on the HCN(1–0) map (grey scale). The contours and grey scale are the same as in Fig. 4. The set of 56 pc size circular apertures used to extract the fluxes of Paα and HCN(1–0) are colour-coded differently, depending on whether they are centred on HCN peaks (in blue), on Paα peaks (in red), or on the 10% of the apertures showing the highest values of SFEdense (in green). The black dots identify the centres of all the Nyquist-sampled 56 pc size apertures where there are reliable estimates of SFEdense across the SB ring. Right panel: distribution of SFEdense values (in Myr−1-units) as a function of LHCN (in K km s−1 pc2-units) in the SB ring of NGC 1068. Symbols are the same as in left panel. Vertical and horizontal error bars at the lower-right corner of the panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

Although Figs. 8 and 9 indicate that there is an overall relationship between the HCN luminosity and recent SF in the SB ring, we explore below the existence of systematic trends in SFE_dense. With this aim, we selected within our initial grid a number of non-overlapping 56 pc size apertures centred on local maxima either in the HCN or the Paα maps, following the same procedure used by Querejeta et al. (2019) in their analysis of M 51 data. The local maxima were identified as the pixels with peak intensities within circular regions (‘clumps’), obtained from successive cuts in the maps. We therefore started from a high threshold (a large multiple of the noise level) and iteratively explored lower values from isolated circular regions. To collect a similar number of points selected from HCN and Paα peaks we continued to identify new clumps down to different threshold values of 12σ and 44σ, respectively. We use different colours to identify in the two panels of Fig. 9 the apertures centred on HCN peaks (blue colour), on Paα peaks (red colour), and also on the 10% of the apertures showing the highest values of SFE_dense (green colour).

Figure 9 shows that the apertures showing the top 10% SFE_dense values are not uniformly distributed throughout the SB ring. High SFE_dense values are preferentially located at the south-west and north-east sections of the ring, namely where the ring is connected to the ends of the stellar bar. We also identify high SFE_dense values further out at the north-east extreme of the SB ring. As expected, the clumps with the strongest Paα emission are also located in the regions with the highest SFE_dense. Furthermore, although the overall distribution of the highest (dense) gas column density apertures also corresponds to the ends of the stellar bar, on small scales there is no one-to-one correspondence between Paα and HCN peaks. This reflects the significant variance in the evolutionary states of the dense gas clumps in these regions. In particular, the gas ridge traced by the HCN maxima tends to appear ‘upstream’ relative to the Paα maxima at the south-west section of the SB ring⁸.

In order to minimise the bias introduced due to the presence of HCN in the two axes of the right panel of Fig. 9, we replaced HCN by HCO⁺ luminosities along the x axis by selecting from our grid of points those positions satisfying $L_{{\text{HCO}}^{\text{+}}}^{\prime }>3\sigma$. Figure 10 shows the results obtained following the same steps leading to Fig. 9. As expected, Fig. 10 eliminates the slight anti-correlation trend identified in the SFE_dense − $L_{\text{HCN}}^{\prime }$ plot of Fig. 9, but the main results described above remain virtually unchanged.

Fig. 10. Same as Fig. 9 but replacing HCN(1–0) with HCO+(1–0) along the x axis of the right panel.

6. Star-formation efficiency and dense gas fraction

The overall SFE of molecular gas (SFE_mol) can be expressed as the product of SFE_dense and the dense gas fraction (F_dense), namely: SFE_mol ≡ SFE_dense × F_dense. There is mounting evidence supporting the existence of significant variations in SFE_dense and F_dense as a function of the galactic environment, based on observations of molecular clouds in the centre of our Galaxy (Longmore et al. 2013; Kruijssen et al. 2014) and in nearby galaxies for a range of spatial scales (Usero et al. 2015; Bigiel et al. 2015, 2016; Gallagher et al. 2018; Querejeta et al. 2019; Jiménez-Donaire et al. 2019; Bešlić et al. 2021). In the following sections we study the trends in SFE_dense and F_dense as a function of a number of physical variables with the aim of resolving the degeneracy of SF laws in the SB ring of NGC 1068. Table 3 lists the Spearman rank parameters for the different combinations of variables and spatial scales explored below.

6.1. Trends as a function of the boundedness of the gas

In this section we use an alternative prescription for the SF relations of the SB ring, which includes explicitly the dependence of SFE_dense on a combination of Σ_dense and the gas velocity dispersion (σ), in an attempt to resolve the degeneracy associated with the scatter in the SFE_dense − L′(HCN) plot of Fig. 9. This approach was first used by Leroy et al. (2017) in their analysis of SF relations in M 51 and also adopted by Kreckel et al. (2018) in a similar study carried out in NGC 628. Leroy et al. (2017) and Kreckel et al. (2018) used CO(1–0) and CO(2–1) to trace the bulk of the molecular gas on 40−50 pc-scales in M 51 and NGC 628, respectively. In the following analysis we use HCN(1–0) to specifically trace the dense molecular gas in the SB ring of NGC 1068 on spatial scales (∼56 pc) comparable to those explored by Leroy et al. (2017) and Kreckel et al. (2018)⁹.

The degree of self-gravity or boundedness of molecular gas clouds determines to a large extent their ability at forming stars. The virial parameter, defined as α_vir ≈ 2KE/UE, captures the balance of gravitational potential (UE) and kinetic energy (KE) and is commonly used in turbulent models of SF as a predictor of the efficiency of SF (Krumholz & McKee 2005; Padoan et al. 2012, 2017).

We define similarly to Leroy et al. (2017) the boundedness parameter, here particularised for the dense molecular gas phase traced by HCN, as $b_{56\mathrm{pc}}\equiv {\mathrm{\Sigma }}_{\mathrm{dense}}/{\sigma }^2\propto \mathrm{UE}/\mathrm{KE}\propto {\alpha }_{\mathrm{vir}}^{-1}$, where σ is the velocity dispersion and Σ_dense is the column density of dense gas measured both at the initial resolution of 56 pc. We derived intensity-weighted averages of b_56 pc using two different apertures ΔA = 100 pc and 400 pc over the Nyquist sampled grid, defined as follows:

$$\begin{array}{c}\hfill {⟨b⟩}_{\mathrm{\Delta }\mathrm{A}}(x_0,y_0)={\displaystyle \frac{\mathrm{\Sigma }w(x,y)I_{56\mathrm{pc}}(x,y)b_{56\mathrm{pc}}(x,y)}{\mathrm{\Sigma }w(x,y)I_{56\mathrm{pc}}(x,y)}}·\end{array}$$(10)

The Gaussian weight w(x, y) is defined as

$$\begin{array}{c}\hfill w(x,y)=exp\left(\frac{-{(\theta (x,y,x_0,y_0))}^2}{2{\sigma }_{\mathrm{\Delta }\mathrm{A}}^2}\right),\end{array}$$(11)

where θ(x, y, x₀, y₀) is the angular distance from the measurement point (x₀, y₀), σ_ΔA is the σ-width of the Gaussian averaging beam, and ΔA corresponds to the adopted averaging scale (100 pc and 400 pc in our case). In Eq. (10), b_56 pc and I_56 pc are the boundedness parameter and the HCN integrated intensity, respectively, measured at a generic position (x, y) of the grid¹⁰.

Figure 11 shows how the depletion time of the dense molecular gas estimated from HCN, $T_{\mathrm{dep}}^{\mathrm{dense}}\equiv {\mathrm{SFE}}_{\mathrm{dense}}^{-1}$, changes as a function of the self-gravity of the gas, measured by the b parameter, at three spatial scales: 56 pc (the initial resolution), 100 pc, and 400 pc. The values of ⟨b⟩_ΔA for the 100 pc and 400 pc apertures were derived from Eq. (10). Furthermore, similarly to the rest of the physical parameters analysed in Sects. 6.2–6.4, the average value of $T_{\mathrm{dep}}^{\mathrm{dense}}$ inside the two values of ΔA were derived using the Gaussian weighting function of Eq. (11).

Fig. 11. Depletion time of the dense molecular gas estimated from HCN, $T_{\mathrm{dep}}^{\mathrm{dense}}\equiv {\mathrm{\Sigma }}_{\mathrm{dense}}$/ΣSFR, as a function of the self-gravity of the gas, measured by the parameter b ≡ Σdense/σ2, at different spatial scales (magenta circles and isodensity contours): 56 pc, namely the initial resolution of ALMA observations (left panel), 100 pc (middle panel), and 400 pc (right panel). A significant trend in $T_{\mathrm{dep}}^{\mathrm{dense}}$ as a function of b, identified in the middle and right panels, indicates a higher rate of SF per unit mass of dense gas for regions with stronger self-gravity. Values of ⟨b⟩ for the 100 pc and 400 pc apertures were derived from an intensity-weighted average of b evaluated at the Nyquist-sampled grid of points for each aperture. Values of $T_{\mathrm{dep}}^{\mathrm{dense}}$ were derived using a Gaussian-weighted average. Black circles show median $T_{\mathrm{dep}}^{\mathrm{dense}}$, and error bars indicate the rms scatter in the bins. The green lines in the right panel identify the two regimes in the $T_{\mathrm{dep}}^{\mathrm{dense}}-⟨b⟩$ plot found by the MARS algorithm to fit the data using 400 pc as averaging scale. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

Figure 11 shows a significant monotonic decrease in $T_{\mathrm{dep}}^{\mathrm{dense}}$ as a function of ⟨b⟩_ΔA when we represent these parameters averaged over the two apertures. Specifically, we obtain (anti) correlation Spearman rank parameters ρ_sp = −0.43 and −0.66 for ΔA = 100 pc and 400 pc, respectively, and associated two-sided p-values < 1%. Overall, this is indicative of a higher rate of SF per unit mass of dense gas (lower $T_{\mathrm{dep}}^{\mathrm{dense}}$) for regions characterised by a stronger self-gravity (higher b or lower α_vir values).

The anti-correlation reported above is more pronounced in the intensity-weighted version of the plot derived at 400 pc scales. We identify a turnover in the T_dep − b scatter plot for ΔA = 400 pc located around log₁₀(⟨b⟩) ≃ 0.2 M_⊙ pc⁻² (km s⁻¹)⁻², as shown in the right panel of Fig. 11. The observed change of tendency around this point defines two regimes in the T_dep − b parameter space. To quantify the two-regime trend and the turnover, we used the Multivariate Adaptive Regression Splines (MARS) fit routine from the Rstudio package¹¹. We performed the MARS fit on 100 realisations of the T_dep versus b relation at ΔA = 400 pc taking into account the data uncertainties. From this Monte Carlo simulation, we obtained a turnover at log₁₀(⟨b⟩) = 0.23 ± 0.02 M_⊙ pc⁻² (km s⁻¹)⁻² which is indicated in the right panel of Fig. 11. Furthermore, for ⟨b⟩ values below the turnover the trend is approximately flat (slope = 0.02 ± 0.11). In contrast, the MARS routine fits a slope = −1.59 ± 0.15 beyond the turnover. This slope is larger than the one obtained by Leroy et al. (2017) (≃ − 0.9) in their analysis of the M 51 CO(1–0) data, derived using the same averaging spatial scales (400 pc). Kreckel et al. (2018) found nevertheless no significant correlation in their analysis of the CO(2–1) data of NGC 628, which used an averaging scale of 500 pc. Compared to the shallower or nonexistent trends identified in M 51 and NGC 628, the steeper decline of $T_{\mathrm{dep}}^{\mathrm{dense}}$ with ⟨b⟩ in the SB ring of NGC 1068 reflects the tighter link between the boundedness of dense molecular gas and SFE.

Figure 12 shows the spatial distribution of $T_{\mathrm{dep}}^{\mathrm{dense}}$ and ⟨b⟩ derived for ΔA = 400 pc in the SB ring of NGC 1068. The $T_{\mathrm{dep}}^{\mathrm{dense}}$ map of Fig. 12 confirms the picture drawn from the analysis of Sect. 5.3: the regions showing comparatively higher (lower) SFE_dense ($T_{\mathrm{dep}}^{\mathrm{dense}}$) values are preferentially located closer to the region where the SB ring is connected to the ends of the stellar bar around PA ∼ 15° −75° ( ± 180°). We also see a similar spatial segregation in the ⟨b⟩ map, which shows higher ⟨b⟩ values closer to the bar-ring interface region. Besides the azimuthal dependence of $T_{\mathrm{dep}}^{\mathrm{dense}}$ and ⟨b⟩ within the SB ring, we also identify a radial dependence for both parameters especially in the southern section of the bar-ring interface: in particular, the highest (lowest) $T_{\mathrm{dep}}^{\mathrm{dense}}$ (⟨b⟩) values tend to appear upstream (smaller radii) along the gas circulation lines if we assume that the sense of rotation of the gas in the disk is counterclockwise. The two branches in the $T_{\mathrm{dep}}^{\mathrm{dense}}-⟨b⟩$ plot of Fig. 11 correspond to a large extent to the two regions of the SB ring identified in Fig. 12. Similar results are found in Appendix E when we consider the HCO⁺ as dense gas tracer. The trends and statistical parameters derived are practically identical to the ones obtained using HCN.

Fig. 12. Maps of the depletion time and the b parameter in the SB ring of NGC 1068. Left panel: map of the depletion time of the dense molecular gas ($T_{\mathrm{dep}}^{\mathrm{dense}}$) obtained after using 400 pc as an averaging scale. Contour levels for log10($T_{\mathrm{dep}}^{\mathrm{dense}}$) go from 1.50 to 2.22 in steps of 0.045 (in Myr). Right panel: map of the intensity-weighted average of the b parameter (⟨b⟩) in the SB ring of NGC 1068 obtained after using 400 pc as an averaging scale. Contour levels for log10(⟨b⟩) go from 0.075 to 0.45 in steps of 0.025 (in M⊙ pc−2 (km s−1)−2). The Δα and Δδ offsets in arc seconds are relative to the phase tracking centre of ALMA. The ellipse identifies the extent of the nuclear stellar bar in the disk. We highlight the approximate extent of the bar-ring interface region in both panels.

We can speculate if the reported trends of $T_{\mathrm{dep}}^{\mathrm{dense}}$ as a function of ⟨b⟩ for ΔA = 400 pc, visualised in Figs. 11 and 12, could be entirely attributed to potential variations in the α_HCN conversion factor on these spatial scales. In this context it is worth noting that the range explored by $T_{\mathrm{dep}}^{\mathrm{dense}}$ as a function of ⟨b⟩, which amounts to 0.6 dex (a factor of four), is a significant factor of three larger than the differences found between the global kiloparsec-scale and the ‘small’ parsec-scale HCN conversion factors reported by Wu et al. (2005), as mentioned in Sect. 3.1. Although we have no way of confirming or refuting the existence of significantly larger changes of α_HCN based on our data, we note that adopting a lower value of the conversion factor for HCN would result in lower molecular gas surface densities, particularly in the regions of the SB ring of NGC 1068 that happen to show comparatively higher (lower) SFE_dense ($T_{\mathrm{dep}}^{\mathrm{dense}}$) values. These regions very likely comprise a collection of hot core-like clouds akin to the Galactic SF cores studied by Wu et al. (2005) for which ${\alpha }_{\mathrm{HCN}}^{\mathrm{cores}}\simeq 0.7\times {\alpha }_{\mathrm{HCN}}^{\mathrm{global}}$. As a direct consequence, the trends shown in Figs. 11 and 12 would be further enhanced rather than being suppressed.

6.2. Trends as a function of the dense gas fraction

Observations of molecular gas in our Galaxy and in a number of nearby galaxies have found clear evidence of trends in SFE_dense as a function of F_dense, which are indicative of anti-correlation (e.g. Longmore et al. 2013; Chen et al. 2015; Murphy et al. 2015; Usero et al. 2015; Bigiel et al. 2016; Gallagher et al. 2018; Querejeta et al. 2019; Jiménez-Donaire et al. 2019). In particular, Usero et al. (2015) found that SFE_dense is about six to eight times lower near the galaxy centres than in the outer regions of the galaxy disks analysed in their IRAM 30-m survey. This radial trend is reversed for F_dense in their sources. Furthermore, Longmore et al. (2013) found anomalously low SFR values in a large fraction of high-density molecular clouds in the centre of the Milky Way, suggestive of an anti-correlation between SFE_dense and F_dense. A similar trend has been found recently by Querejeta et al. (2019) in their study of M 51, and by Jiménez-Donaire et al. (2019) from an analysis of the data obtained by the EMPIRE survey in nine spiral galaxies. We examine in this section the existence of a trend in SFE_dense as a function of F_dense in the SB ring of NGC 1068.

We define two proxies for the ‘dense gas fraction’ in the SB ring. First, the ratio of the dense molecular gas surface density derived from HCN (1–0) to the bulk molecular gas surface density derived from CO(1–0), that is, F_dense = Σ_dense/Σ_mol ∝ I_{HCN(1 − 0)}/I_{CO(1 − 0)}. We consider that CO line traces the bulk of the molecular gas in the galaxy (n_H2 ≳ 10^2 − 3 cm⁻³), while HCN traces material with n_H2 ≳ 10^4 − 5 cm⁻³. We therefore use the HCN(1–0)/CO(1–0) line ratio (hereafter R_HCN/CO) as a proxy for the dense gas fraction in the SB ring. Secondly, we also use the CO(3–2)/CO(1–0) ratio (hereafter R_{3 − 2/1 − 0}) in the SB ring derived by García-Burillo et al. (2014) as an alternative proxy for F_dense.

Figures 13 and 14 show R_HCN/CO and R_{3 − 2/1 − 0} in T_mb units. To derive the brightness temperature ratios, we degraded the HCN(1–0) and CO(3–2) maps to the spatial resolution of the CO(1–0) observations of Schinnerer et al. (2000). The R_HCN/CO ratio changes significantly across the disk of NGC 1068. The highest values of R_HCN/CO, > 0.4 − 1.2, correspond to the CND. Although in the SB ring R_HCN/CO shows a wide range of values, this ratio is higher in the bar-ring interface region (R_HCN/CO ∼ 0.2 − 0.4) than elsewhere in the ring (R_HCN/CO ∼ 0.05 − 0.2). The R_{3 − 2/1 − 0} ratio, shown in Fig. 14, changes also significantly across the disk. The R_{3 − 2/1 − 0} ratio in the CND (∼2 − 6) is higher than in the SB ring, in agreement with previous estimates by Krips et al. (2011) and Tsai et al. (2012). While the average ratio is ∼1.2 ± 0.02 in the SB ring, R_{3 − 2/1 − 0} is higher in the bar-ring interface region (R_{3 − 2/1 − 0} ∼ 0.6 − 3) and comparatively lower elsewhere in the ring (R_{3 − 2/1 − 0} ∼ 0.1 − 0.8).

Fig. 13. HCN(1–0)/CO(1–0) brightness temperature ratio map (RHCN/CO; colour scale and contours) derived at the common (lower) spatial resolution of the CO(1–0) observations of Schinnerer et al. (2000) (≃100 pc). The brown circle of 8″ in radius (≃560 pc) locates the inner region excluded from our analysis of SF relations.

Fig. 14. Same as Fig. 13 but showing the CO(3–2)/CO(1–0) brightness temperature ratio map (R3 − 2/1 − 0; colour scale and contours) derived at the common (lower) spatial resolution of the CO(1–0) observations.

The R_HCN/CO ratio is hardly sensitive to kinetic temperature (T_k), as the energy levels giving rise to both rotational transitions are similar (E_k[J = 1, CO] ∼ 5.5 K, E_k[J = 1, HCN] ∼ 4.3 K). However, as the excitation of both CO lines are sensitive to both n_H2 and T_k, R_{3 − 2/1 − 0} is comparatively a more indirect and less straightforward tracer of F_dense relative to R_HCN/CO. Leaving aside the uncertainties on the value of α_HCN in SF regions, which may reflect a peculiar hot-core like chemistry, R_HCN/CO is therefore the most reliable proxy for the dense gas fraction and as such is widely used in extragalactic studies, and we therefore adopt it in the following analysis. In either case, we note that both line ratios suggest a higher excitation of HCN(1–0) and CO(3–2) lines relative to CO(1–0) in the bar-ring interface region.

Figure 15 represents SFE_dense as a function of F_dense estimated from R_HCN/CO in the SB ring. We adopted the same approach followed in Sect. 6.1 to obtain estimates of both variables over two averaging scales: ΔA = 100 pc and 400 pc. Figure 15 shows that there is no significant trend in SFE_dense as a function of F_dense for ΔA = 100 pc: we estimate a Spearman rank parameter ρ_sp = +0.03 with a two-sided p-value > 1%. For ΔA = 400 pc there is nevertheless a more significant positive correlation (ρ_sp = +0.34, with a p-value < 1%). This trend is a direct consequence of the spatial distribution of SFE${}_{\mathrm{dense}}^{-1}$ and F_dense shown, respectively, in Figs. 12 and 13. The monotonic increase shown in the right panel of Fig. 15 suggests that the comparatively denser molecular gas of the bar-ring interface region inside the SB ring forms stars at higher rate per unit dense gas mass.

Fig. 15. Star-formation efficiency of the dense gas (SFEdense) as a function of dense gas fraction (Fdense: estimated from the HCN(1–0)-to-CO(1–0) line ratio of Fig. 13) in NGC 1068 using two averaging scales: 100 pc (left panel) and 400 pc (right panel). Magenta circles correspond to data measurements evaluated at the Nyquist-sampled grid of points after convolution with the appropriate Gaussian for each aperture. Black circles show median SFEdense values, and error bars indicate the rms scatter in the bins. Vertical and horizontal error bars at the lower-right corner of both panels account for the typical uncertainties, which amount to ±0.13 dex.

This result seems to be in contradiction with the anti-correlation trends found between SFE_dense and F_dense in previous works in other galaxies. Specifically, Usero et al. (2015) found an index N ∼ −1.6 for the SFE${}_{\mathrm{dense}}\propto F_{\mathrm{dense}}^N$ power-law fitting their data. A similar index (N ∼ −1.5) can be estimated from the data compiled by Querejeta et al. (2019).

6.3. Trends as a function of the stellar mass surface density

Several works studying the existence of trends in SFE_dense as a function of Σ_star on kiloparsec scales in a number of nearby galaxies have found that SFE_dense is seen to decrease in the central parts of galaxy disks, namely in regions characterised by high Σ_star. These observations also showed that F_dense tends to increase systematically with Σ_star (Chen et al. 2015; Usero et al. 2015; Bigiel et al. 2016; Gallagher et al. 2018; Jiménez-Donaire et al. 2019). Querejeta et al. (2019) studied different regions in the disk of M 51 and found that similar correlations are recovered at 100 pc-scales. We investigate below the trends in SFE_dense and F_dense as a function of Σ_star in the SB ring of NGC 1068.

Near-infrared emission is commonly used to trace the stellar mass in nearby galaxies, since light at these wavelengths mainly comes from old stars and is less affected by extinction (Quillen et al. 1994). With this aim, we used the HST NIR continuum narrow-band image of NGC 1068 at 1.9 μm, obtained with the F190N filter on the NICMOS 3 camera, to trace the distribution of the stellar mass in the disk of the galaxy. We assumed a constant mass-to-light ratio M/L_1.9 μm = 0.2 M_⊙/L_⊙ (e.g. see Querejeta et al. 2015) to obtain the local stellar mass surface density. Figure 16 overlays the HCN(1–0) integrated intensity contours on the HST/NICMOS F190N image. The stellar bar feature is easily identified in the NIR image. Judging from the overall distribution of Σ_star in the disk, it appears that the bar-ring interface is characterised by higher Σ_star values compared to the regions located elsewhere in the SB ring.

Fig. 16. Overlay of the ALMA HCN(1–0) map (contours) on the HST/NICMOS F190N continuum image at 1.9 μm HST (colour scale as shown in Jy). The contour spacing is 3σ, 7σ, 12σ, 24σ, 50σ, and 120σ, where 1σ = 0.028 Jy beam−1 km s−1. The white circle of r = 8″ locates the inner region of the galaxy disk, which was excluded of our analysis. The HST/NICMOS F190N image has been degraded to the spatial resolution of HCN(1–0) (1${}^{″}\times 0\stackrel{″}{.}6\approx 50$ pc).

Figure 17 represents SFE_dense as a function of Σ_star across the SB ring of NGC 1068 for the three working apertures adopted in Sect. 6.1. SFE_dense does not show any significant correlation with Σ_star at the initial resolution (ρ_sp = +0.10, with a p-value > 1%), or at the 100 pc averaging scale (ρ_sp = +0.15, with a p-value > 1%). On the contrary, SFE_dense shows a statistically significant positive correlation with Σ_star in the SB ring for ΔA = 400 pc (ρ_sp = +0.59, with a p-value < 1%). We compare the location of NGC 1068 in the SFE_dense − Σ_star parameter space with the position occupied by the galaxies studied by the references listed in Table 2 in the middle panel of Fig. 17. NGC 1068 clearly deviates from the overall anti-correlation trend followed by other galaxies. The positive correlation observed in the NGC 1068 SB ring is in stark contrast with the anti-correlation trends identified in the galaxies studied by Usero et al. (2015) (ρ_sp = −0.50), Querejeta et al. (2019) (ρ_sp = −0.57), and Gallagher et al. (2018) (ρ_sp = −0.66). In this context it is noteworthy that Querejeta et al. (2019) noticed that for a given value of Σ_star the data points of M 51 span a significant 1 dex range of SFE_dense, namely similar to the spread of values seen in the NGC 1068 SB ring. This is an indication of the wide range of dynamical environments probed in both galaxies.

Fig. 17. Star-formation efficiency of the dense gas (SFEdense) as a function of the stellar mass surface density (Σstar) in the SB ring of NGC 1068 at the initial resolution of 56 pc (left panel) and at two averaging scales: 100 pc (middle panel) and 400 pc (right panel). Magenta circles and isodensity contours correspond to data measurements evaluated at the Nyquist-sampled grid of points after convolution with the appropriate Gaussian for each aperture. Black circles show median SFEdense values, and error bars indicate the rms scatter in the bins. Symbols and isodensity contours in the middle panel identify the data obtained in other galaxies. The dashed orange line represents the fit of Usero et al. (2015). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

Figure 18 represents F_dense as a function of Σ_star across the SB ring of NGC 1068 for the two averaging scales adopted in Sect. 6.1. F_dense shows already a significant correlation with Σ_star at ΔA = 100 pc (ρ_sp = +0.40, with a p-value < 1%). The correlation is further reinforced when we use 400 pc as averaging scale (ρ_sp = +0.67, with a p-value < 1%). This is qualitatively and quantitatively similar to the trends identified in the data of the galaxies shown in the left panel of Fig. 18 at scales of 100 pc, for which Querejeta et al. (2019) derived a Spearman rank parameter of ρ_sp ∼ +0.45, and ρ_sp ∼ +0.7 when they added the data points from Chen et al. (2017) and Gallagher et al. (2018).

Fig. 18. Fraction of dense gas (Fdense) as a function of the stellar mass surface density (Σstar) in the SB ring of NGC 1068 at the two averaging scales: 100 pc (left panel) and 400 pc (right panel). Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

6.4. Trends as a function of velocity dispersion

The role of velocity dispersion in setting the efficiency of molecular gas at forming stars is central to current models of turbulence-driven SF via their dependence on the Mach number (e.g. Krumholz & McKee 2005; Hennebelle & Falgarone 2012; Federrath 2015). From the observational point of view, and under the hypothesis that the velocity dispersion measured at a given scale reflects mainly turbulent motions in molecular gas, identifying trends in SFE as a function of local velocity dispersion may indicate whether SF tends to be either enhanced or suppressed by turbulence. The predictions of models on the expected trends of SFE_dense as a function of the ‘observed’ velocity dispersion (σ) differ depending on the role of turbulence versus large-scale bulk motions (e.g. streaming motions) that are driven by the background gravitational potential. The degree of coupling to the large-scale gravitational potential implies that the critical threshold for SF is not universal and that SFE_dense is therefore not constant. Specifically, the recent model published by Meidt et al. (2020) foresees that SFE_dense should decrease with σ at fixed values of F_dense and Σ_mol.

Leroy et al. (2017) found a decreasing trend between the efficiency of the bulk molecular gas (SFE_mol) and σ measured from CO(1–0) in M 51. Querejeta et al. (2019) also found a significant anti-correlation between SFE_dense and σ in M 51 (ρ_sp = −0.62). A comparison of the observations of M 31 (Viaene et al. 2018) and M 51 (Querejeta et al. 2019) with the predictions of the Meidt et al. (2020) model shows a fair agreement for apertures of 100 pc.

Figure 19 examines the trends in SFE_dense as a function of σ across the SB ring of NGC 1068 for the three apertures adopted in Sect. 6.1. We note that the value of σ derived on scales of 56 pc is not a good proxy for the ‘internal’ turbulence of the dense cores probed by HCN, which in all likelihood have sizes of around a few parsec. Instead, σ likely encapsulates a mix of the ‘macroscopic’ turbulence between the cores and the residual gradient of large-scale bulk motions within the ALMA beam. The correlation shown in Fig. 19 is not statistically significant for the two extreme values of the averaging apertures (56 pc and 400 pc). However, there is a weak yet significant trend at 100 pc scales, for which we derive a ρ_sp = −0.28, with a p-value = 0.3%. Overall, the results obtained in NGC 1068 show a very marginal agreement with the results of Viaene et al. (2018) and Querejeta et al. (2019) and, therefore, with the predictions of the Meidt et al. (2020) model. As illustrated in Fig. 20, a fraction of the data corresponding to the SB ring in NGC 1068 lies in the region predicted by Meidt et al. (2020) for F_dense ∼ 0.25 and a range for Σ_mol ∼ 500 − 5000 M_⊙ pc⁻². Although these values are in rough agreement with the ones estimated for a high fraction of the SB ring positions, the decreasing trend in SFE_dense as a function of σ is seen to be much shallower than the one predicted by Meidt et al. (2020). However, as noted by Meidt et al. (2020), a gravitational potential linked to a strong density wave (bar or spiral) may tend to degrade the strength of the expected decreasing trend, which is estimated in their model using a purely axisymmetric disk.

Fig. 19. Star-formation efficiency of the dense gas (SFEdense) as a function of the velocity dispersion of HCN (σ) in NGC 1068 for three averaging scales: 56 pc (left panel), 100 pc (middle panel), and 400 pc (right panel). Measurements at 100 pc scales are compared to literature data for other galaxies. Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.06 dex, respectively.

Fig. 20. Same as the middle panel of Fig. 19 but adding a comparison with the trends of SFEdense as a function of the velocity dispersion of the gas (σ) predicted by the model of Meidt et al. (2020). This model foresees different locations of galaxies in the SFEdense − σ parameter space depending on the fraction of dense gas (Fdense) and, also, on the total molecular gas column densities (Σmol). Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner account for the typical uncertainties, which amount to ±0.13 dex and ±0.06 dex, respectively.

Figure 21 explores the trends of F_dense as a function of σ in the SB ring for the two averaging scales adopted in Sect. 6.1. In contrast with the trends identified in M 31 and M 51 by Viaene et al. (2018) and Querejeta et al. (2019) on 100 pc scales, F_dense fails to show any significant correlation with σ in the SB ring of NGC 1068 for ΔA = 100 pc (ρ_sp = +0.10, with a p-value > 1%), However, a significant positive trend appears when we use 400 pc as averaging scale (ρ_sp = +0.50, with a p-value < 1%).

Fig. 21. Dense gas fraction (Fdense) as a function of velocity dispersion in the SB ring of NGC 1068 at the two averaging scales: 100 pc (left panel) and 400 pc (right panel). Measurements at 100 pc scales are compared to literature data for other galaxies. Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.06 dex, respectively.

6.5. Relation between the star-formation rate and the dense gas fraction

Viaene et al. (2018) observed a significant correlation (ρ_sp = 0.63) between the SFR and F_dense in a number of regions in the disk of M 31 observed with a spatial resolution of 100 pc. Querejeta et al. (2019) explored the correlation between Σ_SFR and F_dense in M 51 after incorporating the data obtained by Chen et al. (2017) in M 51 on similar spatial scales and derived a similar positive trend: ρ_sp = 0.63, with a p-value < 1%. This correlation improved when Querejeta et al. (2019) included the M 31 and NGC 3627 data of Viaene et al. (2018) and Murphy et al. (2015), respectively (ρ_sp = 0.66).

Figure 22 shows the region occupied by the NGC 1068 data in the Σ_SFR − F_dense plot for an adopted averaging scale ΔA = 100 pc. We do not find a strong correlation between Σ_SFR and F_dense (ρ_sp = 0.22) when we consider the NGC 1068 data alone. In particular, this correlation is weaker than the one found between Σ_SFR and Σ_dense for the same spatial scales (ρ_sp = 0.35; see Table 1). This result is an indication that the dense gas fraction is not a better predictor of SF than the dense gas content in the SB ring. However, if we include the data published for the galaxies displayed in Fig. 22, the overall correlation improves significantly (ρ_sp = 0.72). Querejeta et al. (2019) used all the available data at 100 pc scales and similarly concluded that the dense gas fraction does not seem to be a better predictor of the SFR surface density than the dense gas surface density.

Fig. 22. Star formation rate surface density (ΣSFR) as a function of dense gas fraction (Fdense) in NGC 1068 measured at 100 pc scales compared to the data obtained in other galaxies (symbols as in Fig. 8). Vertical and horizontal error bars at the lower-right corner account for the typical uncertainties, which amount to ±0.09 dex and ±0.13 dex, respectively.

7. A scenario for star formation in the SB ring

The results of our work support the relevance of dynamical environment in setting the efficiency of SF of the dense molecular gas in galaxy disks. Specifically, we find that SFE_dense is comparatively boosted by up to a factor of three to four in the bar-ring interface of NGC 1068 relative to the regions located elsewhere in the ring.

The velocity dispersion of the dense gas as derived from the HCN(1–0) line at the initial resolution of ALMA shows little variations over the SB ring. For an averaging scale of 400 pc, Fig. 23 shows that σ changes from ∼14 to ∼19 km s⁻¹ with an estimated mean value of about 16 km s⁻¹¹². Specifically, σ shows only a moderate ∼30% increase in the southern region of the bar-ring interface (see Fig. 23). This is in agreement with the picture drawn from the analysis of the velocity dispersion estimated at 42 pc resolution from the CO(3–2) line, as illustrated by Fig. 10 of García-Burillo et al. (2014). Furthermore, the right panel of Fig. 21 shows that, while F_dense increases by ∼0.6 dex, the corresponding increase in σ (∼0.12 dex) is a factor of three smaller. Although F_dense and Σ_dense are not identical quantities, the super-linear trend of F_dense as a function σ implies that the trends in the boundedness parameter described in Sect. 6.1 are mostly driven by the observed boost in F_dense in the bar-ring interface. All in all, these results suggest that molecular gas undergoes an efficient compression likely as a result of an enhanced rate of cloud-cloud collisions in the bar-ring interface. However, the shallow trends in σ suggest that cloud-cloud collisions have not increased to any significant level the ‘macroscopic’ turbulence between the cores of dense gas probed on scales of 56 pc in this region.

Fig. 23. Same as Fig. 12 but showing the map of the velocity dispersion (σ) obtained using 400 pc as an averaging scale. Contour levels for σ go from 14.0 to 19.0 in steps of 0.5 km s−1.

The SB ring is formed by a tightly wound two-arm spiral structure. The gas accumulates in a region where two density-wave resonances are thought to be overlapping: the inner Lindblad resonance (ILR) of the ∼17 kpc size outer stellar oval and the corotation of the ∼2.6 kiloparsec-diameter stellar bar (Bland-Hawthorn et al. 1997; Schinnerer et al. 2000; Emsellem et al. 2006). The accumulation of gas at the SB ring suggests that the barrier imposed by the corotation of the (nuclear) stellar bar has been overcome¹³. This could reflect either the influence of the outer stellar oval, which produces gas inflow down to its ILR, or the influence of a decoupled (lower pattern speed) spiral mode, which would also induce gas inflow inside its own corotation. Based on a Fourier decomposition of streaming motions derived from CO(3–2), García-Burillo et al. (2014) found the signature of systematic gas inflow across the SB ring, an indication that the two-arm spiral constitutes an independent wave feature characterised by a lower pattern speed. The inward motions detected are particularly strong at the region connecting the bar with the spiral (see Fig. 14 of García-Burillo et al. 2014). In either case, orbital crowding and the intersection of molecular gas on the orbits of the bar and the spiral are seen to enhance SFE_dense in their interface region. Rico-Villas et al. (2021) recently analysed the continuum emissions at 147 and 350 GHz observed with ALMA, together with the Paα HST/NICMOS image of NGC 1068, and identified 14 super star clusters in the SB ring. In agreement with the scenario described above, most of the super star clusters (11 out of 14) are seen to be located at the bar ends.

NGC 1068 is not by any means an isolated case that illustrates how the complex molecular gas dynamics in a bar-spiral arm interface can trigger SF activity in a galaxy disk. In particular, Beuther et al. (2012) found evidence that the particular dynamical environment of the W43 mini-SB complex in the Milky Way has likely increased cloud interactions and the subsequent SF activity in this region, which is located near the Galactic bar-spiral arm interface. Furthermore, Beuther et al. (2017), and more recently, Bešlić et al. (2021) found an increase in SF activity in the bar-spiral arm interface region of the strongly barred galaxy NGC 3627, as a result of the accumulation of gas where the two orbit families (related to the bar and the spiral) intersect, triggering cloud-cloud collisions. Specifically, the work of Bešlić et al. (2021) reveals that, compared to other regions in the disk of NGC 3627, the efficiency of the dense molecular gas derived from Hα/HCN ratios is significantly enhanced in the bar-spiral arm interface region. Similarly to the scenario described by Beuther et al. (2017) for NGC 3627, we conclude that a configuration where the bar and spiral may rotate at two different pattern speeds can explain the intense SF at the bar-ring interface of NGC 1068.

Several hydrodynamical numerical simulations have been used to investigate the role of large-scale bars in AGN fuelling and in triggering SF in galaxy disks (e.g. Renaud et al. 2013, 2015; Emsellem et al. 2015). In particular, Renaud et al. (2015) showed that the SFE and the formation of massive stellar associations are enhanced at the extremities of the bar to a level comparable to the one observed in galaxy-galaxy interactions as a result of cloud-cloud collisions. This scenario seems to account for the observed pattern of SFE_dense in NGC 1068.

8. Summary and conclusions

We used ALMA to image the emission of dense molecular gas in the r ∼ 1.3 kpc SB ring of the Seyfert 2 galaxy NGC 1068 with a resolution of ∼56 pc using the 1–0 transitions of HCN and HCO⁺. We also used ancillary data of CO(1–0) as well as CO(3–2) and its underlying continuum emission at the resolutions of ∼100 pc and ∼40 pc, respectively. These observations allow us to probe a wide range of molecular gas densities (n_H2 ∼ 10^3 − 5 cm⁻³). The SFR is derived from Paα line emission imaged by HST/NICMOS. We analysed the influence of the dynamical environment on different formulations of SF relations in the SB ring and compared our results with the general predictions of density-threshold and turbulent SF models.

The main results of this paper are summarised as follows:

• We derived spatially resolved KS laws for a set of eight spatial resolutions, ranging from ∼40 pc up to ∼700 pc. We studied how these relations change depending on the adopted aperture sizes and on the choice of molecular gas tracer. For a given spatial resolution the correlation parameters derived from the high-density tracers (CO(3–2), HCN(1–0) and HCO⁺(1–0)) are about a factor of two to three larger than that derived from CO(1–0).

• The KS correlations lose statistical significance below the critical spatial scale of ≈300-400 pc, common for all gas tracers. For spatial scales ≥300 − 400 pc the correlation improves monotonically as a function of the aperture size for all gas tracers, with the exception of CO(1–0). While dust continuum emission shows a behaviour similar to other high-density tracers, the dust-based KS correlation is significant already at 100 pc scales.

• NGC 1068 lies within the general Σ_SFR − Σ_dense linear KS relationship derived from a compilation of data obtained in other galaxies. In particular, the location of NGC 1068 in the KS plot is intermediate between that of normal galaxies and that of (U)LIRGs, a result that underlines the relatively extreme conditions in the SB ring.

• The SFE of the dense molecular gas, defined as SFE_dense ≡ Σ_SFR/Σ_dense, shows a scattered distribution as a function of the HCN luminosity at the initial resolution of ALMA (∼56 pc) around a mean value of ≃0.01 Myr⁻¹. However, we find evidence of a significant environmental dependence of SFE_dense reflected in the existence of systematic trends across the different regions of the SB ring, which are inconsistent with the predictions of density-threshold models.

• With the aim of resolving the degeneracy associated with the SFE_dense − L′(HCN) plot, we explored an alternative prescription for SF relations, which includes the dependence of SFE_dense on the boundedness of the gas that is measured by the parameter b defined as $b\equiv {\mathrm{\Sigma }}_{\mathrm{dense}}/{\sigma }^2\propto {\alpha }_{\mathrm{vir}}^{-1}$. We identified two branches in the version of the SFE_dense − b plot derived for an averaging scale of 400 pc. The two branches correspond to two dynamical environments defined by their proximity to the region where the SB ring is connected to the stellar bar of NGC 1068.

• We studied the trends in SFE_dense as a function of the dense gas fraction (F_dense), the stellar mass surface density (Σ_star), and the velocity dispersion (σ) in the SB ring using different averaging scales. We find that SFE_dense correlates both with F_dense and Σ_star. However, we find no significant correlation of SFE_dense with σ. These results differ to a large extent from those derived from previous kiloparsec-scale studies of galaxy disks and high-resolution (∼100 pc) observations of M 51.

• The trends in the boundedness parameter in the SB ring are mostly driven by the observed boost in F_dense in the bar-ring interface region. The results reflect a significant compression of molecular gas as a result of an enhanced rate of cloud-cloud collisions in the bar-ring interface on the one hand and an efficient dissipation of turbulence on the other.

The outcome of our work emphasises the relevance of different dynamical environments to modulating SFE_dense in galaxy disks. Taken at face value, this result is inconsistent with the predictions of density-threshold models, which foresee the existence of a nearly constant canonical value for SFE_dense. There is mounting observational evidence that SFE_dense can show systematic trends based on studies carried out in both our Galaxy and nearby galaxies. However, the theoretical scenarios advanced to explain the variations in SFE_dense reported in the literature are diverse. This diversity reflects the variety of trends observed in SFE_dense as a function of key physical variables such as F_dense, Σ_star, and σ in galaxy disks. Similarly to the case of the barred galaxy NGC 3627 discussed by Beuther et al. (2017), the enhanced SF activity in the bar-ring interface of NGC 1068 could be explained by a dynamical configuration in which the bar and the two-arm spiral feature are independent m = 2 modes rotating at different pattern speeds. Further studies are required to explore the frequency of this type of dynamical decoupling between bars and spirals and their role at setting the SFE of dense molecular gas in galaxy disks in general and in galactic rings in particular.

Acknowledgments

We thank the referee for the useful comments and suggestions. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2013.1.00055.S and #2011.0.00083.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan), 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 used observations made with the NASA/ESA Hubble Space Telescope, and obtained from the Hubble Legacy Archive, which is a collaboration between the Space Telescope Science Institute (STScI/NASA), the Space Telescope European Coordinating Facility (ST-ECF/ESA), and the Canadian Astronomy Data centre (CADC/NRC/CSA). M.S.G. acknowledges support from the Spanish Ministerio de Economía y Competitividad through the grants BES-2016-078922 and ESP2015-68964-P. S.G.B. and M.Q. acknowledge support from the research project PID2019-106027GA-C44 of the Spanish Ministerio de Ciencia e Innovación. A.F. acknowledges support from the research project PID2019-106235GB-I00. A.A.H., S.G.B. and A.U. work was funded by grant PGC2018-094671-B-I00 funded by MCIN/AEI/ 10.13039/501100011033 and by ERDF A way of making Europe. A.A.H’s work was done under project No. MDM-2017-0737 Unidad de Excelencia “María de Maeztu”- Centro de Astrobiología (INTA-CSIC). M.P.S. acknowledges support from the Comunidad de Madrid through the Atracción de Talento Investigador Grant 2018-T1/TIC-11035 and PID2019-105423GA-I00 (MCIU/AEI/FEDER,UE). L.C. and M.S.G. acknowledge support from the research project PID2019-106280GB-100.

References

Appendix A: Determination of the missing flux in HCN and CO maps

We estimated the missing flux in the HCN(1–0) map obtained by ALMA and in the CO(1–0) map of the IRAM array (PdBI) (Schinnerer et al. 2000). To this aim, we measured the spatially integrated fluxes in both lines using a set of four apertures covering the region occupied by the SB ring of NGC 1068. The interferometer fluxes are compared with the corresponding single-dish fluxes measured from observations obtained by the IRAM 30-m telescope (Usero et al. 2021, priv. comm.). Figure A.1 illustrates the apertures chosen overlaid on the HCN (upper panel) and CO (lower panel) interferometer images of the galaxy. The sizes of the apertures used to derive the spatially integrated fluxes of the interferometer maps correspond to the single-dish beams of the IRAM 30-m telescope, which are ∼30″ and ∼24″ for HCN and CO, respectively.

Fig. A.1. Selection of four apertures to study the missing flux in the HCN and CO maps of NGC 1068. Upper panel: Overlay of the four ∼30″ size apertures (black circles, identified by the numerical labels 1-4) selected to cover the SB ring of NGC 1068 on the ALMA HCN(1–0) map (colour scale and contours as in Fig. 1), used in Appendix A to estimate the missing flux correction. Lower panel: Same as the upper panel but showing the overlay of the four ∼12″ size apertures (black circles, identified by the numerical labels 1-4) selected to cover the SB ring of the galaxy on the PdBI CO(1–0) map of Schinnerer et al. (2000). The colour scale is given in Jy km s−1 beam−1 units and the position of the AGN is identified by the black star markers in both maps. The magenta filled ellipses at the bottom-left in the two panels represent the beam sizes.

Figures A.2 and A.3 compare the spectra of HCN and CO obtained by integrating the emission of the lines inside the four selected regions shown in Fig. A.1 with the IRAM 30-meter spectra. Based on this comparison we estimate that on scales of ∼30″ (∼2 kpc) approximately 20 − 25% of the total flux gets filtered in the HCN ALMA map. A similar estimate for the CO map yields a flux filtering of about 40-45% on scales ∼24″ (∼1.7 kpc). This is about a factor of two larger than the value estimated by Schinnerer et al. (2000), who compared the fluxes measured with the PdBI with those measured with the combined BIMA array and Kitt Peak single-dish telescope datasets on scales ∼55″ (3.8 kpc).

Fig. A.2. Comparison of the HCN(1–0) emission observed by ALMA (blue histograms) and the IRAM 30-m telescope (black histograms) inside the four selected regions identified in Fig. A.1 (labelled as 1, 2, 3, and 4 in each panel).

Fig. A.3. Same as Fig. A.2, but showing the comparison of the CO(1–0) emission observed by the IRAM-PdBI (blue histograms) and the IRAM 30-m telescope (black histograms) inside the four selected regions identified in Fig. A.1 (labelled as 1, 2, 3, and 4 in each panel).

As in this paper we analyse spatial scales that are much smaller than the critical scales reported above (≈ 2-4 kpc), we can foresee that the missing flux factors will be significantly less than 25-45% for the scales that are relevant to our study.

Appendix B: Paα analysis as SFR tracer

To validate our working assumption, namely that the Paα image used in this work to derive the distribution of recent SF in NGC 1068 does not require a significant correction by dust extinction, we compared the Paα fluxes measured in the HST/NICMOS image of the galaxy with the Hα fluxes measured using the ground-based image of Díaz et al. (2000) over a representative number of hot spots of the SB ring. Figure B.1 shows the position of the eight regions considered in this comparison selected from Table 4 of Díaz et al. (2000). We derived the Paα spatially-integrated fluxes inside the eight selected regions using the sizes of Hα SF knots listed in Table 4 of Díaz et al. (2000) and obtained the Hα/Paα flux ratios for each aperture (listed in Table B.1). To estimate the colour excess (E(B − V)) and extinction (A_Paα) values we used the prescriptions of Calzetti (2001) (see also Calzetti et al. 1994, 2000) and the standard equation:

$$\begin{array}{c}\hfill F_{\mathrm{obs}}(\lambda )=F_{\mathrm{int}}(\lambda ){10}^{-0.4k(\lambda )E(B-V)},\end{array}$$(B.1)

Fig. B.1. HST/NICMOS Paα image of NGC 1068. The white labels identify the regions selected to estimate the extinction of Paα, according to the nomenclature used by Díaz et al. (2000).

where F_obs and F_int are the observed and intrinsic fluxes for the lines, respectively, and E(B − V) is the colour excess. We assumed that the differential extinction between the two lines is determined by k(Hα)−k(Paα) = 2.104 (Calzetti 2001). Table B.1 lists the Hα/Paα line ratios as well as the values estimated for E(B − V) and A_Paα for the eight selected regions.

The median value for the distribution of Hα/Paα ratios listed in Table B.1 is ∼6.8. This is close to the expected value corresponding to the canonical case B recombination: ∼7.8. The implied mean extinction correction at 1.875 μm is A_Paα = 0.03 mag, which translates into very low dust opacities for the ensemble of the studied SF knots. Although in regions 18, 19 and 20, A_Paα is seems to reach values ∼0.1 − 0.2, the implied opacities are still compatible with optically thin emission for Paα. Moreover, assuming that the flux scale is uncertain to ±20% due to absolute calibration errors, we can therefore conclude that Paα emission is not significantly affected by dust extinction in the SB ring of NGC 1068.

Appendix C: Kennicutt–Schmidt plots

We represent the KS plots derived from the different tracers for the range of scales used in this work (Figs. C.1, C.2, C.3 and C.4).

Fig. C.1. Same as Fig. 5, but showing all the spatial scales used in this work for HCO+(1-0). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.

Fig. C.2. Same as Fig. 5, but showing all the spatial scales used in this work for CO(3-2). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.

Fig. C.3. Same as Fig. 5, but showing all the spatial scales used in this work for dust continuum emission. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.

Fig. C.4. Same as Fig. 5, but showing all the spatial scales used in this work for CO(1-0). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex and ±0.10 dex, respectively.

Appendix D: Alternative formulation of the intensity-weighted average of the b parameter

Leroy et al. (2017) used a definition of the intensity-weighted average of the b parameter different to the one used in Sect. 6.1. In particular, instead of deriving the averages of b as defined in Eq. (10), they performed the averages of Σ_gas and σ² separately and defined ⟨b⟩≡⟨Σ_gas⟩/⟨σ²⟩. We derived a new version of Eq. (11) using the definition of Leroy et al. (2017) particularised for Σ_dense and our set of two apertures (100 pc and 400 pc). The results obtained with this definition, shown in Fig. D.1 are virtually identical to those shown in Fig. 11. In particular, we obtain (anti) correlation Spearman rank parameters ρ_sp = −0.42 and −0.60 for ΔA = 100 pc and 400 pc, respectively, and associated two-sided p-values < 1%. Furthermore, after applying the MARS routine to the 400 pc-scale scatter plot we find a turnover at log₁₀(⟨b⟩) = 0.18 M_⊙pc⁻²(km s⁻¹)⁻² and a slope ≃ − 1.4 for ⟨b⟩ values beyond the turnover.

Fig. D.1. Same as Fig. 11 but derived using the alternative formulation of the intensity-weighted average of the b parameter used by Leroy et al. (2017).

Appendix E: Trends as a function of the boundedness of the gas from HCO⁺

We have studied in Sect. 6.1 how the depletion time of the dense molecular gas estimated from HCN changes as a function of the self-gravity of the gas (b paramater). We now analyse this with the other dense gas tracer that we use in this work, HCO⁺. Figure E.1 represents the depletion time as a function of the b parameter at 56 pc, 100 pc and 400 pc. We observe a similar trend to Fig. 11, obtaining (anti) correlation Spearman rank parameters ρ_sp = -0.07 (p-value = 0.06), -0.32 (< 1%) and -0.65 (< 1%) for 56 pc, 100 pc and 400 pc, respectively. At 400 pc scales, we identify a turnover in the plot located around log₁₀(⟨b⟩) ≃ 0.2 M_⊙pc⁻²(km s⁻¹)⁻² and a slope ≃ -0.91 for ⟨b⟩ values beyond the turnover.

Fig. E.1. Same as Fig. 11 but estimated from HCO+.

Fig. E.2. Same as Fig. 12 but estimated from HCO+. Left panel: contour levels for log10($T_{\mathit{dep}}^{\mathit{dense}}$) from 1.90 to 2.50 in steps of 0.045 (in Myr). Right panel: contour levels for log10(⟨b⟩) from 0.02 to 0.64 to 0.025 (in M⊙pc−2(km s−1)−2).

Figure E.1 shows the spatial distribution of ${\mathrm{T}}_{\mathrm{dep}}^{\mathrm{dense}}$ and ⟨b⟩ derived for ΔA = 400 pc in the SB ring of NGC 1068 from HCO⁺. We observe that higher (lower) SFE_dense (${\mathrm{T}}_{\mathrm{dep}}^{\mathrm{dense}}$) values are located closer to the bar-ring interface. We find that the HCN and HCO⁺ show similar correlations in the SB ring of NGC 1068.

All Tables

All Figures

Fig. 1. HCN(1–0) and HCO+(1–0) integrated intensity maps of NGC 1068. Left panel: velocity-integrated intensity map of HCN(1–0) of NGC 1068 obtained with ALMA. The map is shown in a logarithmic colour scale with contour levels of 3σ, 5σ, 7σ, 12σ, 24σ, 50σ, 120σ, and 200σ, where 1σ = 0.028 Jy beam−1 km s−1. We highlight the location of several representative regions of the emission, namely the CND, the bar, and the SB ring. Right panel: same as the left panel but for the HCO+(1–0) line map, with the same list of contour levels but here truncated at 120σ. The position of the AGN is identified by the dark grey star marker. The black bar in the top-right corner of each panel shows the spatial scale (in pc). The magenta filled ellipses at the bottom-left corners of the panels represent the ALMA beam size (1$\stackrel{″}{.}$0 × 0$\stackrel{″}{.}$6 at PA = 69°).

In the text

Fig. 2. HCN(1–0)/HCO+(1–0) brightness temperature ratio map (Rdense) derived at the common spatial resolution of ≃56 pc. The brown circle of 8″ in radius (≃560 pc) locates the inner region of the galaxy disk where the molecular outflow signature has been identified in the kinematics of molecular gas (García-Burillo et al. 2014, 2019). This region, purposely excluded from our analysis of SF relations, shows significantly higher line ratios (Rdense ≃ 1.5 − 3.5) compared to the SB ring (Rdense ≃ 0.5 − 1.5). Other symbols are the same as in Fig. 1.

In the text

Fig. 3. HST/NICMOS continuum-subtracted Paα image of NGC 1068. The map is shown in colour scale with contour levels of 3σ, 9σ, 15σ, 32σ, 50σ, 120σ, and 700σ, where 1σ = 0.212 × 10−16 erg s−1 cm−2 pixel−1. The colour scale range is shown in units of 10−16 erg s−1 cm−2 pixel−1. The (truncated) square region identifies the outer edge of the HST/NICMOS FOV. The angular resolution of the image is 0$\stackrel{″}{.}$26 × 0$\stackrel{″}{.}$26. Other symbols are the same as in Fig. 2.

In the text

Fig. 4. Overlay of the HST/NICMOS Paα image on the molecular gas tracer maps. Left panel: overlay of the HST/NICMOS Paα image (contours) on the ALMA HCN(1–0) map (colour scale). Contours have a logarithmic spacing from 3σ to 1536σ, in steps of 0.3 dex, where 1σ = 0.180−16 erg s−1 cm−2 pixel−1. The colour scale is shown in units of Jy beam−1 km s−1. Middle panel: same as the left panel but showing the overlay of the HST/NICMOS Paα contours on the ALMA HCO+(1–0) map (colour scale). Right panel: same as the left panel but showing the overlay of the HST/NICMOS Paα contours on the PdBI CO(1–0) map of Schinnerer et al. (2000) (colour scale). Contour spacing is the same as in the left and right panels, but with 1σ = 0.162 × 10−16 erg s−1 cm−2 pixel−1. The magenta filled ellipses at the bottom-left corners in all panels represent the beam sizes of the molecular line maps: 1$\stackrel{″}{.}$0 × 0$\stackrel{″}{.}$6 at PA = 69° (left and middle panels) and 1$\stackrel{″}{.}$84 × 1$\stackrel{″}{.}$09 at PA = 24° (right panel). The HST/NICMOS Paα images in each panel have been degraded to the corresponding spatial resolutions of the molecular gas tracers used in this comparison. Other symbols are the same as in Fig. 2.

In the text

Fig. 5. SFR surface density (ΣSFR) as a function of the dense molecular gas surface density (Σdense) derived from HCN(1–0) for the different resolutions explored in this work: 56 pc, 100 pc, 200 pc, 300 pc, 400 pc, 500 pc, 600 pc and 700 pc. The red lines show the ODR fits to the data. The dashed lines identify the correlations that are not found to be statistically significant. We only show on both axes the data points above 3σ that are considered in the fits. We indicate the Spearman’s rank correlation coefficients (ρsp) and the power-law indexes (N) of the best fits. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.

In the text

Fig. 6. Pearson and Spearman correlation coefficients (ρps: left panel; ρsp: right panel) of the KS laws derived for the different gas tracers used in this work as a function of the spatial resolution. Error bars represent the 67% confidence interval around the mean values. For spatial scales larger than the range identified by the grey-coloured region (≤350 pc) all the gas tracers (leaving aside CO(1–0)) show statistically significant correlations in their KS laws, defined by p-values < 1% and ρps, ρsp > 0.4.

In the text

	Fig. 7. Power-law indexes (N) of the KS relations obtained from the ODR fit to the data for the different gas tracers used in this work as a function of the spatial resolution adopted. Symbols and colour codes are the same as in Fig. 6.
In the text

Fig. 8. Comparison of the KS laws derived from HCN(1–0) in NGC 1068 with those obtained in different populations of galaxies for a range of spatial scales: global measurements (around a few kiloparsec) of the disks of SFGs and (U)LIRGs (García-Burillo et al. 2012), ≃0.3−1 kiloparsec-scale regions of the disks of nearby SFGs (Murphy et al. 2015; Usero et al. 2015; Bigiel et al. 2016), and 0.1−0.2 kpc size regions in the disk of M 51 (Chen et al. 2017; Viaene et al. 2018; Querejeta et al. 2019). For the sake of a better comparison with the reference sample, the scatter plots for NGC 1068 (magenta circles and isodensity contours) are displayed for three spatial resolutions: 56 pc (left panel), 100 pc (middle panel), and 700 pc (right panel). The scatter plot for M 51 at 100 pc resolution is displayed by the green circles and isodensity contours. The open squares stand for the median values of the distributions in NGC 1068 and M 51.

In the text

Fig. 9. Distribution of the HCN and Paα peaks and the apertures showing the highest SFEdense in NGC 1068. Left panel: overlay of the HST/NICMOS Paα image (contours) on the HCN(1–0) map (grey scale). The contours and grey scale are the same as in Fig. 4. The set of 56 pc size circular apertures used to extract the fluxes of Paα and HCN(1–0) are colour-coded differently, depending on whether they are centred on HCN peaks (in blue), on Paα peaks (in red), or on the 10% of the apertures showing the highest values of SFEdense (in green). The black dots identify the centres of all the Nyquist-sampled 56 pc size apertures where there are reliable estimates of SFEdense across the SB ring. Right panel: distribution of SFEdense values (in Myr−1-units) as a function of LHCN (in K km s−1 pc2-units) in the SB ring of NGC 1068. Symbols are the same as in left panel. Vertical and horizontal error bars at the lower-right corner of the panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

In the text

	Fig. 10. Same as Fig. 9 but replacing HCN(1–0) with HCO+(1–0) along the x axis of the right panel.
In the text

Fig. 11. Depletion time of the dense molecular gas estimated from HCN, $T_{\mathrm{dep}}^{\mathrm{dense}}\equiv {\mathrm{\Sigma }}_{\mathrm{dense}}$/ΣSFR, as a function of the self-gravity of the gas, measured by the parameter b ≡ Σdense/σ2, at different spatial scales (magenta circles and isodensity contours): 56 pc, namely the initial resolution of ALMA observations (left panel), 100 pc (middle panel), and 400 pc (right panel). A significant trend in $T_{\mathrm{dep}}^{\mathrm{dense}}$ as a function of b, identified in the middle and right panels, indicates a higher rate of SF per unit mass of dense gas for regions with stronger self-gravity. Values of ⟨b⟩ for the 100 pc and 400 pc apertures were derived from an intensity-weighted average of b evaluated at the Nyquist-sampled grid of points for each aperture. Values of $T_{\mathrm{dep}}^{\mathrm{dense}}$ were derived using a Gaussian-weighted average. Black circles show median $T_{\mathrm{dep}}^{\mathrm{dense}}$, and error bars indicate the rms scatter in the bins. The green lines in the right panel identify the two regimes in the $T_{\mathrm{dep}}^{\mathrm{dense}}-⟨b⟩$ plot found by the MARS algorithm to fit the data using 400 pc as averaging scale. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

In the text

Fig. 12. Maps of the depletion time and the b parameter in the SB ring of NGC 1068. Left panel: map of the depletion time of the dense molecular gas ($T_{\mathrm{dep}}^{\mathrm{dense}}$) obtained after using 400 pc as an averaging scale. Contour levels for log10($T_{\mathrm{dep}}^{\mathrm{dense}}$) go from 1.50 to 2.22 in steps of 0.045 (in Myr). Right panel: map of the intensity-weighted average of the b parameter (⟨b⟩) in the SB ring of NGC 1068 obtained after using 400 pc as an averaging scale. Contour levels for log10(⟨b⟩) go from 0.075 to 0.45 in steps of 0.025 (in M⊙ pc−2 (km s−1)−2). The Δα and Δδ offsets in arc seconds are relative to the phase tracking centre of ALMA. The ellipse identifies the extent of the nuclear stellar bar in the disk. We highlight the approximate extent of the bar-ring interface region in both panels.

In the text

	Fig. 13. HCN(1–0)/CO(1–0) brightness temperature ratio map (RHCN/CO; colour scale and contours) derived at the common (lower) spatial resolution of the CO(1–0) observations of Schinnerer et al. (2000) (≃100 pc). The brown circle of 8″ in radius (≃560 pc) locates the inner region excluded from our analysis of SF relations.
In the text

	Fig. 14. Same as Fig. 13 but showing the CO(3–2)/CO(1–0) brightness temperature ratio map (R3 − 2/1 − 0; colour scale and contours) derived at the common (lower) spatial resolution of the CO(1–0) observations.
In the text

Fig. 15. Star-formation efficiency of the dense gas (SFEdense) as a function of dense gas fraction (Fdense: estimated from the HCN(1–0)-to-CO(1–0) line ratio of Fig. 13) in NGC 1068 using two averaging scales: 100 pc (left panel) and 400 pc (right panel). Magenta circles correspond to data measurements evaluated at the Nyquist-sampled grid of points after convolution with the appropriate Gaussian for each aperture. Black circles show median SFEdense values, and error bars indicate the rms scatter in the bins. Vertical and horizontal error bars at the lower-right corner of both panels account for the typical uncertainties, which amount to ±0.13 dex.

In the text

Fig. 16. Overlay of the ALMA HCN(1–0) map (contours) on the HST/NICMOS F190N continuum image at 1.9 μm HST (colour scale as shown in Jy). The contour spacing is 3σ, 7σ, 12σ, 24σ, 50σ, and 120σ, where 1σ = 0.028 Jy beam−1 km s−1. The white circle of r = 8″ locates the inner region of the galaxy disk, which was excluded of our analysis. The HST/NICMOS F190N image has been degraded to the spatial resolution of HCN(1–0) (1${}^{″}\times 0\stackrel{″}{.}6\approx 50$ pc).

In the text

Fig. 17. Star-formation efficiency of the dense gas (SFEdense) as a function of the stellar mass surface density (Σstar) in the SB ring of NGC 1068 at the initial resolution of 56 pc (left panel) and at two averaging scales: 100 pc (middle panel) and 400 pc (right panel). Magenta circles and isodensity contours correspond to data measurements evaluated at the Nyquist-sampled grid of points after convolution with the appropriate Gaussian for each aperture. Black circles show median SFEdense values, and error bars indicate the rms scatter in the bins. Symbols and isodensity contours in the middle panel identify the data obtained in other galaxies. The dashed orange line represents the fit of Usero et al. (2015). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

In the text

Fig. 18. Fraction of dense gas (Fdense) as a function of the stellar mass surface density (Σstar) in the SB ring of NGC 1068 at the two averaging scales: 100 pc (left panel) and 400 pc (right panel). Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.09 dex, respectively.

In the text

Fig. 19. Star-formation efficiency of the dense gas (SFEdense) as a function of the velocity dispersion of HCN (σ) in NGC 1068 for three averaging scales: 56 pc (left panel), 100 pc (middle panel), and 400 pc (right panel). Measurements at 100 pc scales are compared to literature data for other galaxies. Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.06 dex, respectively.

In the text

Fig. 20. Same as the middle panel of Fig. 19 but adding a comparison with the trends of SFEdense as a function of the velocity dispersion of the gas (σ) predicted by the model of Meidt et al. (2020). This model foresees different locations of galaxies in the SFEdense − σ parameter space depending on the fraction of dense gas (Fdense) and, also, on the total molecular gas column densities (Σmol). Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner account for the typical uncertainties, which amount to ±0.13 dex and ±0.06 dex, respectively.

In the text

Fig. 21. Dense gas fraction (Fdense) as a function of velocity dispersion in the SB ring of NGC 1068 at the two averaging scales: 100 pc (left panel) and 400 pc (right panel). Measurements at 100 pc scales are compared to literature data for other galaxies. Symbols are the same as in Fig. 17. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.13 dex and ±0.06 dex, respectively.

In the text

	Fig. 22. Star formation rate surface density (ΣSFR) as a function of dense gas fraction (Fdense) in NGC 1068 measured at 100 pc scales compared to the data obtained in other galaxies (symbols as in Fig. 8). Vertical and horizontal error bars at the lower-right corner account for the typical uncertainties, which amount to ±0.09 dex and ±0.13 dex, respectively.
In the text

	Fig. 23. Same as Fig. 12 but showing the map of the velocity dispersion (σ) obtained using 400 pc as an averaging scale. Contour levels for σ go from 14.0 to 19.0 in steps of 0.5 km s−1.
In the text

Fig. A.1. Selection of four apertures to study the missing flux in the HCN and CO maps of NGC 1068. Upper panel: Overlay of the four ∼30″ size apertures (black circles, identified by the numerical labels 1-4) selected to cover the SB ring of NGC 1068 on the ALMA HCN(1–0) map (colour scale and contours as in Fig. 1), used in Appendix A to estimate the missing flux correction. Lower panel: Same as the upper panel but showing the overlay of the four ∼12″ size apertures (black circles, identified by the numerical labels 1-4) selected to cover the SB ring of the galaxy on the PdBI CO(1–0) map of Schinnerer et al. (2000). The colour scale is given in Jy km s−1 beam−1 units and the position of the AGN is identified by the black star markers in both maps. The magenta filled ellipses at the bottom-left in the two panels represent the beam sizes.

In the text

	Fig. A.2. Comparison of the HCN(1–0) emission observed by ALMA (blue histograms) and the IRAM 30-m telescope (black histograms) inside the four selected regions identified in Fig. A.1 (labelled as 1, 2, 3, and 4 in each panel).
In the text

	Fig. A.3. Same as Fig. A.2, but showing the comparison of the CO(1–0) emission observed by the IRAM-PdBI (blue histograms) and the IRAM 30-m telescope (black histograms) inside the four selected regions identified in Fig. A.1 (labelled as 1, 2, 3, and 4 in each panel).
In the text

	Fig. B.1. HST/NICMOS Paα image of NGC 1068. The white labels identify the regions selected to estimate the extinction of Paα, according to the nomenclature used by Díaz et al. (2000).
In the text

	Fig. C.1. Same as Fig. 5, but showing all the spatial scales used in this work for HCO+(1-0). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.
In the text

	Fig. C.2. Same as Fig. 5, but showing all the spatial scales used in this work for CO(3-2). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.
In the text

	Fig. C.3. Same as Fig. 5, but showing all the spatial scales used in this work for dust continuum emission. Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex on both axes.
In the text

	Fig. C.4. Same as Fig. 5, but showing all the spatial scales used in this work for CO(1-0). Vertical and horizontal error bars at the lower-right corner of each panel account for the typical uncertainties, which amount to ±0.09 dex and ±0.10 dex, respectively.
In the text

	Fig. D.1. Same as Fig. 11 but derived using the alternative formulation of the intensity-weighted average of the b parameter used by Leroy et al. (2017).
In the text

	Fig. E.1. Same as Fig. 11 but estimated from HCO+.
In the text

	Fig. E.2. Same as Fig. 12 but estimated from HCO+. Left panel: contour levels for log10($T_{\mathit{dep}}^{\mathit{dense}}$) from 1.90 to 2.50 in steps of 0.045 (in Myr). Right panel: contour levels for log10(⟨b⟩) from 0.02 to 0.64 to 0.025 (in M⊙pc−2(km s−1)−2).
In the text