A diversity of starburst-triggering mechanisms in interacting galaxies and their signatures in CO emission

TL;DR

This study uses high-resolution simulations of interacting galaxies to explore how various physical processes influence star formation, CO emission, and gas properties, revealing complex spatial and temporal variations that impact observational interpretations.

Contribution

It provides new insights into the diverse mechanisms triggering starbursts and their signatures in CO emission, emphasizing the importance of spatially resolved analysis in galaxy evolution studies.

Findings

Physical processes vary over time and space, affecting observables.

Excess dense gas leads to reduced depletion times and deviations from classical star formation.

Spatial variations in the CO-to-H2 conversion factor are significant.

Abstract

The physical origin of enhanced star formation activity in interacting galaxies remains an open question. Knowing whether starbursts are triggered by an increase of the quantity of dense gas or an increase of the star formation efficiency would improve our understanding of galaxy evolution and allow to transpose the results obtained in the local Universe to high redshift galaxies. In this paper, we analyze a parsec-resolution simulation of an Antennae-like model of interacting galaxies. We find that the interplay of physical processes has complex and important variations in time and space, through different combinations of mechanisms like tides, shear and turbulence. These can have similar imprints on observables like depletion time and CO emission. The densest gas only constitutes the tail of the density distribution of some clouds, but exists in large excess in others. The super-linearity of the star formation rate dependence on gas density implies that this excess translates into a reduction of depletion times, and thus leads to a deviation from the classical star formation regime, visible up to galactic scales. These clouds are found in all parts of the galaxies, but their number density varies from one region to the next, due to different cloud assembly mechanisms. Therefore, the dependence of cloud and star formation-related quantities (like CO flux and depletion time) on the scale at which they are measured also varies across the galaxies. We find that the $\alpha_{\rm CO}$ conversion factor between the CO luminosity and molecular gas mass has even stronger spatial than temporal variations in a system like the Antennae. These results raise a number of cautionary notes for the interpretation of observations of unresolved star-forming regions, but also predict that the diversity of environments for star formation will be better captured by the future generations of instruments.

Equations

$${\sigma }^2=\ln \left(1+b^2{\mathcal{M}}^2\right),$$

Full text

1 Introduction

In the local Universe, galaxy interactions and mergers are often associated with an enhancement of the star formation activity, particularly visible in (ultra) luminous infrared galaxies (ULIRGs, Sanders & Mirabel 1996; Scudder et al. 2012; Knapen et al. 2015). Although this is now routinely reproduced in simulations (e.g. Karl et al., 2010; Teyssier et al., 2010; Moreno et al., 2013; Renaud et al., 2014), the underlying physical reasons and the main driver(s) of such starbursts remain uncertain. Understanding the physics of enhanced star formation is key to establish whether and how the known empirical and theoretical relations between star formation and its triggers derived in the local Universe can be transposed to the extreme conditions at high redshift, during the first phases of galaxy formation.

To date, it is not clear whether the bursts of star formation rate (SFR) result from an enhancement of the quantity of dense gas due to a compression from galactic scale mechanisms (e.g. Combes et al., 1994; Casasola et al., 2004; Kaneko et al., 2017), or from an increase of the efficiency of the formation process itself111Not to be confused with the dimension-less small-scale star formation efficiency (SFE) describing the stellar mass produced from a given amount of gas, typically at sub-cloud scales ($\lesssim 1\mathrm{\leavevmode\nobreak\ pc}$) and often implemented in sub-grid recipes of star formation. For the sake of clarity, we do not use the term “efficiency” to describe the galactic-scale quantity but rather use the depletion time, which is directly its inverse. The small scale SFE and the galactic scale depletion time are not directly related (see Semenov et al., 2018). (see Solomon & Sage 1988; Graciá-Carpio et al. 2008; Michiyama et al. 2016). The former is often quantified by estimating the molecular gas mass, while the latter relates to the depletion time of the gas reservoir ($t_{\mathrm{dep}}=M_{\rm gas}/{\rm SFR}$). A number of observations of local interacting systems have quantified the variations of these two quantities to determine the main driver of starbursts, often as functions of the interaction parameters. The emergent consensus depicts a mixed influence of both effects, with the increase of dense gas content being the dominant one (Pan et al., 2018).

In the Schmidt-Kennicutt plane (i.e. the relation between the surface densities of gas and of SFR), an increase of the dense gas mass leads to the expected increase of SFR along the canonical relation ($\Sigma_{\mathrm{gas}}\propto\Sigma_{\mathrm{SFR}}^{N}$, with $N\approx 1.4\textrm{\,--\,}2.3$, Kennicutt 1998; Kennicutt & Evans 2012; Kravtsov 2003; Elmegreen & Elmegreen 2005; Semenov et al. 2018), i.e. at a roughly uniform depletion time ($\sim 1\mathrm{\leavevmode\nobreak\ Gyr}$, Bigiel et al. 2010; Saintonge et al. 2011; Leroy et al. 2013). Such a high dense gas mass is what is observed in e.g. gas-rich turbulent disk galaxies at redshifts $z\gtrsim 1$ (e.g. Tacconi et al., 2010). However, a decrease of the depletion time as observed, e.g. in LIRGs (Klaas et al., 2010; Martinez-Badenes et al., 2012; Nehlig et al., 2016) translates into a deviation from this relation. A duality of star formation regimes has been proposed by Daddi et al. (2010b) to distinguish galaxies with long and short depletion times, respectively disks and mergers, but independently of their SFR. To reach this conclusion, these authors had to assume a conversion factor $\alpha_{\mathrm{CO}}$ between the CO luminosity and the molecular mass (due to the difficulty to detect H2, see Bolatto et al. 2013). Uncertainties hinder the choice of this factor (Zhu et al., 2003; Graciá-Carpio et al., 2008; Pereira-Santaella et al., 2016). The value of $0.8\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}\mathrm{\leavevmode\nobreak\ K}^{-1}\mathrm{\leavevmode\nobreak\ km}^{-1}\mathrm{\leavevmode\nobreak\ s}\mathrm{\leavevmode\nobreak\ pc}^{-2}$ adopted by Daddi et al. for the mergers is in sharp contrast with that of the disks (4.3 for the Milky Way) and has been highly debated, since an under-estimated $\alpha_{\mathrm{CO}}$ would artificially accentuate or even create the differences between the regimes in depletion times. However, the decrease of $t_{\mathrm{dep}}$ is retrieved in parsec-resolution simulations of interacting galaxies in which the transitions between the disk and the merger regimes (and back) can be resolved, without any assumption on $\alpha_{\mathrm{CO}}$ (Renaud et al., 2014; Kraljic, 2014; Fensch et al., 2017; Renaud et al., 2018). Furthermore, by post-processing similar simulations to derive $\alpha_{\mathrm{CO}}$, Renaud et al. (2019) showed that the choice a low $\alpha_{\mathrm{CO}}$ for ULIRG systems is likely justified (see also Narayanan et al., 2011), for a few $10\mathrm{\leavevmode\nobreak\ Myr}$ after the peak(s) of the their SFR. The positions of starbursting systems (at least in the local Universe) in the Schmidt-Kennicutt plane thus confirm both an elevation of the dense gas content, and a decrease of the depletion time. (Note that this effect gets hidden in variants of this diagram when $\Sigma_{\mathrm{gas}}$ is normalized by the free-fall time, Krumholz et al. 2012; Federrath 2013; Salim et al. 2015.)

Variations of the depletion time have been reported observationally across galactic disks (e.g. Meidt et al. 2013; Bigiel et al. 2015; Usero et al. 2015; Pereira-Santaella et al. 2016; Saito et al. 2016; Leroy et al. 2017; Tomičić et al. 2018) and along the course of interactions either theoretically (e.g. Renaud et al., 2014), or observationally by considering statistical samples of galaxies at different evolutionary stages (e.g. Feldmann et al., 2012; Violino et al., 2018; Pan et al., 2018).

For instance, the evolution of modeled interacting systems in the Schmidt-Kennicutt plane reveals that the successive encounters do not induce the same effects (see e.g. Figure 3 of Renaud et al. 2014 for an Antennae-like system, and Figure 11 of Renaud et al. 2018 for a Cartwheel-like galaxy). Distant passages (usually at an early stage of the interaction, before strong dynamical friction has reduced the orbital energy) trigger a violent increase of the star formation activity over large volumes, while the global surface density of gas remains almost unchanged (measured at a scale of $\sim 100\mathrm{\leavevmode\nobreak\ pc}$). However, closer, penetrating encounters (usually at later stages, moments before the coalescence) yield both an increase of the gas and the star formation rate densities, due to most of the activity being concentrated in the galactic nuclei. Such differences are also visible in the star cluster formation rates, where clusters are preferentially formed in early passages instead of at coalescence (Renaud et al., 2015a).

