© ESO 2019

1. Introduction

The physical properties of the dense phase of the interstellar medium (ISM) are often estimated using various tracers (CO, HCH, HNC, HCO⁺, see e.g. Gao & Solomon 2004; Talbi et al. 1996; Graciá-Carpio et al. 2006) that each probe different density and temperature ranges. However, to convert the observed luminosities into dynamical quantities such as the mass of molecular gas, we must rely on conversion factors like the ratio of molecular gas mass to CO luminosity (α_co, expressed in M_⊙K⁻¹km⁻¹ s pc⁻² in the rest of the paper), which have been shown to vary significantly with galactic environment (Bolatto et al. 2013). Uncertainties are particularly important in interacting galaxies (see, e.g., Graciá-Carpio et al. 2008). Zhu et al. (2003) empirically found CO-to-dense gas conversion factors that were lower than average in actively star-forming regions (see also Graciá-Carpio et al. 2008; Sliwa et al. 2012; Genzel et al. 2015), while they remain similar to the Milky Way value in high-redshift disks (Daddi et al. 2010a). The physical reason(s) for these variations are not fully identified, which further prevents the derivation of models. In particular, it is still poorly known whether all or only the most extreme major mergers have low α_co.

Numerical simulations and (semi-)analytical models have been used to derive the CO emission and provide estimates of the variations of α_co under different physical conditions. Cloud-scale studies have highlighted the importance of metallicity, dust, turbulence, and the local star formation rate (SFR) in changing α_co (e.g., Glover & Mac Low 2011; Shetty et al. 2011; Clark & Glover 2015; Seifried et al. 2017). In parallel, galaxy-scale works focused on the role of kiloparsec-scale (hydro)dynamics in shaping the dense gas regions and the corresponding α_co (e.g., Narayanan et al. 2011, 2012; Feldmann et al. 2012; Narayanan & Krumholz 2014; Bournaud et al. 2015; Vollmer et al. 2017; Gong et al. 2018; Kamenetzky et al. 2018). Such studies typically conduct statistical analyses over a sample of galaxies, but must compromise on the resolution, typically 50–100 pc, which means that they do not resolve the cold and dense phase of the ISM. This then calls for sub-grid recipes (either live or in post-processing) to describe the molecular gas and its turbulence, calibrated on Milky Way-like galaxies. However, the temperature and the turbulence are both critical in setting the intrinsic CO emission and the optical depth of the surrounding medium. These models thus might not describe the structure and properties of the star-forming sites (≲10 pc) accurately enough to capture the variations in α_co in a diversity of galaxies.

Using simulations of galaxies with high SFRs (≈50 M_⊙yr⁻¹, local mergers in their starburst phase and high-redshift galaxies), Bournaud et al. (2015) reported a regime of low α_co (≈2) in starbursting mergers and significantly higher (≈4) in disks. Here, we examine the time evolution of an interacting and merging galactic pair that experiences starburst phases throughout the entire course of its interaction, and explore the transitions between these regimes using isolated galaxies as references. We conduct our study at parsec-scale resolution, thus explicitly resolving the high densities and cold temperatures (∼10⁶ cm⁻³10 K) of molecular clouds and their inner turbulent structure, without having to resort to sub-grid post-processing techniques. However, this comes at the price of a small sample size, rather than a statistical analysis. In this paper, we identify three regimes of α_co along the course of the interaction, and provide fits for these regimes as functions of observable quantities.

2. Method

The analysis was performed on a hydrodynamical simulation of the Antennae merger presented Renaud et al. (2015) at maximum resolution, that is, 1.5 pc in the densest regions. Our simulation (and other comparable works, Karl et al. 2010; Teyssier et al. 2010) predicts two pericenter passages (with a separation of the progenitor galaxies in between), and a third encounter leading to final coalescence. The simulation treats the entire ISM at solar metallicity and includes heating from UV background, atomic and molecular cooling, star formation, and stellar feedback in the form of photoionization, radiative pressure, and type-II supernovae. Supermassive black holes and active galactic nucleus feedback are not implemented.