These variations in both space and time relate to the extended star formation, and in particular the off-nuclear activity (Barnes, 2004; Wang et al., 2004; Cullen et al., 2006; Smith et al., 2008; Hancock et al., 2009; Chien & Barnes, 2010), and suggest that different physical mechanisms trigger starburst events, through an increase of the dense gas mass, a decrease of the depletion time, or both, and that their interplay and relative importance vary in space and time, as a function of the properties of the galactic encounter. For instance, the importance of nuclear inflows induced by tidal torques from one galaxy onto the other (Keel et al., 1985; Barnes & Hernquist, 1991) increases with decreasing pericenter distance. Similarly, shocks have a strong impact during close passages, when the dense regions of the ISM overlap (Jog & Solomon, 1992). Conversely, tidal and turbulent compression found in simulations (Renaud et al., 2014) are of gravitational origin and are thus remote effects, occurring even during fly-bys and spanning large volumes.

Understanding these phenomena, their interplay, and their possible observable signatures is thus an important key to pin down the physics of the physics of star formation, and to determine which quantities should be probed in future observations.

To address this problem with simulations, parsec-resolution is necessary to capture the cold, dense ($\sim 10\mathrm{\leavevmode\nobreak\ K}$, $\gtrsim 10^{3}\mathrm{\leavevmode\nobreak\ cm^{-3}}$) phase of the ISM, and the details of its turbulent nature. The high cost of such simulations caps the number of cases that can be studied. Here we use an Antennae-like simulation which encompasses a wide diversity of star-forming environments (nuclei, shocks, extended formation). We explore these spatial variations, attempt to connect them to underlying physical process(es), and seek possible differences of their signatures in the CO emission that could be observed in real galaxies. Such differences could be particularly insightful when numbers of observational diagnostics cannot be resolved and only integrated fluxes are available.

2 Method

2.1 Numerical method

The numerical method followed is identical to that of Renaud et al. (2019), and summarized here. The simulation considered is an hydrodynamical model of an Antennae-like major merger (Renaud et al., 2015a). The simulation has a resolution of $1.5\mathrm{\leavevmode\nobreak\ pc}$ and includes, background UV heating, atomic and molecular cooling, star formation (at an efficiency of 2% per free-fall time above $50\mathrm{\leavevmode\nobreak\ cm^{-3}}$, i.e. where the molecular gas fraction is expected to be close to unity at solar metallicity, Krumholz & Gnedin 2011), photo-ionization, radiative pressure and type-II supernova feedbacks.

Our simulation only comprises a limited number of feedback effects. In particular, AGN feedback is not included. Although this could affect our conclusions, in particular in the nuclear regions of the galaxies where X-ray photons could heat up the gas and excite CO (e.g. Krips et al., 2011), the presence of an AGN in the Antennae has not been unambiguously established (see Brandl et al. 2009, and also Ueda et al. 2012 which mention the possibility of an hidden AGN causing high line ratios in the nucleus of NGC 4039).

The results related to star formation depend on the sub-grid recipe used and details would change when using different methods. The implementation and parameters we use have been designed, at our working resolution, to reproduce observations of integrated quantities like the SFR and the depletion time of the progenitor galaxies run in isolation, typical of main sequence galaxies. The simulation also retrieves the observed star formation activity of the Antennae during the interaction (see below). Uncertainties remain on the validity of non-integrated quantities like the number, spatial distribution and masses of star-forming clouds that would change, e.g., by varying the density threshold and SFE, even if keeping the SFR unchanged. However, the morphology of the star-forming areas and the properties of the young massive cluster formed in the simulation are in line with the observations (Renaud et al., 2015a), such that our models can be used to help the interpretation of observation data down to the scale of $\sim 10\mathrm{\leavevmode\nobreak\ pc}$.

2.2 CO emission

The large velocity gradient (LVG) method is applied in post-processing to compute the intensity of the CO (J$\to$J-1) emission lines based on the local density and temperature of the ISM and lookup data (Weiß et al., 2005), and velocity gradients along the line of sight to infer absorption (see Bournaud et al., 2015, for details). The analysis is conducted at $1.5\mathrm{\leavevmode\nobreak\ pc}$ resolution for almost all the lines-of-sight, with rare exceptions at $3\mathrm{\leavevmode\nobreak\ pc}$ where the highest densities are not met. Data presented here are velocity-integrated Rayleigh-Jeans brightness temperatures, in units of $\mathrm{\leavevmode\nobreak\ K\ km\ s^{-1}}$. For convenience, we also provide the conversion to flux integral in units of $\mathrm{\leavevmode\nobreak\ Jy\ km\ s^{-1}}$. The CO-to-H2 conversion factor $\alpha_{\mathrm{CO}}$ (expressed in $\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}\mathrm{\leavevmode\nobreak\ K}^{-1}\mathrm{\leavevmode\nobreak\ km}^{-1}\mathrm{\leavevmode\nobreak\ s}\mathrm{\leavevmode\nobreak\ pc}^{-2}$ in the rest of the paper) is inferred using the mass of the gas denser than $50\mathrm{\leavevmode\nobreak\ cm^{-3}}$ (i.e. matching the star formation density threshold), which leads to a good agreement with empirical values for a Milky Way-like galaxy (Renaud et al., 2019). The dispersions presented are the root mean square values of a sample of 5 lines-of-sights, obtained by tilting the system by $\pm 15^{\circ}$ along the axes of the plane of the sky.

2.3 Regions of interest

After exploring the time evolution of the system in Renaud et al. (2014, 2019), we now focus on the snapshot representing the best morphological match to observational data. This instant is found $5\mathrm{\leavevmode\nobreak\ Myr}$ after the second pericenter passage, and $19\mathrm{\leavevmode\nobreak\ Myr}$ before the third one which marks the onset of final coalescence. At this time, the SFR of the entire system reaches $15.7\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}$ (i.e. in the range of the observed values, $7\textrm{\,--\,}20\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}$, see Zhang et al. 2001; Brandl et al. 2009; Klaas et al. 2010) and the global $\alpha_{\mathrm{CO}}$ is 2.9 (which puts the system in the post-burst regime, see Renaud et al. 2019).

For the sake of simplicity, we adopt a line-of-sight and an orientation comparable to that of the real system. In this snapshot, we select four star-forming regions of $1\mathrm{\leavevmode\nobreak\ kpc}\times 1\mathrm{\leavevmode\nobreak\ kpc}$, as shown in Figure 1. The selection comprises the two galactic nuclei, the overlap region (on the Eastern side, left in Figure 1, where the disk remnants mix) and the arc on the opposite side of the overlap in the Northern galaxy. In the real Antennae, the corresponding regions yield strong CO(1-0) emission (Wilson et al., 2003), and host or have recently hosted intense star formation activities, in particular in the form of massive star clusters ($\approx 5\textrm{\,--\,}10\mathrm{\leavevmode\nobreak\ Myr}$ old, see e.g. Whitmore et al., 1999; Bastian et al., 2009; Herrera et al., 2012). Our previous explorations of the system as a whole suggest that a variety of physical processes are responsible for this activity (Renaud et al., 2014, 2015a). We seek here signatures of these processes in observable properties of the ISM, and in the CO emission.

For comparison, we also use another snapshot of the same simulation, during the starburst episode triggered by the coalescence, when the SFR peaks at $97.2\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}$. We refer to this point as “ULIRG” in the following. We also include the Milky Way-like case ($1\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}$) and the gas-rich clumpy galaxy ($45\mathrm{\leavevmode\nobreak\ M}_{\odot}\mathrm{\leavevmode\nobreak\ yr}^{-1}$, representing a disk at redshift $\sim 2$) presented in Bournaud et al. (2015).

3 Results

3.1 The triggers of star formation and their delays

We start by seeking spatial (anti-)correlations between the star formation activity and processes known to trigger (quench) it. The top-left panel of Figure 2 overlays the gas density map of our Antennae model with the depletion times measured at the scale of $125\mathrm{\leavevmode\nobreak\ pc}$, using the stars younger than $10\mathrm{\leavevmode\nobreak\ Myr}$.222Therefore, only the regions having hosted star formation in the last $10\mathrm{\leavevmode\nobreak\ Myr}$ appear in this map. A larger age cut would provide better statistics, but would also encompass stars that might have significantly drifted away from the locations (and thus the physical conditions) of their formation.

We first find that, across the galaxies, $t_{\mathrm{dep}}$ has a scatter of several orders of magnitude, with sharp gradients especially around the overlap and the galactic nuclei (up to $\approx 20\mathrm{\leavevmode\nobreak\ Myr}\mathrm{\leavevmode\nobreak\ pc}^{-1}$), as discussed in Section 3.5. Other regions also host a range of depletion times. The structure of the Northern arc (including our right-most region of interest and its gaseous extension toward the top-left) yields a broad distribution of depletion times, as well as dense gas clumps not hosting star formation. The map of shear (top-right panel of the figure) indicates that such clouds found at the edge of the disk remnant lie in areas of strong shear (e.g. at $x\approx 0.9\mathrm{\leavevmode\nobreak\ kpc}$, $y\approx 2.2\mathrm{\leavevmode\nobreak\ kpc}$ on Figure 2), likely caused by the differential tidal stretch that truncates the disks. The shear is thus destroying the over-densities and these structures will not survive long-enough to form stars. Conversely, other clouds (further away from the tidal truncation of the disk) are on the verge of hosting star formation in the next few Myr (e.g. at $x\approx-1.9\mathrm{\leavevmode\nobreak\ kpc}$, $y\approx 2.3\mathrm{\leavevmode\nobreak\ kpc}$).

This indicates that the large scatter of $t_{\mathrm{dep}}$ is partly due to gas structures being at different evolutionary stages, toward collapse or dissolution, and that the influence of the large-scale environment is blended into this “intrinsic” time sequence (see Grisdale et al. 2019 on the diversity of SFEs at cloud scale). The variations of depletion times are further discussed in Section 3.3, in relation with the distribution of dense gas.

The four panels of Figure 2 illustrate the complexity of the picture by displaying the maps of some factors triggering or quenching star formation. In these maps, we retrieved expected general trends in most of the regions: volumes with strong shear have long depletion times, and those with strong compressive turbulence (i.e. high Mach number and high compressive to solenoidal ratio) tend to favor efficient star formation. Notable exceptions are the vicinity of the galactic centers, discussed in Section 3.5. However, a more detailed and quantitative inspection of the relations between the plotted quantities reveals no clearer trends than the general qualitative ideas mentioned above. In that sense, the causal link between these processes and star formation is not un-ambiguously established here, although it appears when considering the time evolution of the system (see Renaud et al. 2014, their figure 1).

The reason for this is likely the time delay between the onset of a physical process and the measure of the resulting star formation activity. In our case, using stars younger than $10\mathrm{\leavevmode\nobreak\ Myr}$ implies that the delay is at least $10\mathrm{\leavevmode\nobreak\ Myr}$. In other words, the process that triggers (quenches) star formation, as early as the formation of a cloud, can be long gone when the depletion time is measured to be dropping (rising). This is particularly important in fast evolving systems like galaxy mergers where dynamical times are short. If causal relations are not found in our Eulerian scheme, they would likely get clearer by using Lagrangian methods (e.g. with tracers in grid codes, see Genel et al. 2013; Cadiou et al. 2019) and by tracking the gas material throughout the phases of compression, cooling, fragmentation and star formation (Semenov et al., 2017).

Apparent in these maps is the diversity of environments leading to a broad dynamical range in the quantities plotted. For instance, the vicinity of the nuclei host strong shear and the depth of the potential implies fast moving gas structures, while the arc in the outskirts of the disk evolves with longer timescales. As a result, the relevant timescales for cloud evolution (and thus star formation) are strong functions of position and time, such that it is impossible to define global timescales. Therefore, observational tracers of any stage of the star formation process (from dense gas to gas-free star cluster) suffer from this diversity of timescales by likely introducing biases in the interpretation of galaxy-wide maps. For the same reason, the diversity of timescales translates into a diversity of scale-lengths which also introduce biases, as discussed in Section 3.4.

According to our model, the on-going collision occurs as the two galaxies are falling back onto each other after $\sim 150\mathrm{\leavevmode\nobreak\ Myr}$ of separation since their first encounter. The Northern galaxy (NGC 4038) is moving southward toward the other progenitor. The overlap of the disks starts on the East side (left on Figure 1) and moves roughly south-westward as the two galaxies fly through each other. The time taken by the overlap to move from one side of the stellar disks to the other is about $10\mathrm{\leavevmode\nobreak\ Myr}$, i.e. comparable to the timescales of the physical processes mentioned above, to the depletion times of the clouds, and to the propagation of feedback. If starbursts were exclusively triggered by shocks (when the gaseous material of the two galaxies overlap, either through cloud-cloud collisions, or between the gas reservoirs of the progenitors, see Jog & Solomon 1992 and Section 4.2), one would expect to find a signature of this evolution with star formation running from East to West of the system. In that case, the galaxies would exhibit gradients of tracers of star formation from the signatures of the latest stages of star formation (e.g. $\mathrm{H}\alpha$, UV, IR) toward the East, behind the shocks themselves, to more early-stage indicators (e.g. dense gas, CO, HCN, HNC, etc.) to the West. However, MUSE observations covering the entire central region reveal no such gradients in the positions or sizes of the Hii regions (Weilbacher et al., 2018), indicating that shocks are not the only triggers of the starburst activity. The observed complexity is in line with our results and we predict that similar conclusions could be reached with other tracers.

In this context, using velocity dispersion only to trace star formation might lead to errors, since an above-average dispersion could be both the cause of star formation (e.g. compressive turbulence, shocks) and its consequence (feedback). Some strong shocks not associated with young stars have been detected in the Antennae (Haas et al., 2005), but exclusively in the overlap, very likely marking the onset of localized shock-triggered star formation.

Due to this diversity in star-forming regions, their triggers and their delays, it is involved to establish clear global relations. To circumvent this, we examine now the properties of the dense ISM to identify features specific to interacting systems.

3.2 Excess of dense gas in the star-forming regions

Figure 3 shows the gas density probability distribution functions (PDFs) in our four regions and in the entire system. All PDFs yield a secondary peak of dense gas ($\gtrsim 10^{4}\mathrm{\leavevmode\nobreak\ cm^{-3}}$), but this feature is most pronounced in the nuclei and, to a lower extent in the overlap. Together, the four kpc-scale regions gather 87% of this dense gas across the galactic system, and only a few other clouds account for the rest (on the Eastern side of the Northern arc, and at the South-most end of the southern galaxy, recall Figure 1). We note that the average density of the secondary peak varies from one region to the next. Although gas less dense than this secondary peak can also form stars, the density dependence of the star formation law ($\propto\rho^{3/2}$) makes the dense peak the main contributor ($85\%$) to the total star formation activity.

As shown in Figure 3 (bottom panel), the gas clouds of the nuclei themselves (i.e. the innermost $\sim 100\mathrm{\leavevmode\nobreak\ pc}$) constitute a significant fraction of the mass in the secondary peak of their regions, but are not the sole contributors: the peak is still present in the PDF of the kpc-scale regions without the nuclei and its maximum is found at the same density with and without the nuclei. Star formation in the central regions is thus not restricted to the sole innermost gas structures but span numbers of other clouds with comparably high densities, as visible in the density maps (Figure 1). These star-forming clouds in the immediate vicinity of the nuclei could result from the fragmentation of infalling material, but the excess of dense gas and the similarities of the PDFs with and without the nuclei suggests a more efficient process.

In addition to the presence of the peak in the overlap and the arc, this indicates that the physical origin of the secondary peak affects large-scale volumes and is not (only) related to nuclear inflows. Furthermore, Renaud et al. (2014, their figure 2) showed that the secondary peaks appear as early as the first passage in the interaction, and are thus not only limited to penetrating encounters, in-line with the global increase of compressive turbulence. The link between the shape of these PDF and turbulence is further discussed in Section 4.1.

It is important to note that these PDFs do not result from a shift nor a widening of those found in isolated galaxies: the presence of a secondary peak is unique to interactions and reflects that star formation in such systems is not a simple extrapolation of that of isolated galaxies. Therefore, the physical reason for the enhanced star formation, both in term of rate and efficiency (i.e. SFR and $t_{\mathrm{dep}}$), lies in the process(es) that generate this excess of dense gas. In the next section, we connect this excess to reduced depletion times.