Fig. 1. Top panel: evolution of the SFR measured over the last 5, 20, and 50 Myr. Vertical lines indicate the pericenter passages. Symbols mark the instants selected for analysis, with their color-coding time. Their shapes and the shaded areas correspond to the regimes identified below: white, dark gray, and light gray for the disk, starburst, and post-burst regimes, respectively. The black square indicates the instant of best morphological match with the observed Antennae merger. The gray line is the SFR of the progenitor galaxies, run in isolation (i.e., an almost constant SFR of ≈1.5 M⊙yr−1 for each galaxy). Bottom panel: evolution of αco.

Fig. 2. αco as a function of the SFR (measured over 5 Myr) in the merger and in simulations of other galaxies. (The system starts in the top left corner, and at first order, moves clockwise). The symbols are the same as in Fig. 1. The dashed lines indicate linear fits (αco = aSFR + b) to three regimes: disk galaxies, starbursts, and post-burst phases (see text). The typical dispersions of αco (estimated by varying the line of sight, see text) is shown in the bottom left corner.

We performed a large velocity gradient (LVG) analysis to model CO emission, as presented in Bournaud et al. (2015). In short, the intrinsic emission was estimated based on the gas density and the temperature in each cell of the simulation and by referring to lookup emission data (Weiß et al. 2005). Then, if the velocity difference between the source and a gas element along the line of sight was smaller than the intrinsic width of the emitted line, the flux is absorbed. Conversely, if the velocity shift was large enough, the emission flux is unaffected. This process was repeated for all the grid cells over all columns along the line of sight to compute the total CO flux. To compute the CO luminosity to molecular gas mass conversion factor α_co, we considered the mass of the gas denser than 50 cm⁻³. (Adopting a threshold of 10 cm⁻³ would induce variations of 10–20% in α_co, with no systematic trend. This is on the order of the dispersions of our α_co values, as explained below, and does not affect our conclusions.) This criterion gives a mass of molecular gas of 0.2 × 10⁹ M_⊙ in each of the isolated disks, and a decline from 1.9 × 10⁹ to 1.3 × 10⁹ M_⊙ for the pair along the interaction due to the enhanced fragmentation of clouds (see also Fig. 2 of Renaud et al. 2014). Our method gives α_co = 4.6 for a Milky Way-like galaxy (Bournaud et al. 2015), that is, a value well within the uncertainties of that usually adopted observationally: α_co = 4.3 (±0.1 dex, Bolatto et al. 2013). To estimate the dispersion of the measured quantities, we tilted the system by ±15° along two axes, and repeated the LVG analysis. The dispersion was evaluated as the root mean square of these five lines of sight.

3. Results and interpretation

Figure 1 shows the evolution of the SFR and α_co throughout the interaction. The SFR is measured over 5, 20, and 50 Myr, corresponding to the instantaneous SFR and that from observational tracers of different timescales such as H_α and far-infrared. Figure 2 shows the relation between α_co and the instantaneous SFR. Additional points correspond to the progenitor galaxies run in isolation, and to the median values of α_co from simulations of a Milky Way-like galaxy and a z ∼ 2 gas-rich clumpy disk (both in isolation, see Bournaud et al. 2015). In this diagram, we identify three regimes:

• Disk regime: high α_co, regardless of the SFR. This corresponds to the isolated disks (the Milky Way, the high-redshift disk, and our Antennae progenitors), and the earliest instants of the interaction (our first three points), including when the SFR increases but is already enhanced (20–40 M_⊙yr⁻¹).

• Starburst regime: low α_co, high SFR. These points corresponds to the interacting system during the peaks of star formation.

• Post-burst regime: intermediate α_co, intermediate SFR. This sample exclusively gathers the snapshots made during the decline of the SFR after the first burst and after final coalescence, but without a clear distinction between them.

At the moment of best match with the observations, our simulation gives α_co = 2.9.

Our model of the Antennae starts in the disk regime and remains there ∼20 Myr after the first pericenter passage, that is, 25 Myr after the earliest rise of the SFR. Our sample comprises a point at 40 M_⊙yr⁻¹, but still with a α_co comparable to the high-redshift disk. The first reason for this delay in reaching the starburst regime is likely the non-instantaneous propagation of the triggering mechanism that leads to the burst. When the progenitors first interact, only a fraction of their gas disks is immediately affected by the boost in star formation. It takes several Myr for the enhancement mechanism to propagate across the disks, depending on the orbital configuration of the encounter. The second reason is related to the local modification of the ISM properties (density, temperature, and velocity field) by stellar feedback over timescales of ∼1–10 Myr. Before enhanced star formation spans the majority of the disks’ volume, both these phenomena are local and only affect a fraction of the CO-emitting dense gas.

This is similar to the situation found in our high-redshift clumpy disk, where the star formation activity is restricted to an handful of massive clumps. In the interaction, however, this compact and sparse configuration ceases when a larger portion of the ISM is compressed by turbulence to high densities (Renaud et al. 2014), resulting in wide emission lines and shallower density contrasts between the star-forming and non-star-forming regions than in isolated disks (at low and high redshifts, see Bournaud et al. 2015). The boost in star formation then covers larger volumes such that the system moves to the starburst regime. Why the transition between the two regimes takes place in only a few Myr is unclear, however.

When it reaches starburst mode, the system moves back and forth between the burst and post-burst regimes, without reaching the disk regime again. This evolution is a complex interplay of tidal disruption, debris falling back onto the disks, global decrease of the SFR, cooling of the feedback bubbles, and other mechanisms that drive the system toward lower SFRs and higher α_co, and back when the late galactic encounters occur. The cost of the simulation prevented us from running it for a long time after coalescence, but it is very likely that the merger remnant will return to the isolation phase (via the post-burst regime) to reach a low SFR (≲10 M_⊙yr⁻¹) and high α_co (≈4–5).