3.3 Variations of depletion time

At the instant considered, the entire galactic system has a depletion time of $100\mathrm{\leavevmode\nobreak\ Myr}$ ($180\mathrm{\leavevmode\nobreak\ Myr}$, if considering the molecular and atomic gas), i.e. about one order of magnitude shorter than that of normal star-forming disk galaxies, which places it in the regime of LIRGs/ULIRGs defined by Daddi et al. (2010b).

Our kpc-sized regions have depletion times of $65\mathrm{\leavevmode\nobreak\ Myr}$ (resp. 120) in the overlap, $80\mathrm{\leavevmode\nobreak\ Myr}$ (resp. 120) in the nuclei and $170\mathrm{\leavevmode\nobreak\ Myr}$ (resp. 250) in the arc (when accounting for the molecular phase only, and the total gas mass respectively). These measurements are compatible with the estimates of Bigiel et al. (2015) who reported shorter depletion times when using tracers of dense gas ($\sim 40\textrm{\,--\,}160\mathrm{\leavevmode\nobreak\ Myr}$ in total infrared) and longer times with CO, but at slightly larger scales ($\approx 2\textrm{\,--\,}3\mathrm{\leavevmode\nobreak\ kpc}$) and by adopting disk-like luminosity-to-mass conversion factors (while the $\alpha_{\mathrm{CO}}$ model from Renaud et al. 2018 would lead to depletion times about 1.5 shorter). Already by selecting the regions hosting the most intense star formation activity, we report significant variations of $t_{\mathrm{dep}}$, but these 4 regions are not representative of all the star-forming areas, and thus we also explore the variation of $t_{\mathrm{dep}}$ in a more systematic way, and with a larger statistical sample, on a cloud-basis.

The top panel of Figure 4 shows the positions of star-forming clouds in the Schmidt-Kennicutt diagram. The surface densities are computed at the scale of each cloud ($\sim 30\mathrm{\leavevmode\nobreak\ pc}$), independently of each others. The overall relation is super-linear as the points are best-fitted by power-laws of indexes 1.4 when including the 2 centermost clouds, and 2.0 when excluding them, but with large deviations (69.2 and 19.2 respectively, in the units of the figure).

Instead of interpreting this distribution as a unique power-law with a large scatter, we highlight that a line of constant depletion time splits the clouds into two families: those with short depletion times ($t_{\mathrm{dep}}\lesssim 40\mathrm{\leavevmode\nobreak\ Myr}$) are found in the dense-end of the distribution while the others are systematically less dense. The dichotomy directly relates to the PDFs of these clouds, shown in the bottom panel of Figure 4. Clouds with short $t_{\mathrm{dep}}$ systematically yield an excess of dense gas ($\gtrsim 10^{4}\mathrm{\leavevmode\nobreak\ cm^{-3}}$) with respect to the others. More than the mere expected relation between short $t_{\mathrm{dep}}$ and high densities, this shows that all the clouds with $t_{\mathrm{dep}}<40\mathrm{\leavevmode\nobreak\ Myr}$ contain an excess of dense gas, while none of the others do. Furthermore, not all star-forming clouds but only those with a short $t_{\mathrm{dep}}$ contribute to the secondary peak noted in Figure 3, demonstrating that the excess of dense gas noted before as the presence of a secondary peak is associated with short depletion times, in addition to high SFRs.

At small scale, the density of star formation rate depends on the gas density in a super-linear fashion (e.g. $\rho_{\rm SFR}\propto\rho^{1.5}$ in the case of constant efficiency per free-fall time, as adopted here). This implies that, for a given gas mass, an excess of dense gas (our secondary peak) leads to a more efficient star formation (i.e. a decrease in depletion time).

Both the Schmidt-Kennicutt diagram and the distribution of $t_{\mathrm{dep}}$ show that the interacting galaxies can contain simultaneously clouds with unimodal PDFs and long depletion times (in line with those of isolated galaxies), and clouds with an excess of dense gas and significantly shorter depletion times. The mix of both types sets the properties of the galaxies at kpc-scale which, in our case at this instant, puts the Antennae among the galaxies with shortened depletion time, in the starburst regime.

Comparable variations of $t_{\mathrm{dep}}$ have also been found observationally across the Antennae (Bigiel et al., 2015), IC 4687 (Pereira-Santaella et al., 2016) and NGC 2276 (Tomičić et al., 2018). In the Antennae (the observation and our model), the complexity of the geometry due to the deformation from the previous interaction and the partial overlap of two galaxies lead to complex galactic-scale gradients, neither radial, nor from one side of the galaxy to the other (unlike in NGC 2276 where ram pressure and tides imprint clear variations of $t_{\mathrm{dep}}$, Tomičić et al. 2018). The absence of such large-scale, smooth variations indicates that the triggers of star formation are a complex interplay of physical processes, and/or of processes with timescales or the order of the dynamical times of the clouds, as hinted in Figure 2.

The exact physical origin of this dense gas is difficult to assess in our models, as discussed in Section 3.1. However, some additional clues can be obtained by examining the star formation histories of the clouds. Clouds with long depletion times ($t_{\mathrm{dep}}>40\mathrm{\leavevmode\nobreak\ Myr}$) are forming less stars now than $5\mathrm{\leavevmode\nobreak\ Myr}$ ago. This indicates that their formation activity is slowing down, and will probably stop within a few Myr, due to the consumption of dense gas and/or its dispersion by feedback. The picture is less clear for the clouds with shorter $t_{\mathrm{dep}}$. Some show an increasing star formation activity while others behave like the longer $t_{\mathrm{dep}}$ objects mentioned previously. This points toward a mix between evolutionary and environmental effects shaping the star formation activity of these clouds (see Grisdale et al. 2019, in the context of the Milky Way.).

Interestingly, the locii of the clouds in the Schmidt-Kennicutt diagram (Figure 4) is independent of their host region, and more generally of their position in the galaxies. In particular, we find dense clouds with short depletion times in all the regions, which suggests that differences in the kpc-scale environment play a minor role on the properties of the clouds themselves. The other possible explanation is that the diversity of processes leading to short depletion times across the galaxies are indistinguishable in a Schmidt-Kennicutt diagram (with the exceptions of the two nuclei which reach extreme densities, but comparable depletion times as many other clouds).

3.4 Scale dependence

In Figure 5, we consider different areas (or “beam” size), centered on the star formation peaks within our four regions of interest. Due to the spatially extended enhanced star formation activity, any of our regions yields a depletion time an order of magnitude shorter than that of equivalent volumes in a normal star-forming disk ($\sim 1\mathrm{\leavevmode\nobreak\ Gyr}$, as expected). Not surprisingly, increasing the beam size leads to, on average, including a larger fraction of non-star-forming gas, which increases the depletion time (see e.g. Schruba et al., 2010). Details of this relation however depend on the size of individual star-forming structures and the separation between them within the beams, and the latter varies across the galaxies due to the diversity of processes involved (recall Figure 2). This is illustrated by the differences between our four regions of interest in Figure 5, on top of the expected general trend of increase of $t_{\mathrm{dep}}$ with the beam size, as also seen in the isolated galaxy.

For instance, star formation in the arc occurs in a handful of clouds, surrounded by non-star-forming, diffuse material. Therefore, the depletion time in this area rapidly increases when it is measured at scales larger than the size of these clouds. In the other regions, this effect is still visible but without a particular scale identifiable due to the larger filling factor of clouds in the area333The mean separation between clouds is $270\mathrm{\leavevmode\nobreak\ pc}$ in the arc and $120\textrm{\,--\,}140\mathrm{\leavevmode\nobreak\ pc}$ in the other regions. The clouds have different shapes (for instance, some are elongated due to shear), but their volume is approximately constant across the galaxy, corresponding to spheres of radius $30\textrm{\,--\,}35\mathrm{\leavevmode\nobreak\ pc}$., which means that increasing the beam size leads to including more clouds, i.e. a lower fraction of non-star-forming material than in the arc.

Figure 6 shows the distribution of depletion times when measured at cloud (30 pc) and galactic (500 pc) scales. At large scale, the regions of shortest depletion times are blended with less and no star-forming material, which truncates the distribution on the short $t_{\mathrm{dep}}$ side. For the same reason, regions not forming stars ($t_{\mathrm{dep}}=\infty$) are binned with regions with long (but finite) depletion times ($\gtrsim 5\mathrm{\leavevmode\nobreak\ Gyr}$), which extends the distribution to its long side with respect to that at small scale. As expected, using beam sizes larger than star-forming regions dilute them towards longer $t_{\mathrm{dep}}$.

The galactic-scale distribution reveals the two modes of depletion times: the main one at $\sim 100\mathrm{\leavevmode\nobreak\ Myr}$ and a secondary at $\sim 2\mathrm{\leavevmode\nobreak\ Gyr}$. These average $t_{\mathrm{dep}}$’s approximately correspond to that of the disk and starburst regimes found in Schmidt-Kennicutt diagrams (e.g. Daddi et al., 2010b). At cloud scale, however, the distribution is narrower (although it appears wider due to the log-scale of Figure 6) because of the lack of the dilution effect mentioned above. However, it still spans several orders of magnitude in depletion times, despite the uniform star formation recipe at constant SFE imposed (0.02 at $1.5\mathrm{\leavevmode\nobreak\ pc}$). Details on the physical origin of this scatter are discussed in Grisdale et al. (2019) in the context of more quiescent disk galaxies, and Kraljic et al. (in preparation) in the case of interacting systems.

In summary, all regions contain clouds with short and long depletion times, with an excess of dense gas or not. This implies that the galactic-scale dynamics has either no influence on the depletion times of individual clouds and thus that the scatter is a purely intrinsic evolution, or that the different dynamical triggers at stake have the same signatures on the depletion times. However, the details in the scale dependence of $t_{\mathrm{dep}}$ reveal a different organization of the star-forming material from one region to the next (in particular in the arc), linked with the kpc-scale triggering of cloud formation.

3.5 Long depletion times in the nuclei

Figure 5 shows that the depletion time yields a steeper dependence with scale for beam sizes $\lesssim 100\mathrm{\leavevmode\nobreak\ pc}$ in the nuclei than in the overlap and the arc. The reason for this originates from the peculiarity of the nuclei in the gravitational potential of the galaxies. Being natural convergence points of flows and minima in the potential, they induce strong shearing and tidal effects along the fueling channels toward them (recall Figure 2). These effects smooth out most (but not all) of the gas clumps, and in turn prevent their collapse on over-dense seeds. Consequently and despite large amounts of dense gas, star formation is considerably slowed down a few $\sim 100\mathrm{\leavevmode\nobreak\ pc}$ away from the centers. However, when the gas reaches the nuclei themselves ($\sim 50\textrm{\,--\,}100\mathrm{\leavevmode\nobreak\ pc}$ in our simulation), it accumulates, leading to important over-densities and the observed star formation activity. As a result, efficient star formation is concentrated in the main structure of the nuclei, while it is more evenly distributed in the overlap, due to the large shocked and compressed volumes. This situation is qualitatively comparable to that of barred galaxies where the bar(s) generate the torques necessary to fuel large amounts of gas, but where shear increases $t_{\mathrm{dep}}$ of this gas, or even quenches star formation (see Emsellem et al., 2015; Renaud et al., 2015b, in the context of the Milky Way), with the notable difference that the excess of dense gas in the Antennae allows some clouds to survive destruction due to their higher binding energy (recall Figures 1 and 2 and see Section 4.2).

This point explains why observations focusing on the nuclei can report longer depletion times there than in the overlap, despite higher SFRs (see e.g. Bigiel et al. 2015). Our results predict that observing these regions with beams comparable to the sizes of the star-forming area only (i.e. encapsulating as little surrounding non-star-forming material as possible) would lead to the expected shorter $t_{\mathrm{dep}}$ in the nuclei. We also note that the scale dependence of the depletion time varies with the clumpy nature of star formation, which itself depends on the tracer used to measure the SFR. For instance, a tracer of molecular gas is expected to reflect larger objects and thus milder variation of $t_{\mathrm{dep}}$ with the scale than, e.g., the $\mathrm{H}\alpha$ luminosity which is more concentrated around the star-forming sites.

3.6 CO emission

3.6.1 Spectral line energy distributions

Figure 7 shows the spectral line energy distributions (SLEDs), in term of integrated intensity and flux integral, of our regions ($1\mathrm{\leavevmode\nobreak\ kpc}$) compared to entire galaxies ($10\textrm{\,--\,}30\mathrm{\leavevmode\nobreak\ kpc}$): the Antennae at the same time (named “global”), at the time of maximum starburst (named “ULIRG”), a Milky Way-like galaxy, and a high redshift clumpy disk.

Our regions of interest encompass a range of densities, from the peaks in the clouds, to the more diffuse warm medium in between them. For this reason, we expect our modeled CO emission to differ from that observed when exclusively focusing on clouds. Our results would thus better match observations with beams not resolving the individual clouds and thus blending the peaks of cold emission with warmer media, i.e. on scales of the order of and larger than $\sim 1\mathrm{\leavevmode\nobreak\ kpc}$. Indeed, while our results for the Antennae are in qualitative agreement with the observations of Bigiel et al. (2015, $\sim 0.6\mathrm{\leavevmode\nobreak\ kpc}$ beams, centered on peaks), they better match those of Schirm et al. (2014, $\sim 4\mathrm{\leavevmode\nobreak\ kpc}$ beams), as discussed in Appendix A. All our 4 regions and the Antennae system (global and ULIRG in Figure 7) yield a strong CO(3-2) emission, scattered around CO(3-2)/CO(1-0) $\approx 0.6$, which is similar to expected values in (U)LIRGs (Yao et al., 2003; Papadopoulos et al., 2012; Bolatto et al., 2015).

The high CO(3-2) intensity (in particular the ratio of brightness temperatures CO(3-2)/CO(2-1)) in the kpc-scale regions is not retrieved in any of the large scale measurements, regardless of the star formation activity. A similar trend might be seen in the observational SLEDs from Bigiel et al. (2015) who used a beam size smaller than our regions ($\approx 0.6\mathrm{\leavevmode\nobreak\ kpc}$). This suggests that such emission is localized in small volume(s) and thus gets blended with regions of weaker emission at larger scales. This may correspond to the CO(3-2) detection in the real system of a bright super giant molecular cloud (the so-called “firecracker”) supposedly being the future formation site of a massive star cluster (Whitmore et al. 2014, see also Johnson et al. 2015).

At higher excitation levels ($J_{\mathrm{upper}}\geq 4\textrm{\,--\,}5$), the nuclei and the overlap regions have SLEDs comparable to that of the ULIRG phase, but the arc (and the galactic-scale “global” measurement) resembles more the high-redshift clumpy galaxy, with a significant bend in their SLEDs (in integral intensity). We propose an interpretation of these features in Section 4.3.

Figure 8 shows the contribution of the regions (1 kpc scale) to the global CO(J$\to$J-1) intensity (10 kpc scale). As expected, each region contributes more than average to the global intensity, for all transitions. Their combined contribution ranges from $\approx 40\%$ (CO(1-0)) to $75\%$ (CO(6-5)) of the total intensity, and 85% of the SFR, although they only span $4\%$ of the area considered. The remaining, yet significant, fraction of the flux originates from less concentrated regions.

The relative contribution of the individual regions increases with $J_{\mathrm{upper}}$ up to 6 for the nuclei and the overlap. However, that of the arc peaks at CO(4-3), which further tells apart the arc. As shown in Figures 7 and 8, the Southern nucleus is more excited and has a higher CO intensity than the northern one at all $J_{\mathrm{upper}}$, despite hosting lower amounts of dense gas and with lower maximum densities (Figure 3), and a comparable distribution of temperatures. One would then expect a stronger emission from the Northern nucleus, contrarily to our measurements. However, it hosts a stronger shear and higher Mach number (Figure 2) than its southern counterpart, i.e. an overall larger gas velocity dispersion, which accounts for a stronger absorption of the CO flux captured by the LVG method. The exact relative roles of the gas density, the temperature and the velocity dispersion are difficult to assess due to sharp variations of these quantities on small-scales and the limited resolution of our analysis. Our results on the differences between the regions presented in the previous sections highlight that the relative importance of each quantity varies across the galaxy (and thus likely in time too), in a complex manner (see also Section 4.3).

3.6.2 Spatial variations of $\alpha_{\mathrm{CO}}$