When the enhancement in SFR is interpreted as an additive process to the pre-interaction activity, linear relations of the form α_co = a SFR + b are searched for. The parameters of the best fits of these relations (plotted in Fig. 2 and the standard deviation σ are (for α_co and the SFR expressed in M_⊙K⁻¹ km⁻¹ s pc⁻² and M_⊙yr⁻¹):

$$\begin{array}{llll}\mathrm{diskregime}:\hfill & a=-9.5-3\hfill & b=4.35\hfill & (\sigma =0.13)\hfill \\ \mathrm{starburstregime}:\hfill & a=-9.3-3\hfill & b=2.60\hfill & (\sigma =0.15)\hfill \\ \mathrm{post-burstregime}:\hfill & a=-70.2-3\hfill & b=4.52\hfill & (\sigma =0.26).\hfill \end{array}$$

However, the enhancement of star formation can also be seen as a multiplicative boost of the SFR with respect to that expected for main-sequence galaxies. In that case, it is more natural to use power-law relations ${\alpha }_{\mathrm{CO}}=d{\mathrm{SFR}}_{\mathrm{boost}}^c$ to fit our data points, as shown in Fig. 3. Because of the reaction time between the variations of SFR and that of α_co (Fig. 1), measuring the SFR over long timescales (e.g., in the UV or IR) smoothes out its peaks and thus blurs the differences between the starbursts and post-burst regimes. When we consider all phases of the interaction (dashed-dotted lines in Fig. 3), the SFR measured over the last 20 Myr provides the least dispersed relation to α_co, corresponding to the timescale of the propagation of feedback to galactic scales.

Fig. 3. αco as a function of the SFR measured over the last 5 (top panel), 20 (middle panel) and 50 Myr (bottom panel), and normalized to its value in the isolated galaxies. The symbols are the same as in Fig. 1. Dotted, dashed, and dash-dotted lines are fits of the high and low αco regimes (without distinctions between starburst and post-burst points) and the entire sample, respectively. The fitted relation is ${\alpha }_{\mathrm{CO}}=d{\mathrm{SFR}}_{\mathrm{boost}}^c$, with the values c and d indicated in the legend, with the standard deviations (σ).

This relation yields several modes in which the SFR boost is measured over timescales shorter than its variations (≲50 Myr). However, when smoothing it out over longer timescales, the different regimes blend together toward a unimodal relation (but with a larger scatter), in qualitative agreement with that of Sargent et al. (2014, their Fig. 15).

The decrease in α_co in starbursting galaxies can be interpreted as resulting from the observed increase in velocity dispersion of the ISM (Irwin 1994; Ueda et al. 2012), which would lower the absorption along the line of (e.g.,sight Narayanan et al. 2011). Figure 4 shows a more complex picture, however, as no correlation exists between α_co and the velocity dispersion (measured here at the scale of 40 pc and mass-weight-averaged over the entire system). A given velocity dispersion can correspond to a wide range of α_co, and no simple relation can be found, even when the points are split into different regimes. The reason for this is the dual origin of the increased velocity dispersion: turbulence and feedback. The velocity dispersion is enhanced both before and after a peak in SFR by turbulence (which triggers the burst) and the feedback that results from it, respectively (Renaud et al. 2014). However, these epochs of enhanced dispersion correspond to different physical states, and as shown in Fig. 2, to different α_co. Therefore α_co is not a function of the velocity dispersion alone.

Fig. 4. αco as a function of the velocity dispersion of the system. The symbols are the same as in Fig. 1. There is no unequivocal relation(s) between αco and the velocity dispersion.

This argument also suggests a tighter link with feedback, as is indicated also by the delay between the increase in SFR and the drop in α_co noted above. We can consider the SFR normalized by the dense gas mass, that is, the inverse of the depletion time t_dep, as a proxy for the specific energy injected by stellar feedback. Observations report significantly shorter depletion times in interacting systems than in disks (e.g., Daddi et al. 2010b), but these results suffer from uncertainties on α_co. However, simulations do reproduce these differences, without having to choose any α_co. These differences are interpreted as a change in nature of the turbulence that induces efficient compression of the gas (Renaud et al. 2012, 2014; Kraljic 2014). Accounting for the depletion time could thus help to tighten the models of α_co by bringing disks and mergers together.

Figure 5 shows the relation between t_dep and α_co, still with the presence of several modes. The post-burst regime highlighted above (2.5 ≲ α_co,≲ 3.5) is split into two branches of approximately constant t_dep: that at t_dep ≈ 100 Myr corresponds to the separation phase, after the starburst induced by the first passage but before the second passage (see Fig. 1), while that at t_dep ≈ 40 Myr is found after final coalescence. The degeneracy in the α_co – t_dep relation in each of these branches illustrates that α_co is not a function of the feedback injected in the ISM alone (see also Bournaud et al. 2015).

It further shows that successive encounters do not induce the same effects (see, e.g., Renaud et al. 2014 for an Antennae-like system, and Renaud et al. 2018 for a Cartwheel-like galaxy). In the separation phase, only distant, remote triggers such as tidal compression can enhance the SFR, over large scales but at a moderate efficiency, that is, a mildly enhanced SFR and intermediate t_dep. At coalescence, however, nuclear inflows induced by tidal torques are the main trigger of star formation. This local and strong activity efficiently fuels gas and convert it into stars, but over a small volume, which leads to a comparably mild SFR but a shorter t_dep than before. While these differences are visible in Fig. 5, α_co is found to better correlate with t_dep than the SFR.

The timescale for the remnant to reach the disk regime again is uncertain. Lower resolution simulations have shown that the period of reduced t_dep lasts several 100 Myr after the last starburst, even though the SFR recovers its pre-interaction value within only a few 10 Myr (see Kraljic 2014). The reason for this is still debated, and it is possible that a similar effect applies to α_co. Low-resolution simulations like this do not capture the diversity of turbulence in molecular clouds, however, and are thus biased toward universal values of t_dep and α_co such that they cannot be used to answer this point.

Fig. 5. αco as a function of the depletion time from our simulation sample. The symbols are the same as in Fig. 1. The dashed line is the fit αco = 1.33 log(tdep/Myr) + 0.13 (with a standard deviation of 0.34), in the units indicated in the figure.

4. Discussion and conclusions

By applying the LVG method to hydrodynamical simulations resolving molecular clouds, we derived the evolution of the α_co parameter throughout a galactic interaction. Our results confirm the existence of distinct α_co in disks and mergers, even for comparable SFRs (see, e.g., Bournaud et al. 2015). By following a galactic system throughout the course of its interaction and probing the full variation of its SFR, we further identified a third regime that corresponds to the post-burst phase, which can last at least ≈100 Myr. The transitions between the regimes are very rapid, however (∼2–10 Myr). Our main conclusions are as follows.

• The variations of α_co are not a smooth, unimodal function of the SFR because of the modification of the structure and properties of the ISM induced by the interaction, and by star formation and the associated feedback. This leads to three distinct regimes.

• The first regime gathers galaxies for which no interaction has significantly modified the ISM properties over large scales. It thus comprises isolated galaxies and the earliest phases of interactions (while the SFR still increases, but can already be significantly enhanced).

• The starburst regime represents the phases of most active star formation activity, typically ≈20–50 Myr after a close passage, and during coalescence. In this regime, the compression of gas triggers star formation over large volumes, from which the feedback induces a drop of α_co.

• The post-burst regime is found during the separation phase of the galaxies (i.e., between the encounters) and after the final coalescence, while the SFR decreases. α_co increases again, toward the disk regime.

• The transitions from one regime to the next are set by large-scale hydrodynamics effects conveyed to smaller scales by turbulence. The system crosses the starburst regime in several 10 Myr, while it can stay for several 100 Myr in the post-starburst mode. The timescale for the transition is on the order of that of the propagation of stellar feedback, which alters the structure of the ISM.

• As a consequence, the latter two regimes are blended together when estimating the SFR over longer timescales, which results in estimates over 20 Myr providing the tightest α_co -SFR relation. The disks and the earliest phases of the interaction still stand out as a distinct mode, however.

• The depletion time, which varies with the nature of turbulence setting the structure of the ISM, provides a better correlation with α_co than the SFR.

• α_co is not only set by the velocity dispersion or the feedback energy, but is a function of both.

Previous works advocated for a smooth transition of α_co with galaxy properties (Narayanan et al. 2012; Bolatto et al. 2013; French et al. 2015). Our results reveal indeed a more complex picture than a bimodal relation. The post-burst regime we highlight is the signature of this smooth transition. However, we note that the initial transition from the disk to the starburst regime is extremely fast (∼2 Myr), likely because of the interplay of processes acting on short timescales (dynamics of the interaction, increase in SFR, and feedback), but with details that are yet to be identified and understood. This could thus be detected and interpreted as a non-smooth or even discontinuous transition, favoring the idea of bimodality, while in reality, any other transition between regimes is smoother and occurs over longer timescales.

Observational derivations of molecular gas masses based on a universal or bimodal α_co are thus likely to be affected by large errors. By using the relations we provide here on observational data and accounting for the gas fraction and the merger phase rather than the sole SFR or velocity dispersion, the empirical differences in ISM properties and star formation activities in different environments would likely be accentuated. In the case of post-starburst galaxies, for instance, our results suggest that intermediate values of α_co should be adopted for at least several 10 Myr after each burst, and in particular, after final coalescence. In their sample of post-starburst galaxies, French et al. (2015) measured low SFRs (∼0.1 M_⊙ yr, i.e., even lower values than those probed by our simulations¹ which suggests that values of α_co should be adopted that are comparable to or even higher to that of disks. Our findings thus support the high values of α_co they used (∼4–6, instead of choosing an ultraluminous infrared galaxy value of α_co = 0.8, see their Fig. 12), and thus their conclusion that star formation efficiencies are lower in their post-starburst galaxies than in disks despite large amounts of CO gas, for a reason that is yet to be determined (French et al. 2018, but see also Alatalo et al. 2016 for the opposite conclusion from different galaxy selection criteria).

Our finding of several regimes of α_co complements the work of Narayanan et al. (2011) and Narayanan & Krumholz (2014), who adopted a universal sub-grid model to describe turbulence and advocated for a smooth, unimodal transition of the CO-to-H₂ factors between disks and mergers. The distinct regimes we find likely originate from the diversity of turbulence that our simulations capture. For instance, these variations explain that a given SFR surface density can be found in very different media (e.g., a low-redshift merger and a high-redshift disk) because of the variation in depletion time induced by compressive turbulence (Renaud et al. 2014). Thus, our analysis distinguishes galaxies based on the physical trigger of star formation, rather than on the star formation activity itself. As mentioned in Bournaud et al. (2015), the structure of the ISM, including the non-star-forming material, and in particular its clumpy nature, influences the propagation of feedback effects that affect the CO emission (and absorption). Although more investigations are required, this likely explains that t_dep provides a better correlation with α_co than the SFR alone.

The differences between our regimes demonstrate that the turbulent cascade needs to be captured that sets the structure of the ISM down to its dense, cold phase (∼10 K, ∼1–10 pc). However, this implies a high cost that so far prevents modeling many interacting systems. Our conclusions are thus only based on the Antennae-like major merger presented here. While the response of galaxies to interactions strongly depends on orbital and intrinsic parameters, the case presented here includes all the possible phases of interaction and thus provides a wide sampling of the physical conditions that can be expected in mergers.

Our simulations only include a few aspects of feedback, but it is likely that the other mechanisms affect the density structure and the temperature of the star-forming regions and thus the emission CO. Our method provides a value of α_co for a Milky Way-like galaxy that is close to the observational estimate, however, which suggests that including the missing physics would not significantly change our conclusions. In the merger, however, the lack of AGN feedback might play a more important role in regulating the emission from the nuclear regions. Similarly, accounting for different coolants, including dust, would alter the temperature of the ISM and thus the excitation of CO.

The variations found in α_co suggest a tight but complex link between the physical conditions that shape the ISM and drive star formation and CO emission. Understanding them could allow us to invert the analysis and infer the underlying physics from the detection of CO emission lines. For instance, it might be possible to distinguish clumpy disks from mergers using several CO emission lines. Spatial variations of the CO emission, as suggested by Sandstrom et al. (2013), could provide interesting clues, in particular in mergers (Renaud et al., in prep.).

Acknowledgements

We thank the referee for a constructive report. FR acknowledges support from the Knut and Alice Wallenberg Foundation. This work was supported by GENCI (allocations A0030402192 and A0050402192) and PRACE (allocation pr86di) resources.

References

All Figures

Fig. 1. Top panel: evolution of the SFR measured over the last 5, 20, and 50 Myr. Vertical lines indicate the pericenter passages. Symbols mark the instants selected for analysis, with their color-coding time. Their shapes and the shaded areas correspond to the regimes identified below: white, dark gray, and light gray for the disk, starburst, and post-burst regimes, respectively. The black square indicates the instant of best morphological match with the observed Antennae merger. The gray line is the SFR of the progenitor galaxies, run in isolation (i.e., an almost constant SFR of ≈1.5 M⊙yr−1 for each galaxy). Bottom panel: evolution of αco.

In the text

Fig. 2. αco as a function of the SFR (measured over 5 Myr) in the merger and in simulations of other galaxies. (The system starts in the top left corner, and at first order, moves clockwise). The symbols are the same as in Fig. 1. The dashed lines indicate linear fits (αco = aSFR + b) to three regimes: disk galaxies, starbursts, and post-burst phases (see text). The typical dispersions of αco (estimated by varying the line of sight, see text) is shown in the bottom left corner.

In the text

Fig. 3. αco as a function of the SFR measured over the last 5 (top panel), 20 (middle panel) and 50 Myr (bottom panel), and normalized to its value in the isolated galaxies. The symbols are the same as in Fig. 1. Dotted, dashed, and dash-dotted lines are fits of the high and low αco regimes (without distinctions between starburst and post-burst points) and the entire sample, respectively. The fitted relation is ${\alpha }_{\mathrm{CO}}=d{\mathrm{SFR}}_{\mathrm{boost}}^c$, with the values c and d indicated in the legend, with the standard deviations (σ).

In the text

	Fig. 4. αco as a function of the velocity dispersion of the system. The symbols are the same as in Fig. 1. There is no unequivocal relation(s) between αco and the velocity dispersion.
In the text

	Fig. 5. αco as a function of the depletion time from our simulation sample. The symbols are the same as in Fig. 1. The dashed line is the fit αco = 1.33 log(tdep/Myr) + 0.13 (with a standard deviation of 0.34), in the units indicated in the figure.
In the text