Figure 9 shows the CO-to-molecular conversion factor $\alpha_{\mathrm{CO}}$ in our regions. As noted in Renaud et al. (2019), the SFR and its surface density are not good tracers of variations in $\alpha_{\mathrm{CO}}$, and the recent SFH must also be considered, in particular to account for the effects of feedback on temperature and velocity dispersion. When doing so, we find that all regions contain stars formed during the first passage ($\approx 150\mathrm{\leavevmode\nobreak\ Myr}$) and a recent boost of the activity ($\approx 20\mathrm{\leavevmode\nobreak\ Myr}$) due to the second encounter. Tracking backward in time the elements constituting our four selected regions reveals that the arc is made exclusively of material from the outskirts of the disk, while the overlap comprises material from a wider range of radii, including only a couple of kpc from the galactic centers. The difference in gas densities between these volumes explains that, during the separation of the galaxies, star formation slows down significantly more in the future arc material than in the elements ending up in the other regions. This translates into a deficit of $20\textrm{\,--\,}100\mathrm{\leavevmode\nobreak\ Myr}$ old stars in the arc compared to the other regions. Furthermore, tides, shear, shocks, and the resulting turbulence reshape the structures of clouds, particularly in the nuclei and the overlap, which shifts young stars away from the densest gas. As a result, stellar feedback can propagate further and with less resistance to large scales in the nuclei and the overlap (see Section 4.4). Therefore, the ISM of the nuclei and the overlap contain more gas of intermediate densities and warm temperatures ($\sim 10\mathrm{\leavevmode\nobreak\ cm^{-3}}$, $10^{4}\mathrm{\leavevmode\nobreak\ K}$) than the arc. The media in which CO emission propagates are thus different in the arc and the other regions, which translates into a different absorption and in turn, into a different $\alpha_{\mathrm{CO}}$. The same argument explains the slight difference between the two nuclei (recall also Section 3.6.1).

The new onset of star formation triggered by the second passage ($\sim 5\mathrm{\leavevmode\nobreak\ Myr}$ ago) is too recent to have significantly altered the CO emission and propagation over kpc-scales (see Renaud et al. 2019) and therefore only the older activity affects our results. The galaxies thus lie on the post-burst regime and with slightly higher $\alpha_{\mathrm{CO}}$’s in the outskirts of the disks than in the central regions.

As noted in Renaud et al. (2019), $\alpha_{\mathrm{CO}}$ for the entire system runs from 4.27 in the isolated phase to 1.78 at the peak of the starburst triggered by the final coalescence, i.e. a variation by a factor of 2.4 along the course of the merger. This is comparable to the variation between our regions of interest from our selected snapshot (from 1.3 to 2.8, i.e. a factor 2.2, see Figure 9). However, we note that we have focused our analysis on actively star-forming regions, and that $\alpha_{\mathrm{CO}}$ is likely much higher in other areas. Therefore, we expect variations of $\alpha_{\mathrm{CO}}$ across space to be possibly larger than those along time, even when a galaxy goes from the disk regime to the starburst one, which could have important implications on the interpretation of observations resolving the diversity of star-forming environments, even across a given galaxy.

4 Discussion

4.1 Compressive turbulence and the shape of the density PDF

Sub-kpc scale simulations of the ISM are commonly used to study the inner physics of molecular clouds, down to the individual star level (see a review in Hennebelle & Falgarone, 2012). These simulations do not include the galactic scales from which part of the turbulence is generated ($\gtrsim 0.1\textrm{\,--\,}1\mathrm{\leavevmode\nobreak\ kpc}$, see Agertz et al. 2009; Bournaud et al. 2010; Renaud et al. 2013; Falceta-Gonçalves et al. 2015; Grisdale et al. 2017), and they must resort to manually forced turbulence. Such works indicate that the turbulent ISM yields a log-normal density PDF (see e.g. Vazquez-Semadeni, 1994) of which the variance reads:

[TABLE]

with $\mathcal{M}$ the Mach number, and $b$ a dimensionless parameter related to the forcing nature of the turbulence, i.e. the mix between its compressive (i.e. curl-free) and solenoidal (i.e. divergence-free) modes (see e.g. Kritsuk et al., 2007; Federrath et al., 2010). In these simulations, the PDF is found to always have a log-normal shape444A log-normal density PDF may arise artificially from observations due to incompleteness of the lowest column density bins truncated by a limited field of view (Alves et al., 2017). The log-normal shape is however retrieved when considering large enough areas to encompass the low density envelops of clouds in closed contours (Körtgen et al., 2019)., regardless of the turbulence forcing (with a power-law tail at high densities when self-gravity is captured, Ballesteros-Paredes et al. 2011; Elmegreen 2011; Renaud et al. 2013). In other words, for a given Mach number these simulations report a wider, but still log-normal, PDF in compressive turbulence than for equipartition.

As presented in Renaud et al. (2014) and Figure 3, the shape of the density PDF in our interacting galaxy models significantly deviates from the classical log-normal shape during the interactions, with the formation of a secondary peak at high densities (see also Teyssier et al., 2010; Bournaud et al., 2011). The same bimodal feature is also seen in real interacting galaxies, in term of line-width distributions (Sun et al., 2018), further suggesting the superposition of a high $\sigma$ mode to the classical one in efficiently star-forming systems. Renaud et al. (2014) attributed this excess of dense gas to turbulent compression that arises during the interaction: changing the turbulence from energy equipartition between the compressive and solenoidal modes before the interaction to a compressive-dominated turbulence until the post-merger phase increases the gas density contrast and generates this excess of dense gas. Therefore, a discrepancy exists between the log-normal PDF from forced compressive turbulence in small-scale simulations, and the double-peaked one measured in the compressive turbulence of galaxy mergers (as is, i.e. without any manual forcing).

In the context of the Schmidt-Kennicutt relations, the analytical model of Renaud et al. (2012) uses the secondary peak to explain the observed differences in depletion times between disks and starbursting galaxies. In this model (and others of the same nature, see e.g. Padoan & Nordlund 2011; Hennebelle & Chabrier 2011; Federrath & Klessen 2012), using a log-normal PDF shape for all galaxies would lead to a unimodal star formation law, in which the SFR is set by the width of the PDF (Equation 1), and would thus not retrieve the difference between a gas-rich isolated disk with a high SFR (e.g. BzK galaxies, Daddi et al. 2010a) and a local (U)LIRG with the same SFR (Kennicutt, 1998). In other words, increasing the variance of the log-normal PDF (either by increasing $\mathcal{M}$ or $b$ in Equation 1) would only shift a galaxy along the canonical Schmidt-Kennicutt relation, but keeping it on the relation. Such a model could retrieve a diversity of SFRs, but not the diversity of depletion times. Another effect is needed to change the regime, e.g. moving from the relation of disks to that of starbursts by decreasing the depletion time.

As discussed in Section 3.3, the super-linearity of the star formation law linking the local density of SFR to that of gas is the key to change the regime, i.e. to move the galaxy from the relation of disks to that of starbursts. With a wide, unimodal PDF, the densest gas constitutes the tail of the PDF, and increasing its width lead to an higher SFR. Conversely, the presence of a secondary peak in the distribution (i.e. the excess of dense gas we report) and the super-linearity of the star formation law imply that most of the star formation activity takes place at the highest densities, therefore with a higher efficiency than usual (i.e. a shorten depletion time). (See Renaud et al. 2012 for an analytical demonstration, and recall Figure 4.)

We thus report a fundamentally different shape of the PDF between that from manual forcing of compressive turbulence in sub-kpc scale simulations, and that of compressive turbulence consistently arising in our simulations of interacting systems. We do not fully understand the reason of this difference but, assuming it originates from the different setups of the simulations, we postulate it could be related to the injection scale of the turbulence, its forcing mechanism and/or its dissipation. It could also be that galactic simulations contain large gas reservoirs to replenish the PDF when compression occurs, while simulations at cloud scales would suffer from a rarefaction of gas to populate both our two peaks. By nature, galactic-scales simulations like ours are likely better capturing the injection and forcing mechanisms, but the limited resolution implies the dissipation is less accurate than in small-scale studies. Further investigations, with a large dynamical range are required to conclude.

The exact nature of turbulence has important implications on the structure of the ISM and the propagation of feedback, potentially down to the stellar initial mass function (IMF) itself. For instance, theoretical models predict that a compressive turbulence would lead to a bottom-heavy IMF (Chabrier et al., 2014). However, Renaud et al. (2014) have established the existence of compression-dominated turbulence in interactions only at the scale of $\sim 40\mathrm{\leavevmode\nobreak\ pc}$, i.e. several orders of magnitude larger than that of pre-stellar cores. Therefore, there are still many uncertainties on how this type of turbulence, driven by galactic dynamics, cascades and potentially dissipates down to the core scale, in addition to our lack of understanding of the core-to-star mapping (Goodwin et al., 2008; Motte et al., 2018). Small-scale studies with detailed implementations of a range of physical processes conclude on the independence of the IMF with environmental conditions, in particular the turbulence and the metallicity (Bate, 2009, 2019; Myers et al., 2011). However, including self-consistently the full range of galactic-driven effects, in particular turbulence injected at $\gtrsim 100\mathrm{\leavevmode\nobreak\ pc}$ scale, seems essential to capture and predict the possible variations of the IMF. This remains out-of-reach for present-day models and supercomputing resources.

4.2 The galactic nuclei as collision environments

All the aspects of our analysis reveal similarities between the overlap and the nuclei, suggesting that the former could be a scaled-down version of the latter, i.e. with similar physical processes active, but at lower rates. From a morphological point of view, the two types of regions show several dense, cold star-forming clumps, with an average separation of $\approx 130\mathrm{\leavevmode\nobreak\ pc}$, distributed along elongated structures of warm gas at intermediate densities ($\sim 10^{2\textrm{\,--\,}3}\mathrm{\leavevmode\nobreak\ cm^{-3}}$, $\sim 10^{1\textrm{\,--\,}2}\mathrm{\leavevmode\nobreak\ K}$), themselves surrounded by more diffuse medium ($\sim 1\textrm{\,--\,}10\mathrm{\leavevmode\nobreak\ cm^{-3}}$, $\sim 10^{2\textrm{\,--\,}4}\mathrm{\leavevmode\nobreak\ K}$).

The overlap is, by definition, the location of kpc-scale shocks between the gas reservoirs and the clouds (and/or marginally stable over-densities) they contain (Jog & Solomon, 1992). As simulated and observed in other contexts, cloud-cloud collisions trigger star formation (Tan, 2000; Tasker & Tan, 2009; Inoue & Fukui, 2013), in the form of massive clusters (due to the enhanced external pressure Elmegreen & Efremov 1997; Fukui et al. 2014) and possibly massive stars (Motte et al., 2014; Takahira et al., 2014). Although it is difficult in our Eulerian simulation to track individual gas structures and establish collision rates, the presence of thin and dense tidal tails and bridges connecting the density peaks in the overlap confirms the ubiquity of cloud-cloud collisions (Figure 1). It is thus likely that a significant fraction of the star formation activity in the overlap is linked to cloud-cloud or cloud-reservoir collisions.

In the nuclei, the gas is fueled toward the center by negative gravitational torques, mainly from one galaxy on the other (Keel et al., 1985; Hernquist & Katz, 1989; Barnes & Hernquist, 1991), but also from asymmetric structures (see e.g. Emsellem et al., 2015). The clumpy nature of the fueled material (recall Section 3.5) makes the nuclear accretion comparable to a series of cloud-cloud collisions. In that sense, the nuclei share a common mechanism with clouds in the overlap region, but likely at a higher rate due to their preferential location in the potential.

In both the overlap and the nuclei, the trigger(s) of star formation is a short-range effect, of which intensity peaks during penetrating galaxy encounters (but remains relatively high during separation phases due to the high gas density in these inner galactic regions), which then explains that the $\alpha_{\mathrm{CO}}$’s of these areas are closer to the starburst regime highlighted in Renaud et al. (2019), than the transitional post-burst regime as it is clearer for the arc.

In short, the nuclear inflows induced by gravitational torques during close encounters play a similar role as the overlap of the outer disks, by triggering collisions between dense gas structures that lead to a starburst activity lasting $\sim 100\mathrm{\leavevmode\nobreak\ Myr}$ (Renaud et al., 2014). The similar physical processes translate into SLEDs and $\alpha_{\mathrm{CO}}$ factors in the galactic nuclei and the overlap resembling those of starbursting systems, while the activity in the arc appears to have fundamentally different signatures, closer to that of disks with a high but steady SFR.

4.3 The imprint of the ISM structure on the CO emission

By fitting separately the SLEDs from the cold and warm phases, Schirm et al. (2014) established a transition at $J_{\mathrm{upper}}\approx 4\textrm{\,--\,}5$, with the emission from the cold phase dominating the small $J_{\mathrm{upper}}$ regime. All our regions share a comparable distribution of cold ISM ($\lesssim 100\mathrm{\leavevmode\nobreak\ K}$) and have similar SLEDs at small $J_{\mathrm{upper}}$ (Figure 7), which confirms the conclusion of Schirm et al. (2014).

In this regime, the behavior of our SLEDs corresponds to a ULIRG (or starburst) type of activity (as opposed to a disk, steady regime of star formation). However, at higher $J_{\mathrm{upper}}$, the arc deviates from the common ULIRG-like behavior and resembles more the high redshift disk. This is due to a different inter-cloud medium in the arc where the dense gas is highly concentrated into clouds surrounded by relatively diffuse medium ($\sim 10\mathrm{\leavevmode\nobreak\ cm^{-3}}$, recall Figure 3), like the massive clumps in high redshift gas rich disks (see a discussion in Bournaud et al., 2015). Conversely, due to destroying tides, shear and the overall shorter separation between the star-forming regions (Section 3.4), the inter-cloud media of the nuclei and the overlap gather larger amounts of dense and warm gas ($\sim 10^{2\textrm{\,--\,}3}\mathrm{\leavevmode\nobreak\ cm^{-3}}$, $10^{4}\mathrm{\leavevmode\nobreak\ K}$), likely heated by stellar feedback which propagates further from the young stars in such inhomogeneous, porous media (see Ohlin et al. 2019, and the next section). Because these types of regions typically dominate the emission of (U)LIRGs (Figure 8), they are good representatives of our large-scales SLEDs (“global” and “ULIRG” in Figure 7).

The state of the medium surrounding the density peaks is critical to set the external pressure on the clouds, and possibly sets the mass of the star clusters formed (Elmegreen & Efremov, 1997; Ashman & Zepf, 2001; Maji et al., 2017). Although our model predict comparable pressures for the overlap and the nuclei (as estimated observationally by Schirm et al. 2014), a comparison on the mass of the clusters is likely inadequate when considering the extreme cases of nuclear clusters. However, our analysis also leads (at first order) to a higher pressure in the overlap than in the arc, that would translate into more massive clusters in the former region. Yet, comparable massive clusters are found in both regions in our simulation (Renaud et al., 2015a), and in observations (e.g. Bastian et al., 2009), which questions the sole role of pressure in setting the mass of young clusters.

4.4 On the importance of capturing the structure of the cold ISM

Our results on the CO emission and the diversity of depletion times reflect the distribution of dense gas but also the nature of the surrounding, more diffuse ISM. The hierarchical nature of the ISM (Elmegreen & Efremov, 1997; Elmegreen, 2008) is here key to assess how star-forming seeds are fueled with gas, and how fluxes are emitted and absorbed. For instance, we have noted above that the decoupling of the stars from their gas nursery (mainly because of the dissipative nature of the ISM, Renaud et al. 2013, or because of deformation of the clouds by tides or shear) can lead to significant drifts of young stars away from the densest regions of the clouds. In this context, the porosity of the ISM is equally important, particularly to regulate the propagation of feedback effects. This porosity is not only set by the level of turbulence (i.e. the Mach number), but also by the nature of turbulence (compressive vs. solenoidal). While a solenoidal-dominated turbulence tends to homogenize the ISM by smoothing out density contrasts (and thus slowing down or even quenching star formation), a compressive-dominated turbulence favors the formation of high density filaments and clumps, surrounded by low density cavities and chimneys. Such low density volumes then facilitate the propagation of feedback effects (ionization, energy, momentum, chemical enrichment) to large distances, but along preferred directions of least-resistance, as demonstrated by Ohlin et al. (2019). It is thus likely that the regulation of star formation by feedback effects is also a function of the nature of turbulence, and thus, following the results presented here, the depletion time of the cloud itself.

To conclude, the role of turbulence is multiple: by setting the distribution of dense gas (at $\sim 10-100\mathrm{\leavevmode\nobreak\ pc}$ scale), it shapes the formation of the star-forming clouds (which can then fragment and gravitationally collapse Vázquez-Semadeni et al. 2009), it makes the low-density channels in which the propagation of feedback is facilitated, and it sets the optical depth through a series of high and low density structures. It is therefore of prime importance for models and simulations to capture the turbulence and its potential diversity over a large dynamical range.

5 Conclusions

We present an analysis of the spatial variations of the star formation activity in a numerical model of Antennae-like interacting galaxies at parsec resolution. We identify several mechanisms enhancing or slowing down star formation, and study their interplay and observable signatures in different regions of the galaxies. Our main conclusions are as follows.

• •

The interplay of the physical processes acting on star formation exhibits complex variations in space. Their timescales are comparable to that of the dynamics of interacting galaxies (i.e. a few $10\mathrm{\leavevmode\nobreak\ Myr}$), implying that the initial cause(s) for the formation (or destruction) of clouds can be long gone when star formation itself occurs. As a result, no clear (anti-)correlation can be found between causes (inflows, shocks, tidal and turbulent compression, shear, tidal stress) and consequences of star formation. Such correlations are, however, unambiguous when integrating over galactic scales and examining the time-evolution along the interaction and the merger (Renaud et al., 2014).

• •

All star-forming clouds yield enhanced SFRs and reduced depletion times with respect to isolated galaxies. At cloud-scale, a correlation exists between short depletion times ($\lesssim 40\mathrm{\leavevmode\nobreak\ Myr}$) and excess of dense gas visible as a secondary peak in the density PDF ($\gtrsim 10^{4}\mathrm{\leavevmode\nobreak\ cm^{-3}}$). This occurs in all regions of the galaxies, with no apparent dependence on space (i.e. no gradient or symmetries). This translates at galactic scale into the entire system moving from the canonical Schmidt-Kennicutt relation of star-forming disks to the regime of starbursts, i.e. with a depletion times reduced by a factor $\approx 10$.

• •

While all clouds have comparable sizes (radii of $\sim 30\mathrm{\leavevmode\nobreak\ pc}$ at their $50\mathrm{\leavevmode\nobreak\ cm^{-3}}$ contour), the diversity of processes acting on the assembly and destruction of gaseous overdensities leads to different separations between the clouds in different regions of the galaxies ($120\textrm{\,--\,}270\mathrm{\leavevmode\nobreak\ pc}$). This has an imprint on the scale-dependence of measurements of the star formation activity like the depletion time.

• •

The differences in the kpc-scale ISM between the arc and the other regions is visible in the CO emission, at high $J$ transitions. This is likely due to an excess of warm gas in the other regions, originating from the efficient propagation of stellar feedback in the inter-cloud medium ($\sim 10\mathrm{\leavevmode\nobreak\ cm^{-3}}$). It follows that the arc shows CO SLEDs comparable to that of isolated disks (at low and high redshifts). The other regions contribute in larger proportions to the global emission of the system and thus imprint the CO emission of (U)LRIG-like galaxies.

• •

We report significant spatial variations of the $\alpha_{\mathrm{CO}}$ factor across the galaxies, at least as important as the time variations along the course of the interaction i.e. between isolated progenitor phase and the peak of starburst.

• •

The excess of dense gas correlates with the onset of a compression-dominated turbulence, but the shape of the PDF fundamentally differs from that measured in simulations of isolated clouds. The reasons for this discrepancy are still under investigation.

Although the links between the physical processes and their observational signatures are highly degenerated, combining different tracers and indicators of the amounts and distributions of dense gas (CO, HCN, HNC, etc.), its excitation (velocity dispersion, SLEDs including both low and high transitions), and the star formation activity (its rate and efficiency) would help disentangling comparable physical conditions. This could be particularly useful when individual star-forming regions cannot be resolved, e.g. in high redshift galaxies. For instance, at a given (surface) density of gas, a short depletion time traces an excess of dense gas, which calls for an kpc-scale trigger. In such a medium, an highly excited CO at $J_{\mathrm{upper}}\gtrsim 4\textrm{\,--\,}5$ likely corresponds to shocks (or equivalently to an infall toward a unique convergence point). However, a more moderate SLED would indicate compression from larger scales, possibly of tidal and turbulent origin (like in the arc in our study). How these different states are combined and mix across the galaxy could then help inferring the underlying nature of the galaxy and its evolutionary stage in case of an interacting system. It would also provide constraints on the value of $\alpha_{\mathrm{CO}}$ to be adopted (see also Renaud et al., 2019). Note however that such a diversity of mechanisms comes with a diversity of timescales, which can be of the same order as the dynamical time of galaxies, especially in fast evolving systems (e.g. interactions, bars, galactic centers). This blurs the identification of the physical processes at stake in a complex manner, which further calls for the use of as many independent tracers as available.

The identification of the physical trigger of the star formation enhancement in interacting galaxies is a difficult task, even in simulations. Our analysis reveals resemblances between the overlap and the two nuclei, with a comparable influence of shocks between the gas reservoirs and cloud-cloud collisions. This is balanced by a strong shear on the infalling material at the vicinity of the nuclei. In this respect, the arc stands out, as it hosts no cloud-cloud collisions, no shocks and no inflows. The only physical process left to explain the enhancement of star formation (high SFR and short depletion time) is that of tidal and turbulence compression (Renaud et al., 2014). We note that clouds formed there are less numerous than in other regions, but once formed they have comparable properties (sizes, mass, PDF, depletion times) as in the most active regions. The similarities at $1\textrm{\,--\,}30\mathrm{\leavevmode\nobreak\ pc}$ scales is also retrieved in the formation of massive star clusters (Renaud et al., 2015a). The main difference we found is for the CO emission from the warm phase of the ISM ($\sim 10^{4}\mathrm{\leavevmode\nobreak\ K}$, $J_{\mathrm{upper}}\gtrsim 4\textrm{\,--\,}5$). The key in determining observationally the physical trigger of enhanced star formation could thus be hidden in the medium surrounding the clouds themselves.

Acknowledgements.

FR and OA acknowledge support from the Knut and Alice Wallenberg Foundation. OA acknowledges support from the Swedish Research Council (grant 2014-5791). This work was supported by GENCI (allocations A0030402192 and A0050402192) and PRACE (allocation pr86di) resources.

Appendix A Comparison of CO SLEDs with observations

Figure 10 compares the simulated SLEDs with that from observations in the two nuclei and the overlap. The data from Schirm et al. (2014) is obtained for a beam size of $43^{\prime\prime}$, while Bigiel et al. (2015) focuses on the emission peaks at a resolution of $6.5^{\prime\prime}$. Assuming a distance to the Antennae of $19\textrm{\,--\,}22\mathrm{\leavevmode\nobreak\ Mpc}$ (Whitmore et al. 1999; Schweizer et al. 2008), these correspond to $4.0\textrm{\,--\,}4.5\mathrm{\leavevmode\nobreak\ kpc}$ and $0.6\textrm{\,--\,}0.7\mathrm{\leavevmode\nobreak\ kpc}$ respectively.

We report a better agreement between our simulated SLEDs and the values of Schirm et al. than between the two observations, which likely originates from differences in beam sizes. While Bigiel et al. focused on the peaks of emission, the pointings in Schirm et al. cover a much larger area and thus encompass a larger fraction of gas with weaker emission. This is qualitatively closer to the approach followed in our simulations, where our 1 kpc regions gather a wide range of CO-emitting media.

When considering de-normalised values (not shown), our model yields intensities $\sim 11\textrm{\,--\,}18$ times higher than that of Schirm et al. (for $J_{\mathrm{upper}}\leq 4$, see below), which is compatible with a scale factor of $\pi(4\mathrm{\leavevmode\nobreak\ kpc}/2)^{2}/(1\mathrm{\leavevmode\nobreak\ kpc})^{2}\approx 13$ between our surface areas.

We note that our model yields stronger excitations at $J_{\mathrm{upper}}\gtrsim 5$, compared to Schirm et al. Following the arguments presented in Sections 3.6.1 and 4.3, we interpret this difference as a dilution of the contribution from the warm gas ($\sim 1\textrm{\,--\,}10\mathrm{\leavevmode\nobreak\ cm^{-3}}$, $10^{4}\mathrm{\leavevmode\nobreak\ K}$) in the larger beams of Schirm, compared to our regions.