Interacting galaxies on FIRE-2: The connection between enhanced star formation and interstellar gas content

TL;DR

This study uses high-resolution galaxy merger simulations to explore how interactions influence star formation and interstellar gas phases, revealing significant increases in cold-dense and hot gas during galaxy interactions.

Contribution

It provides detailed insights into the connection between galaxy interactions, star formation, and multi-phase interstellar medium evolution using FIRE-2 simulations.

Findings

Star formation rate increases by ~30% during interactions.

Cold-dense gas content rises by ~18% in merging galaxies.

Hot gas increases dramatically (~400%) during interactions.

Abstract

We present a comprehensive suite of high-resolution (parsec-scale), idealised (non-cosmological) galaxy merger simulations (24 runs, stellar mass ratio ~2.5:1) to investigate the connection between interaction-induced star formation and the evolution of the interstellar medium (ISM) in various temperature-density regimes. We use the GIZMO code and the second version of the 'Feedback in Realistic Environments' model (FIRE-2), which captures the multi-phase structure of the ISM. Our simulations are designed to represent galaxy mergers in the local Universe. In this work, we focus on the 'galaxy-pair period' between first and second pericentric passage. We split the ISM into four regimes: hot, warm, cool and cold-dense, motivated by the hot, ionised, atomic and molecular gas phases observed in real galaxies. We find that, on average, interactions enhance the star formation rate of the pair (~30%, merger-suite sample average) and elevate their cold-dense gas content (~18%). This is accompanied by a decrease in warm gas (~11%), a negligible change in cool gas (~4% increase), and a substantial increase in hot gas (~400%). The amount of cold-dense gas with densities above 1000 cm^3 (the cold ultra-dense regime) is elevated significantly (~240%), but only accounts for 0.15% (on average) of the cold-dense gas budget.

Tables

Table 1. Table 1 : Simulation specifications: mass resolution (dark matter, gas and stars), highest density and spatial resolution (gas), and gravitational softening (collisionless components and gas).

Property	Value
Mass resolution (dark matter)	$1.9\times {10}^5$ $M_{\odot }$
Mass resolution (gas)	$1.4\times {10}^4$ $M_{\odot }$
Mass resolution (stars)	$8.4\times {10}^3$ $M_{\odot }$
Highest gas density	$5.8\times {10}^5$ ${\mathrm{cm}}^{-3}$
Highest spatial gas resolution	$1.1$ pc
Gravitational softening (collisionless)	0.01 kpc
Gravitational softening (gas)	0.001 kpc

Table 2. Table 2 : Properties of progenitor galaxies: dark matter halo mass, stellar mass, bulge mass, total gas mass, halo spin parameter, stellar disk radial scale length, thickness of stellar disk in units of radial scale length and bulge scale length in units of disk scale length.

Property	Primary	Secondary
$M_{\mathrm{halo}}$	$7.5\times {10}^{11}{M}_{\odot }$	$3.5\times {10}^{11}{M}_{\odot }$
$M_{\mathrm{stellar}}$	$3.0\times {10}^{10}{M}_{\odot }$	$1.2\times {10}^{10}{M}_{\odot }$
$M_{\mathrm{bulge}}$	$2.5\times {10}^9{M}_{\odot }$	$7.0\times {10}^8{M}_{\odot }$
$M_{\mathrm{gas}}$	$8.0\times {10}^9{M}_{\odot }$	$7.0\times {10}^9{M}_{\odot }$
${\lambda }_{\mathrm{halo}}$	0.05	0.05
$R_{\mathrm{disk}}$	2.85 $\mathrm{kpc}$	1.91 $\mathrm{kpc}$
$h_{\mathrm{disk}\left(\mathrm{bulge}\right)}/R_{\mathrm{disk}}$	0.14 (0.13)	0.2 (0.136)

Table 3. Table 3 : Spin-orbit angles for each of our three orientations (approximately prograde, polar and retrograde) used in this work (and in Patton et al., 2013 ; Moreno et al., 2015 ) , originally drawn from Robertson et al. ( 2006 ) . Here, θ 𝜃 \theta and ϕ italic-ϕ \phi specify the total angular momentum of each gas disc relative to the orbital plane. See Figure 1 of Moreno et al. ( 2015 ) for a schematic.

	Prograde (“e”)	Polar (“f”)	Retrograde (“k”)
Primary
${\varphi }_1$	60∘	60∘	-30∘
${\theta }_1$	30∘	60∘	-109∘
Secondary
${\varphi }_2$	45∘	0∘	-30∘
${\theta }_2$	-30∘	150∘	71∘

Table 4. Table 4 : Temperature-density regime demarcations (hot, warm, cool and cold-dense). The bottom two rows split the cold-dense regime into a cold moderately-dense and a cold ultra-dense component (Section  4.2 , Figure  16 ).

ISM regimes	Temperature-density demarcations
warm
cool
cold-dense
hot	$\left(T>{10}^{6}K\right)$
cold moderately-dense
cold ultra-dense

Equations

$$\mathrm{index}={x_{50}}_{x_{25}-x_{50}}^{x_{75}-x_{50}},\mathrm{where}{x}_q=\mathrm{the}{q}^{\mathrm{th}}\mathrm{percentile}.$$

$${\rm index}\Big{|}_{\rm\small SFR\,enhancement}={1.30}^{+0.36}_{-0.24},$$

$$p_{X,Y}=\frac{\mathrm{cov}(X,Y)}{{\sigma }_X{\sigma }_Y},$$

$$\frac{\mathrm{d}{M}_{\alpha }(t)}{\mathrm{d}t}=\sum _{\beta \ne \alpha }{\mathcal{R}}_{\alpha \leftrightarrow \beta }(t),$$

$${\mathcal{R}}_{\alpha \leftrightarrow \beta }(t)=R_{\alpha \to \beta }(t)-R_{\beta \to \alpha }(t),$$

$$R_{\alpha\rightarrow\beta}(t)=\frac{{\rm d}M(t)}{{\rm d}t}\big{|}_{\alpha\rightarrow\beta}$$

Full text

1 Introduction

Galaxy mergers and interactions are an integral component of any galaxy formation model. In the currently accepted $\Lambda$CDM cosmology, structures grow hierarchically, with small objects merging with each other and with larger objects to form the massive galaxies, groups, and clusters that exist in the Universe today (White & Frenk, 1991). Specifically, observations suggest that galaxy-galaxy interactions enhance star formation rates (Ellison et al., 2008, 2013; Patton et al., 2013; Knapen et al., 2015), decrease nuclear metallicities (Kewley et al., 2006; Rupke et al., 2010; Montuori et al., 2010; Torrey et al., 2012; Scudder et al., 2012), and drive AGN activity (Ellison et al., 2011; Khabiboulline et al., 2014; Satyapal et al., 2014; Weston et al., 2017; Goulding et al., 2018). On the numerical side, idealised (non-cosmological) binary merger simulations continue to play an important role in shaping our understanding of the merger process and its relative importance to galaxy formation (e.g., Barnes & Hernquist, 1991, 1996; Mihos & Hernquist, 1996; Iono et al., 2004; Di Matteo et al., 2007, 2008; Teyssier et al., 2010; Torrey et al., 2012; Hayward et al., 2014; Moreno et al., 2015; Blumenthal & Barnes, 2018).

A number of the previously published merger simulation suites employ fairly simple models for the interstellar medium (ISM). One traditional approach is to treat the multi-phase ISM as a single, pressurised fluid (e.g., Springel & Hernquist, 2003; Schaye & Dalla Vecchia, 2008; Murante et al., 2015). This has the advantages of being well behaved numerically (i.e. convergent), computationally inexpensive, and of being based on known empirical relations (e.g., the Schmidt-Kennicutt relation, Schmidt, 1959; Kennicutt, 1998), which can be imposed with care. The principal disadvantage of this approach is that the predictive power of the numerical models is restricted by requiring strict assumptions to be made. For example, in Torrey et al. (2012), the assumed stiffness of the ISM effective equation of state, as well as the mass loading factor for galactic winds, both impact the obtained results.

Models employing a barotropic equation of state have recently pushed resolution to levels capable of resolving Giant Molecular Clouds (GMCs) and the structure of the ISM (Bournaud et al., 2008; Bournaud et al., 2010; Renaud et al., 2008, 2009; Teyssier et al., 2010; Renaud & Gieles, 2013; Renaud et al., 2015). These works highlight the importance of stellar feedback in regulating the ISM structure. Likewise, increasingly detailed models that capture the multi-phase structure of the ISM and adopt feedback-regulated star formation have been developed (e.g., Hopkins et al., 2011; Agertz et al., 2013; Hopkins et al., 2014, 2018b). A key characteristic of these models is a desire to resolve the scales where feedback from young stellar populations first have a global impact on the pressurization of the interstellar medium, whilst simultaneously adopting realistic physical prescriptions at those scales. Since turbulent pressure is believed to play a key role in preventing run-away ISM gas collapse (Ostriker & Shetty, 2011; Hayward & Hopkins, 2017), resolving the injection of momentum and energy sources that drives the supersonic turbulent ISM pressure is key to obtaining a self-regulating ISM. Moreover, it has been found that including local sources of stellar feedback (photoionisation, radiation pressure, supernovae energy/momentum injection, and stellar winds) naturally leads to the launching of large scale galactic outflows (Hopkins et al., 2013a; Hopkins et al., 2013d). The current disadvantage of these more detailed models is that they are computationally expensive, and so a large and comprehensive suite of galaxy merger simulations has not yet been published.

The goal of this paper is to employ modern, state-of-the-art, simulations to investigate the evolution of the ISM in different temperature-density regimes during the merger. We use the GIZMO simulation code (Hopkins, 2015, 2017) - in conjunction with the second version of the ‘Feedback In Realistic Environments’ (FIRE-2) model (Hopkins et al., 2018b) - to build a comprehensive suite of high resolution galaxy merger simulations. This work expands on a smaller set of simulations by Hopkins et al. (2013a); Hopkins et al. (2013d), which use an earlier version of our model for only a handful of mergers. The suite we present here effectively replaces an earlier suite of GADGET-3 (Springel & Hernquist, 2003) galaxy merger simulations created by our group (Patton et al., 2013; Moreno et al., 2015).

Figure 1 displays the various stages of the interaction sequence for one of our simulations. The top two rows show mock ugr SDSS-band composite images. The bottom two rows show the multi-phase structure of the ISM. This colour projection displays gas in terms of temperature in approximately three bins: $T<1000\,K$ (magenta), $T\sim 10^{3}-10^{4}\,K$ (green), and $T>10^{6}\,K$ (red). For this figure, we employ the same visualisation techniques as in Hopkins et al. (2018b) and in other works by the FIRE collaboration. The interaction sequence unfolds as follows. The first pericentric passage produces tidal tails and a bridge (first and third rows), which are more evident in the gaseous component. The galaxies approach again and merge, displaying disturbed morphology first, and ultimately settling down into a disk galaxy with a prominent bulge. The post-merger remnant exhibits faint stellar features, such as shells and streams, as well as a hot gas halo. This work focuses on the early stages of interaction. In the pre-coalescence period, mergers can be relatively ‘cleanly’ identified in observations by searching for galaxy pairs, whereas identifying and timing coalescencing galaxies is more challenging. Thus, by focusing on the pre-coalescence period (hereafter, the ‘galaxy-pair period’), we can more straightforwardly compare with observational constraints on the ISM structure of galaxy pairs.

Section 2 reviews the FIRE-2 model and describes our suite of galaxy merger simulations. Section 3 describes the evolution of star formation, the various temperature-density gas regimes, and their interplay. We present our discussion and conclusions in Sections 4 and 5.

2 Methods

2.1 FIRE-2: The ‘Feedback In Realistic Environments’ Physics Model (Version 2)

The FIRE-2 model (Hopkins et al., 2018b) is very similar to its predecessor (FIRE, Hopkins et al., 2014), but with updated hydrodynamic solvers and supernova coupling algorithms. Differences between these two versions of the FIRE model only have a minor impact on phase mixing in the circumgalactic medium, which is largely irrelevant for this paper. We provide a brief description of the model here, and refer the reader to the above two papers for details. Unlike most recent papers employing the FIRE-2 model, this work is not the result of zoom-in simulations selected from a cosmological box. Instead, we use FIRE-2 physics to run idealised (non-cosmological) galaxy merger simulations***For information on the FIRE Project, visit https://fire.northwestern.edu. Videos of our galaxy merger simulations are available at research.pomona.edu/galaxymergers..

The treatment of radiative heating and cooling includes free-free, photo-ionisation/recombination, Compton, photoelectric, dust-collisional, cosmic ray, molecular, metal-line and fine-structure processes – and it tracks 11 species independently. The code accounts for UV background (Faucher-Giguère et al., 2009) and locally-driven photo-heating and self-shielding. Star formation is constrained to self-gravitating (at the resolution scale), self-shielded gas denser than 1000 cm*-3* (Hopkins et al., 2013b). Stellar feedback mechanisms include momentum flux from radiation pressure; energy, momentum, mass and metal injection from SNe (Types Ia and II) and stellar mass loss (OB and AGB stars). The location where supernovae go off can potentially have an impact on the ISM (Gatto et al., 2015). In our model, this is governed solely by the location of the star particle when the explosion occurs (Hopkins et al., 2018a). Other stellar feedback channels include photoionisation and photo-electric heating. Single stellar population properties are attributed to each star particle. We tabulate stellar masses, ages, metallicities, feedback event rates, luminosities, energies and mass-loss rates with starburst99, a stellar synthesis model (Leitherer et al., 1999). We assume a Kroupa (2001) initial mass function (IMF). Our simulations do not include feedback caused by supermassive black hole (SMBH) gas accretion because the coupling between this type of feedback and the circumnuclear ISM at the resolutions explored here is not well understood yet (but see Torrey et al., 2017). For works on SMBH feedback with FIRE, see Hopkins et al. (2016) and Anglés-Alcázar et al. (2017). Due to the lack of observational constraints, we do not include hot gas atmospheres at the start of the simulation. In this paper we employ the Meshless Finite-Mass (MFM) version of GIZMO, which is useful to keep track of inter-phase gas transitions.

2.2 Suite of Galaxy Merger Simulations

The suite in this paper is similar to a previous merger suite created by our group (Patton et al., 2013; Moreno et al., 2015), but with fewer runs, substantially higher resolution, and with a new physical model. Namely, GIZMO/FIRE-2 (Hopkins et al., 2018b) replaces an effective equation-of-state model with GADGET-3 (Springel & Hernquist, 2003; Springel, 2005). For a thorough comparison of galaxy merger simulations using GADGET-3 and FIRE, we point the reader to Figure 9 of Hopkins et al. (2013a), which shows the bursty nature of SFR in this new model, in line with our work (Section 3.1). Differences between the two versions of FIRE are largely insignificant in the regimes we explore in this paper. Table 1 displays our adopted specifications. Our MFM particle masses are $1.9\times 10^{5}$ $M_{\odot}$, $1.4\times 10^{4}$ $M_{\odot}$ and $8.4\times 10^{3}$ $M_{\odot}$ (dark matter, gas and stellar components). The highest gas density and spatial resolution achieved in our simulations are $5.8\times 10^{5}$ ${\rm cm}^{-3}$ and $1.1$ pc, respectively. Our gravitational softening are 0.01 kpc for the dark matter and stellar components, and 0.001 kpc for the gaseous component.

2.2.1 Isolated Galaxies

Isolated galaxies in our simulations are set up following analytic work by Mo et al. (1998), employing the procedure outlined by Springel et al. (2005). Simulated bulges and dark matter haloes follow a Hernquist (1990) profile. Table 2 describes our two progenitor galaxies, hereafter labelled ‘primary’ and ‘secondary’ (stellar mass ratio $\sim$2.5:1). As we increase the stellar mass, we adjust the corresponding dark matter halo mass, following abundance matching results (Moster et al., 2013). We also adjust the bulge mass so that the bulge-to-total mass ratio decreases, in line with trends in the Sloan Digital Sky Survey (SDSS) (Mendel et al., 2014). We follow Saintonge et al. (2016) to set up our gas masses. For both galaxies, we set up our gaseous components to have disc lengths of $\sim$10 kpc, which yield diameters of $\sim$25 kpc (defined as the radial extents at which the gas surface density falls beneath 1 $M_{\odot}$ pc*-1*), in line with observations by Broeils & Rhee (1997). We also list values for the halo spin parameter, the stellar disk radial scale length, the thickness of stellar disk in units of radial scale length and bulge scale length in units of disk scale length. For a more thorough description on our adopted radial profiles (dark matter halo, stellar disc and bulge, gaseous disc) and the vertical structure of the gas disc, please read Section 2 of Springel et al. (2005). For reproducibility purposes, we point the interested reader to the following website for a complete list of parameters: research.pomona.edu/galaxymergers/isolated-disks-initial-conditions/.

Our isolated galaxies are stable and well-behaved on Gigayear timescales (Section 3.2), in line with earlier tests using the same initial conditions (Hopkins et al., 2012, 2013b; Torrey et al., 2017). We note that, in our set up, we do not start with a multi-phase ISM. Instead, we initially set our gas to $10^{5}$ K with solar metallicity, from which a turbulent ISM emerges as a result of feedback produced by star formation, and in which subsequent metal enrichment is tracked self-consistently. We also note that our initial conditions do not include hot gas atmospheres, and the hot gas described in this paper arises solely from feedback heating the ISM. We elected not to include this hot component because of the lack of detailed hot gas profile measurements around spiral galaxies. Moster et al. (2011) attempt to include hot atmospheres in their merger simulations by adopting $\beta$-profiles from galaxy cluster observations (Cavaliere & Fusco-Femiano, 1976; Jones & Forman, 1984; Eke et al., 1998). However, it is not clear that such an extrapolation is valid at galactic scales. We refer the reader to Section 2.4, where we expand on the suitability of this approximation.

2.2.2 The Fiducial Run

We study a fiducial run in detail to introduce terminology and to provide intuition. For this run, we adopt a nearly prograde orbit with small impact parameter ($\sim$7 kpc) to maximise the effects of the encounter. This limits the duration of the galaxy-pair period; i.e., the time between first and second pericentric passage. To compensate, we adopt a highly eccentric orbit. This choice is slightly different from the fiducial run in Moreno et al. (2015), which only aims to maximise duration. Figure 2 shows the orbital evolution of this run. The top (bottom) panel shows the galaxy-galaxy separation (relative velocity) as a function of time (separation). The vertical dashed lines in the top panel indicate the times at first and second pericentric passages, plus coalescence (defined here as the time at which the two central supermassive black holes are at a distance of 500 pc for the last time). Hereafter, we shift time so that $t=0$ Gyr corresponds to first pericentric passage, unless stated otherwise. Black lines refer to the entire orbit and blue lines correspond to the ‘galaxy-pair period’, corresponding to the times between first and second pericentric passage, the primary focus of this paper. We mask times outside the galaxy-pair period in light grey. The large blue star and pentagon denote orbital properties at first pericentre and apocentre respectively. These markers will be useful when comparing this fiducial run to other simulations in our merger suite. The jagged nature of the lines in the bottom panel is driven by the small-scale dynamical interaction between each the black hole and its neighbouring circumnuclear medium (we smooth our curves using a Savitzky-Golay filter for display purposes only). In this work we use the locations of the supermassive black holes as proxies for galactic centres. We verified the suitability of this approximation visually.

2.2.3 Galaxy Merger Simulations

To explore the effects caused by variations in orbital merging configuration, we constructed a suite of 24 galaxy merger simulations involving our primary and secondary galaxies (Table 2). Following Moreno et al. (2015), we adopt a choice of stellar masses (and mass ratio) representative of close galaxy pair catalogues drawn from the SDSS (Patton et al., 2013) and samples extracted from cosmological simulations with abundance matching (Moreno, 2012; Moreno et al., 2013). Our suite is split into three families, corresponding to three distinct spin-orbit orientations: nearly ‘prograde’, ‘polar’ and ‘retrograde’ (Table 3), corresponding to the “e”, “f” and “k” orientations from Robertson et al. (2006), also used by Patton et al. (2013) and Moreno et al. (2015).

The top panel of Figure 3 shows separation at apocentre versus time at that location for each merger in our suite. Recall that time is rescaled to zero at first pericentric passage. The separation at apocentre and the time it takes the two galaxies to reach that separation are indicative of the orbit’s spatial extent and duration. Large blue circles denote nearly prograde mergers; medium-sized green circles denote nearly polar mergers; and small orange circles denote nearly retrograde mergers. Different shades correspond to different runs within a given orientation. We cover a broad range of separations and timescales. Our fiducial run, represented by the large blue pentagon, is a “middle-of-the-pack” orbit. In Moreno et al. (2015) and Patton et al. (2013), for each orientation we selected a fixed set of specific eccentricities and impact parameters, consistent with published low-resolution cosmological N-body estimates (Khochfar & Burkert, 2006). With the advent of large cosmological simulations of galaxy formation at higher resolution with baryons, however, it is unclear if those original N-body estimates still hold (Rodriguez-Gomez et al., 2015). In this work, we do not fine-tune these orbital parameters. Instead, we tailor our orbits to have certain properties at first pericentric passage, which govern their spatial separation extent and their duration. The bottom panel of Figure 3 shows relative velocity versus separation at first pericentric passage for all of our mergers (the large blue star represents our fiducial run). By design, our 24 mergers arrange themselves into three clusters, two containing 9 mergers and one containing 6 mergers – each with increasing (decreasing) separation (relative velocity) at first pericentric passage: $\sim$(7 kpc, 560 km sec*-1*), $\sim$(17 kpc, 460 km sec*-1*) and $\sim$(28 kpc, 410 km sec*-1*). Each of these groups can be regarded as a set of strong, moderate and weak encounters (D’Onghia et al., 2010). For each orientation, we initially had 9 mergers, but we dropped those orbits with merging times exceeding 5 Gyr (thus reducing our total number of runs from 3$\times$9 = 27 down to 24). This is because it is unlikely for a pair of interacting galaxies to evolve in isolation from the rest of the Universe for such an extended time (Moreno, 2012; Moreno et al., 2013). Also, our low-redshift initial conditions would no longer be applicable for such long-lasting interactions. This explains why the last cluster at large separations in the bottom panel of Figure 3 only contains 6 mergers.

With our range of orbital choices, we are able to cover a broad range of separations and merging timescales. Figure 4 shows this. The 2D histograms in this Figure incorporate all times at and after first pericentric passage, including times well past coalescence (recall that most of the paper focuses only on the galaxy-pair period). The top panel shows separation versus time, colour coded by relative velocity. Our orbital choices produce separations as high as 300 kpc, consistent with observations of wide galaxy pairs by Patton et al. (2013); Patton et al. (2016). At large separations, galaxies slow down relative to their companions, thus spending long periods of time near their orbital apocentres. The bottom panel shows relative velocity versus separation, colour coded by time. Soon after first pericentric passage, galaxy pairs have small separations and high velocities (dark blue bins). At later times, as galaxies approach their orbital apocentres, they tend to have larger separations and smaller relative velocities (green bins, bottom right). Short-lived orbits approach their second pericentres with high relative velocities (blue and cyan bins formed in a diagonal with separations $\sim$30-150 kpc). Comparing Figure 2 to Figure 4 shows that our fiducial run is an ‘average’ orbit. We warn the reader that the use of fixed line-of-sight velocity cuts to select galaxy pairs (e.g., 300 km sec*-1*, Ellison et al., 2008; Patton et al., 2013) may delete true interacting pairs at small separations (upper-left corner).

In summary, we split our suite of merger simulations as follows: 24 mergers $=$ (3 orientations $\times$ 3 first-pericentre separations $\times$ 3 relative velocities at first pericentric passage) $-$ (3 mergers with duration greater than 5 Gyr).

2.3 ISM Temperature-Density Regimes

To study the structure of the ISM, we split our gas component into four temperature-density regimes: hot, warm, cool and cold-dense. See Figure 5 and Table 4 for demarcations. These four temperature-density regimes are motivated by the phases making up the structure of the ISM in observations: the hot, ionised, atomic and molecular gas phases. In this work, we require hot gas to have temperature above 1 million Kelvin. Whilst arbitrary, this choice is widely adopted in the literature because it corresponds approximately to the virial temperature of Milky-Way type haloes (Kereš et al., 2005). We remind the reader that our initial conditions do not include hot atmospheres, rather our hot component is produced solely as a result of feedback heating up the colder ISM components. Below 1 million Kelvin, we split the ISM into three temperature-density regimes. We expect the warm regime to be dominated by warm-ionized gas, indicated by the bright band above 8,000 Kelvin. The lower right corner (with $n>0.1$ and $T<8000$ K) contains the cool and cold-dense regimes. This region has the following structure: a diffuse valley immediately below the 8000-K cut and a cloud peaking at a few hundred Kelvin. For a similar diagram using our model to simulate Milky-Way like galaxies, see Figure 5 of Guszejnov et al. (2017). We expect this corner to contain a mixture of atomic and molecular gas, with the bottom part of this cloud to host most of the molecular gas (see e.g., Mihalas & Binney, 1981; Kulkarni & Heiles, 1988; Draine, 2011). To split this corner into the cool and cold-dense regimes, we follow Orr et al. (2018) and adopt their conservative cut at 300 K. For alternative cold-dense temperature demarcations, see Bournaud et al. (2015).

In this paper, we do not attempt to employ sophisticated (and computationally intensive) models of the ionised, atomic and molecular gas directly. For the next generation of simulations by our group, we are currently refining realistic radiative-transfer methods (Krumholz et al., 2008; Krumholz & Gnedin, 2011; Narayanan & Krumholz, 2017), coupled with chemical-network solvers, to capture the relative importance of each ISM phase (Lakhlani et al., in prep, Richings et al., in prep). Meanwhile, incipient work aimed to validate our FIRE-2 simulations against the ISM in real galaxies has been conducted already. Orr et al. (2018) demonstrate that their treatment of HI and H2 gas yields reasonable agreement with the observed Kennicutt Law (Kennicutt, 1998; Kennicutt & Evans, 2012). Guszejnov et al. (2017) demonstrates that our model successfully reproduces the GMC mass function in the Milky Way (Rice et al., 2016) and the line-width size relation (e.g., the Larson scaling relationship, Larson, 1981) in our Galaxy (Heyer et al., 2009; Heyer & Dame, 2015) and in nearby galaxies (Bolatto et al., 2008; Fukui et al., 2008; Muraoka et al., 2009; Roman-Duval et al., 2010; Colombo et al., 2014; Tosaki et al., 2017). Lakhlani et al. (in prep) is also conducting a through comparison between properties of Milky-Way GMCs in our models and observations in our Galaxy (Dame et al., 2001; Miville-Deschênes et al., 2017). Exporting these methods to our galaxy mergers is in the plans. Our main goal here, however, is to explore the behaviour of gas in various temperature-density regimes with a minimum set of assumptions.

2.4 Caveats and Limitations

As with any numerical simulation, our work has limitations that impact the extent to which our results should be interpreted and compared with observations. First, it is difficult to compare our various temperature-density ISM regimes directly against observations without employing full radiative-transfer calculations coupled with chemical network solvers. Additionally, our simulations lack feedback from SMBH accretion, which is likely to affect the ISM. For galaxy pairs, this effect is likely to be small (Treister et al., 2012). Likewise, we do not include hot gas atmospheres at the start of our simulations, which might alter the evolution of our hot gas component and other phases. Namely, cooling of hot gas may feed the ionised gas component (and, subsequently, other gas phases) by ‘hot gas accretion’ (Kereš et al., 2005; Moster et al., 2011; Karman et al., 2015). Another limitation in our approach is the lack of cosmological context. Interacting galaxies embedded in the cosmic web can increase their cold-gas supply via ‘cold mode’ accretion (Kereš et al., 2009) and their dynamics can be influenced by third bodies (Moreno, 2012; Moreno et al., 2013). A promising approach to avoid these problems is to extract merging systems from cosmological simulations (Tonnesen & Cen, 2012; Sparre & Springel, 2016; Bustamante et al., 2018; Hani et al., 2018, Patton et al., in prep; Blumenthal et al., in prep). Additionally, placing mergers in a cosmological setting helps us identify false pairs (projected interlopers) or contamination by pairs selected before their first approach (Kitzbichler & White, 2008; Hayward et al., 2013; Robotham et al., 2014).

Nevertheless, idealised binary merger simulations have significant advantages. First, this method allows us to resolve the ISM to extremely high densities. Secondly, this ‘experimental’ set up affords us control on the types of progenitors and orbital parameters we employ (e.g., Cox et al., 2006; Di Matteo et al., 2007, 2008; Moreno et al., 2015), and is ideal for numerical experiments aimed at simulating specific well-known systems (Karl et al., 2010; Karl et al., 2011; Karl et al., 2013; Privon et al., 2013; Renaud et al., 2015; Lahén et al., 2018). In order to maximise the advantages, and mitigate the limitations of the isolated binary-merger method, we take the following steps. First, we avoid long-lived galaxy-galaxy interactions, for which gas accretion from the cosmic web or interactions with external ‘third’ galaxies cannot be ignored. Secondly, in this work we quantify interaction-induced effects by comparing merging systems against isolated ‘control’ galaxies (e.g., by reporting SFR enhancements and mass excesses). This approach mitigates the effects caused by other environmental factors. On average, if matched correctly, we expect real non-interacting and interacting galaxies to experience the same level of external effects – i.e., cosmic web gas fueling, interactions with third bodies, etc. – as they evolve. Thus, even though real galaxy pairs are embedded in the cosmic web, our idealised binary-merger simulations provide a reasonable approximation as long as we report quantities in terms interacting-to-isolated ratios and avoid long-lived interactions.

3 Results

3.1 Fiducial Run: Star Formation

We summarise the change in SFR experienced by our simulated galaxies during their interaction. To compare interacting galaxies to their ‘control’ isolated counterparts, we define SFR ‘enhancement’ as the ratio of SFR in the interacting galaxies to the sum of the SFRs in the two isolated galaxies. This mimics what we do in observations (Ellison et al., 2008; Scudder et al., 2012). The rationale behind this is to tease out the effects caused by interactions from those caused by other environmental effects (Patton et al., 2016). The only difference here is that, instead of looking at individual galaxies, we calculate this ratio for the entire galaxy-pair system. The top panel of Figure 6 shows the evolution of star formation rate (hereafter SFR) in our fiducial merger. The jagged thin grey solid line shows outputs at every 5 Myr, demonstrating the bursty nature of FIRE galaxies (Hopkins et al., 2014; Orr et al., 2017; Sparre et al., 2017; Faucher-Giguère, 2018)– in contrast to previous works (i.e., Mihos & Hernquist, 1996; Torrey et al., 2012; Moreno et al., 2015). Hereafter, we interpolate our curves at every 25 Myr for easy display purposes (thick solid purple line). The dashed lines represent the sum of the SFR in the two corresponding isolated galaxies. The vertical dashed lines represent first and second pericentric passages, plus coalescence. Regions outside the galaxy-pair period are masked in light grey. The bottom panel of Figure 6 shows star formation rate enhancement versus time. Overall, interactions elevate SFR in galaxies, in qualitative agreement with previous work (Mihos & Hernquist, 1996; Di Matteo et al., 2007, 2008; Teyssier et al., 2010; Renaud & Gieles, 2013; Moreno et al., 2015). At first pericentric passage, the SFR exhibits a sudden spike, followed by a prolonged period of enhancement (by factors of $\sim$2-3), especially between t$=$0 and $\sim$1.3 Gyr. Although not the central focus of the paper, we briefly comment on the sudden rise in SFR and SFR enhancement at second pericentric passage. For our fiducial run, these quantities increase by a factor of $\sim$20. After checking, about half of our runs exhibit similar, albeit generally weaker, upturns. The intensity of SFR and SFR enhancement at second passage and coalescence depends on several factors, including internal properties of the colliding galaxies (their ISM structure) and the geometry of collision at second approach.

3.2 Fiducial Run: The Structure of the ISM

One of the advantages of the FIRE-2 model is its ability to capture the multi-phase structure of the ISM. The central focus of this paper is to study the impact of the encounter on various temperature-density regimes. Table 4 and Figure 5 contain our adopted demarcations. Figure 7 shows the evolution of gas mass in our fiducial run, split into four ISM temperature-density regimes: warm (orange), cool (green), cold & dense (blue) and hot (red lines) - displayed from top to bottom in decreasing order of fractional contribution to the entire gas mass content. This approximate ordering is maintained for the isolated galaxies (dashed curves), and for most of the simulation time for the interacting system (solid lines). In all cases, the isolated runs are stable and well-behaved on timescales of several Gigayears. The bottom panel shows mass ‘excesses’ for each regime, defined as the ratio of the masses in the interacting system divided by the sum of the masses in the two isolated galaxies. For the interacting runs, the four regimes exhibit the following behaviour:

• •

**Warm gas (fiducial run):

**Most of the gas mass is in this regime (orange lines). Warm gas is gradually depleted as a function of time in both isolated and interacting galaxies. This depletion is magnified by interactions, especially when the two galaxies are in close proximity (except briefly at $\sim$0.2-0.3 Gyr).

• •

**Cool gas (fiducial run):

**Cool gas (green lines) mass is depleted in both interacting and isolated cases. For the isolated case, this is slow and steady, whereas the interacting case undergoes more dramatic changes: a brief boost, followed by a drop and a long-term steady recovery.

• •

**Cold-dense gas (fiducial run):

**This component (blue lines) experiences depletion over long timescales. The encounter provokes the following behaviour: a brief and sudden spike, followed by a mild and brief suppression. Soon after, a cold-dense gas reservoir is replenished and survives (whilst experiencing slow depletion) until the end of the galaxy-pair period.

• •

**Hot gas (fiducial run):

**This component (red lines) comprises the smallest contribution to the total gas mass budget (recall that we did not include a hot atmosphere in our initial conditions). At both pericentric passages, and at coalescence, the amount of hot gas increases dramatically, possibly due to shock heating. Hot gas excess appears before the actual pericentric passages, indicated by the first two vertical lines, because this process begins as soon as the outer regions of the extended gas component collide with one another. The first image on the third row of Figure 1 shows this. During the period between first and second pericentric passage, an excess of hot gas is maintained. This excess approximately doubles between t$\sim$1.3-1.9 Gyr. This phase experiences the largest levels of interaction-induced enhancement.

3.3 Merger Suite: Star Formation

Figure 8 shows star formation rate enhancement versus time for the entire merger suite. We only show the galaxy-pair period, between first and second pericentric passage, whose duration varies from merger to merger. This is a generalised version of Figure 6 (excluding times masked in light grey). Darker hexagons indicate higher incidence of interacting systems per bin, displayed on a logarithmic colour scale. The thick solid black line indicates the median, and the thin solid lines encompass the top and bottom quartiles. On average, star formation is enhanced across our merger suite, in qualitative agreement with observations (Woods et al., 2006; Ellison et al., 2008; Lambas et al., 2012; Scott & Kaviraj, 2014; Patton et al., 2013; Patton et al., 2016) and earlier simulations (Barnes & Hernquist, 1996; Di Matteo et al., 2007, 2008; Teyssier et al., 2010; Renaud & Gieles, 2013; Patton et al., 2013; Moreno et al., 2015). The level of enhancement, and the scatter around the median, diminishes with time. We report the median SFR enhancement, and the corresponding upper and lower quartiles, across the entire sample (including all times within the galaxy-pair period) in the upper-right corner of the figure. Hereafter, we display this information using an index with this format:

[TABLE]

To illustrate, our median SFR enhancement is represented as

[TABLE]

where $x_{50}$ $=$ 1.30, with upper quartile at $x_{75}$ $=$ 1.66 $=$ $x_{50}$ $+$ 0.36 and lower quartile at $x_{25}$ $=$ 1.06 $=$ $x_{50}$ $-$ 0.24. We report similar indices to indicate gas mass excesses in the ISM (Figures 9 and 17). Our presentation combines several mergers with different time durations between first and second pericentric passage. At late times, these results are dominated by those mergers with longer durations. A few alternative ways to present our results include: (1) using ‘merging stages’ (i.e., as in Veilleux et al., 2002; Haan et al., 2011; Barrera-Ballesteros et al., 2015; Hung et al., 2015; Larson et al., 2016; Smith et al., 2018); (2) rescaling time by dividing by the time between first and second passage (Privon, 2014); and (3) using three-dimensional or randomly-selected projected separations (Patton et al., 2013). We elect to present our results as a function of time to display how these physical processes unfold more clearly.

3.4 Merger Suite: The Structure of the ISM

Figure 9 shows the time evolution of gas mass excess in various temperature-density regimes in our merger suite (galaxy-pair period only). It generalises the bottom panel of Figure 7 (excluding regions masked in light grey). From top to bottom, panels are organised by regime: warm (orange), cool (green), cold-dense (blue) and hot (red). Thick (thin) solid lines refer to the median (top and bottom quartiles). Vertical axes on each panel have different ranges for display purposes. The index in the upper-right corner shows the median and quartiles (equation 1). Each phase exhibits the following behaviour:

• •

**Warm gas (merger suite):

**In general, interactions suppress the warm gas mass content (sample median $=$ 0.89). The intensity and duration of this depletion varies from merger to merger.

• •

**Cool gas (merger suite):

**Cool gas suppression is followed by a slow and steady recovery. By t$\sim$1 Gyr and after, most mergers in our suite exhibit mild cool gas excess ($\lesssim$10% above unity). This depletion plus recovery results in a negligible overall mass excess (sample median $=$ 1.04).

• •

**Cold-dense gas (merger suite):

**Interacting galaxies in our suite exhibit an excess in cold-dense gas at all times (sample median $=$ 1.18). After a brief dip, the establishment of a reservoir of cold-dense gas occurs, which survives for several Gyrs.

• •

**Hot gas (merger suite):

**This component experiences the highest levels of mass excess (sample median $=$ 4.90), especially at pericentric passages. After the large boost at first pericentre, the hot gas excess stabilises to values generally above unity ($\sim$2-4). Dramatic spikes are driven by hot gas boosts at second pericentre for those mergers experiencing second approach within that window of time. (We remind the reader that our initial conditions do not include a hot gas atmosphere.)

3.5 Merger Suite: Star Formation and its Connection to the ISM

This section explores possible correlations between star formation and the evolution of the various temperature-density regimes in the ISM. Our goal is quantify potential correlations between SFR enhancement and gas mass enhancements or suppressions. Figure 10 shows SFR enhancement versus mass excess for each regime (top-to-bottom): warm (orange), cool (green), cold-dense (blue) and hot (red). We display the corresponding Pearson correlation coefficient on the top corner of each panel. This is defined as

[TABLE]

where cov$(X,Y)$ is the covariance between $X$ and $Y$, and $\sigma_{X}$ and $\sigma_{Y}$ are the standard deviations of $X$ and $Y$. Here, $\{p=1,-1,0\}$ refer respectively to perfect correlation, perfect anti-correlation, and zero correlation. We find the following connections to SFR enhancement:

• –

No correlation with warm gas mass suppression (upper-left panel, Pearson coefficient $=$ -0.13).

• –

Weak anti-correlation with cool gas mass excess (upper-right panel, Pearson coefficient $=$ -0.27).

• –

Weak correlation with cool-dense gas mass excess (lower-left panel, Pearson coefficient $=$ +0.31).

• –

No correlation with hot gas mass excess (lower-right panel, Pearson coefficient $=$ +0.13).

On a final note, we comment on the bimodal structure of the 2D histogram in Figure 10 corresponding to the warm gas component (top panel). We find (not shown) that the peak to the right of $M/M_{\rm isolated}$=0.9 is dominated by retrograde mergers, whilst the peak to the left of $M/M_{\rm isolated}$=0.9 is dominated by prograde and polar mergers. In other words, our simulations suggest that spin-orbit orientation governs the effectiveness of warm gas depletion. Older binary merger simulations show that the spin-orbit orientation plays a role in establishing the following: (1) the level of star formation enhancement (Di Matteo et al., 2007; Moreno et al., 2015); (2) the spatial extent of star formation (Moreno et al., 2015); and (3) the likelihood of the remnant to be a fast or a slow rotator (Bois et al., 2011). Recently, cosmological simulations have probed the physical role of orbital parameters in merging systems. Lagos et al. (2018, EAGLE) show that co-rotating mergers are more conducive to the enhancement of the stellar specific angular momentum, whilst Li et al. (2018, Illustris) find that most prolate galaxies are produced by nearly radial orbits. We defer a systematic study on the dependance of orbital parameters to a future paper.

3.6 Inter-Regime Transition Rates

To obtain a better physical picture, we investigate in detail how the various gas regimes feed and drain one another as star formation unfolds. To achieve this, in this section we explore the behaviour of mass transition rates, and how these build up each gas component. Here, we focus only on the fiducial case because we have shown that this run successfully captures the most relevant trends in our merger suite (Figures 6 and 7 versus Figures 8 and 9). Focusing on a single run allows us to provide a detailed narrative of some of the changes each regime experiences; although the details of these will be different for different mergers, the salient characteristics can be found throughout our merger suite.

Figures 11, 12, 13 and 14 show mass transition rates relevant to each of our four temperature-density regimes in the ISM: the warm, cool, cold-dense and hot gas component respectively. Vertical dashed lines indicate first and second pericentric passages, plus coalescence. Periods of time outside the galaxy-pair period are masked in light grey. We use $M_{\alpha}(t)$ to represent mass as a function of time (top panel of Figure 7), where Greek indices refer to our four temperature-density gas regimes, plus new stars. The time derivatives, ${\rm d}M_{\alpha}(t)/{\rm d}t$, are shown in the top (bottom) panels as thick solid (dashed) lines fluctuating about zero for the interacting (isolated) case. These fluctuations produce upturns and downturns in $M_{\alpha}(t)$. Each derivative can expressed as the sum of net inter-regime transition rates:

[TABLE]

where ${\cal R}_{\alpha\leftrightarrow\beta}$ represents the net rate at which two temperature-density regimes “$\alpha$” and “$\beta$” (labelled “$\alpha$” $\leftrightarrow$ “$\beta$” in the Figures) exchange mass. By definition,

[TABLE]

where

[TABLE]

is the mass transition rate from regime “$\alpha$” to regime “$\beta$”. In practice, to calculate these rates, we employ particle IDs across consecutive snapshots (separated by 5 Myr) and match gas elements that covert from one gas regime into another, or into new stellar particles. Here, net transition rates are shown as thin solid (dashed) lines in the top (bottom) panel for the interacting (isolated) case. For display purposes, we smooth over 25-Myr timescales. Different time binning changes the ‘jaggedness’ of the curves, but does not alter our main results. Positive (negative) net transitions correspond to curves above (below) zero. Note that net transition rates are anti-symmetric: ${\cal R}_{\alpha\leftrightarrow\beta}=-{\cal R}_{\beta\leftrightarrow\alpha}$; i.e., vertical mirror images of the inverse transition. For example, the green line in Figure 11 (top panel, labelled “warm $\leftrightarrow$ cool”) is the negative (vertical mirror image) of the orange line in Figure 12 (top panel, labelled “cool $\leftrightarrow$ warm”).

Before describing how each gas regime exchanges mass with others in detail, a few comments are in order. First, we emphasise that our figures only display net transition rates, not regular rates – i.e., ${\cal R}_{\alpha\leftrightarrow\beta}$ in equation (5), as opposed to ${R}_{\alpha\rightarrow\beta}$ in equation (6). Secondly, whilst it is true that flows between two gas regimes occurs in both directions, our simulations show that inter-regime transitions tend to favour a preferred direction (thin lines rarely cross the horizontal zero line). We discuss specific cases where halting and reversal occurs, resulting in interesting behaviour. Thirdly, we note that the main effect caused by encounters is the amplification of net transition rates: solid lines in top panels (corresponding to interacting galaxies) versus dashed lines bottom panels (corresponding to isolated galaxies). We devote most of the discussion to merging galaxies during this period of time. Lastly, we note that the large fluctuations before first pericentric passage (left of the first dashed vertical line, top panel of Figure 7) are caused primarily by the ISM of the two galaxies coming to contact before the two centres (marked by the positions of the two supermassive black holes) reach their minimum separation for the first time. The first two panels on the third row of Figure 1 illustrate this. A secondary effect is the fact that, for the very first few snapshots, the galaxies in our simulation are becoming stable (i.e., spikes on the left-end of the bottom panels). We verified that stabilization-driven fluctuations disappears before the two galaxies come into first contact.

• •

**Warm gas (transition rates, fiducial run):

**Figure 11 shows the evolution of the warm gas component and its interplay with other ISM regimes. In the top panel, the thick solid orange line represents the the time-derivative of mass in this regime as a function of time. This quantity is governed primarily by the competition of two inter-regime transitions: a net gain caused by cold-dense gas converting into warm gas (thin blue line) and a net loss caused by warm gas converting into cool gas (thin green line). In general, the rate at which warm gas turns into cool gas exceeds the rate at which cold-dense gas turns into warm gas, resulting in an overall slow depletion of warm gas (top panel). Deviations from this general behaviour are caused by special circumstances. For instance, the brief recovery at t$\sim$0.3 Gyr (dashed and solid orange lines touching, top panel of Figure 7) is driven by periods in which the net transformation of warm gas into cool gas is halted or reversed (green lines in top panel crossing into the positive side). Namely, this recovery is caused by brief periods where both the cool and cold-dense gas become warm gas – possibly due to intense stellar feedback.

• •

**Cool gas (transition rates, fiducial run):

**Figure 12 shows the evolution of the cool gas component and its interplay with other ISM regimes. In the top panel, the thick solid green line represents the the time-derivative of mass in this regime as a function of time. This quantity is governed primarily by the competition of two inter-regime transitions: a net gain caused by warm gas converting into cool gas (thin orange line) and a net loss caused by cool gas converting into cold-dense gas (thin blue line). In general, the rate at which cool gas turns into cold-dense gas exceeds the rate at which warm gas turns into cool gas, resulting in a slow depletion of cool gas (top panel). Deviations from this general behaviour are caused by special circumstances. The boost in cool gas right after first passage (solid green line, top panel of Figure 7) is caused by a high net influx of warm gas that is not compensated by a high-enough net loss into cold-dense gas (top panel, green spikes at t$\sim$0 Gyr). Also, strong dips in cool gas are caused by periods in which the net transformation of warm gas into cool gas is halted or reversed (top panel, orange lines crossing into the negative side).

• •

**Cold-dense gas (transition rates, fiducial run):

**Figure 13 shows the evolution of the cold-dense gas component and its interplay with other ISM regimes. In the top panel, the thick blue line represents the the time-derivative of mass in this regime as a function of time. This quantity is governed primarily by the competition of three transitions: a net gain caused by cool gas converting into cold-dense gas (thin green line), a net loss caused by cold-dense gas converting into warm gas (thin orange line), and the consumption of cold-dense gas into new stars (thin brown line). In general, the net rate at which cold-dense gas turns into warm gas and new stars exceeds the net rate at which cool gas turns into cold-dense gas, resulting an overall slow depletion of cold-dense gas (top panel). Note that depletion of cold-dense into warm gas tends to dominate over consumption of cold-dense gas into new stars. Deviations from this general behaviour are caused by special circumstances. The boost in cold-dense gas right after first pericentric passage (solid blue line, top panel of Figure 7) is caused by a high net influx of both cool and warm gas, possibly because gas is being compressed by the encounter (green and yellow peaks near the first vertical dashed line, middle panel). Another example occurs at $\sim$0.4 Gyr after first pericentric passage. For a brief period of time, the conversion of cold-dense into warm gas is diminished, whilst the conversion of cool into cold-dense gas increases, causing the build up of a cold-dense gas reservoir (orange lines approaching zero from below, middle panel).

• •

**Hot gas (transition rates, fiducial run):

**Figure 14 shows the evolution of the hot component and its interplay with other ISM regimes. In the middle panel, the thick solid red line represents the the time-derivative of mass in this regime as a function of time. This quantity is primarily governed by the net mass interchange between the hot and warm regimes. In particular, the boosts and drops in hot gas right before and after pericentric passages are caused by large net exchanges with the warm gas component. In other words, the hot gas “follows” the warm gas. Because the mass in hot gas is orders of magnitude smaller than the mass in warm gas, the warm gas component is largely unaffected by the fluctuations experienced by the hot gas component. Exceptions to this behaviour in which warm gas “controls” the hot gas take place in the early stages of interaction. During the first $\sim$0.4 Gyr after first passage, the sudden drop in hot gas (caused by a net transformation into warm gas) is mildly mitigated by small gains caused by the net conversion of some cold-dense and cool gas into hot gas (solid red line, top panel of Figure 7). We speculate that this is driven by feedback associated with enhanced star formation. These net losses do not affect the cold-dense and cool components because the hot gas mass budget is significantly smaller than the budget contributions from those two components. The amplitude of phase transitions involving hot gas are generally small - except near pericentric passages, where hot gas is produced as gas located in the outskirts of galaxies is shock heated by the encounter.

4 Discussion

4.1 An Emerging Picture

The central goal of this paper is to investigate how various temperature-density regimes in the ISM interact with one another, and their connection to star formation. The FIRE-2 model, which has the ability to resolve the multi-phase structure of the ISM, is ideal to achieve this goal. By investigating net transition rates between ISM temperature-density regimes, we learn that the content in each regime is governed by the competition of net transitions between each regime and two other “neighbouring” regimes:

• •

cold-dense $\rightarrow$ warm              $\rightarrow$ cool,

• •

warm             $\rightarrow$ cool                  $\rightarrow$ cold-dense,

• •

cool                 $\rightarrow$ cold-dense $\rightarrow$ warm ($+$ new stars),

with the hot component feeding from the other regimes:

• •

warm ($+$ cool and cold-dense) $\rightarrow$ hot.

We illustrate this with a Sankey flow diagram in Figure 15. The thickness of each pipe approximately represents typical transition rate amplitude (not drawn exactly to scale). At the end of the pipeline, a fraction of the cold-dense gas turns into new stars, but most of it loops back into warm gas. We find that galaxy-galaxy interactions are responsible for the following interesting effects:

• •

Warm gas depletion is amplified.

• •

Cool gas depletion is amplified, especially early.

• •

A reservoir of extra cold-dense gas is enhanced.

Inter-regime net transition rates offer a physical explanation. Namely, whilst the direction of flow in this pipeline is approximately steady and one-directional, interactions have the potential to accelerate, halt, or even reverse the direction of these inter-regime transitions.

4.2 Cold Ultra-Dense Gas

Our simulations show that the gradual depletion of cold-dense gas is driven primarily by a net conversion into warm gas. Occasionally, depletion is halted and reversed, resulting in the replenishment of a cold-dense gas reservoir, which is depleted on much longer timescales. This build up of cold-dense gas is driven by a brief period in which the transition from cool to cold-dense gas dominates. Figure 10 suggests that cold-dense gas excess and SFR enhancement in interacting galaxies are mildly correlated. By construction, SFR in the FIRE-2 model occurs in only the very densest portions of the cold-dense gas budget, as expected by observations (Lada et al., 2012). The focus of this section is to discuss the behaviour of cold ultra-dense gas (with $n>1000\,{\rm cm}^{-3}$) in our simulations.

Figure 16 (top panel) shows the evolution of cold-dense gas (blue) in our fiducial run, split into cold moderately-dense gas ($n=10-1000\,{\rm cm}^{-3}$, cyan) and cold ultra-dense gas ($n>1000\,\,{\rm cm}^{-3}$, magenta). The density cut at $n=1000\,{\rm cm}^{-3}$ captures all the star-forming cold-dense gas (by construction). Note that the ultra-dense gas component between apocentre and second passage is close to our mass resolution limit ($1.4\times 10^{4}$ $M_{\odot}$, Table 1), resulting in broken lines. The middle panel shows the fraction of cold-dense gas in the cold moderately-dense (cyan) and cold ultra-dense (magenta) components. Overall, during the galaxy-pair period, cold ultra-dense gas only accounts for, at most, a few percent of the entire budget. After apocentre, it drops below a thousandth. The bottom panel shows interaction-induced mass excesses. Whilst the cold-dense and cold moderately-dense budgets are elevated by factors of only $\sim$1-2, the cold ultra-dense gas mass increases by a factor higher than $\sim$10.

Figure 17 explores the behaviour of the cold ultra-dense gas content across our entire merger suite. The top panel shows the fraction of cold-dense gas with densities higher than $n=1000\,\,{\rm cm}^{-3}$. We find that, on average, only $\sim$0.15% of the cold-dense gas mass in our simulations achieves such high densities. The middle panel shows mass excess. We find that interactions strongly elevate the amount of ultra-dense gas in galaxies, by a factor of $\sim$3.41 on average. The bottom panel explores the possibility of a correlation between SFR enhancement and excess of ultra-dense cold gas in interacting galaxies. By construction, and motivated by observations (Lada et al., 2012), cold ultra-dense gas, as selected here, corresponds to star-forming gas. However, we only find a mild correlation (Pearson coefficient $=$ 0.2) when ‘enhanced’/‘excess’ quantities are considered. We speculate that this scatter is driven by the fact that we are exploring the high-density tail of the gas density function, in a regime where our simulated galaxies contain only a handful (if any) of gas elements at those densities. Higher resolution is required to accurately investigate this correlation.

4.3 Connection to Observations

Our suite of simulations shows that close interactions enhance star formation in galaxies by $\sim$30% on average (Figure 8), confirming predictions by earlier simulations (Mihos & Hernquist, 1996; Di Matteo et al., 2007, 2008; Teyssier et al., 2010; Renaud & Gieles, 2013; Hayward et al., 2014; Moreno et al., 2015), which are qualitatively consistent with observations (Woods et al., 2006; Ellison et al., 2008; Lambas et al., 2012; Ellison et al., 2008, 2013; Patton et al., 2013; Scott & Kaviraj, 2014). Beyond this, the central goal of this paper is to connect star formation enhancement with the behaviour of the various temperature-density ISM regimes in interacting galaxies. Because of its ability to capture the multi-phase structure of the ISM, FIRE-2 is well positioned for this purpose (Hopkins et al., 2014, 2018b). Below we discuss connections between our simulation work and observations.

• •

**Cool gas (connection to observations):

**We predict that close interactions do not elevate or suppress cool-gas mass levels in galaxies significantly during the galaxy-pair period (median $\sim$4% increase – Figure 9, upper-right panel). This median value is the combination of a brief period of mass suppression during the first Gyr after first pericentric passage, followed by a recovery period in the cool gas mass content. Beyond $\sim$1.5 Gyr after first passage, the majority of our mergers exhibit elevated cool gas mass content. To connect the evolution of this component with observations, HI is the standard tracer of cool-gas in galaxies.

There are currently conflicting indications in the observational literature concerning the HI fraction in merging galaxies. Some studies (e.g., Braine & Combes, 1993; Ellison et al., 2015; Zuo et al., 2018) have found no difference in the HI fractions of merging and control galaxies, whereas others (Casasola et al., 2004; Huchtmeier et al., 2008; Jaskot et al., 2015; Janowiecki et al., 2017, Ellison et al. 2018) find enhanced gas fractions. One of the main observational challenges in determining the HI content in observations is that the vast majority of data are obtained with single dish telescopes. Such facilities have very large beams, often several arcminutes in diameter. Resolving the two components in an interacting pair is therefore often not possible. Although corrections can be attempted for this blending (e.g., Zuo et al., 2018), it is nonetheless an uncertainty in the analysis. Along this vein, cosmological simulations show that this effect can overestimate the neutral-hydrogen content in galaxies, especially in satellites (Stevens et al., 2018, IllustrisTNG). Alternatively, the blending problem can be mitigated with high spatial resolution interferometric observations. Some individual galaxies have been mapped in detail this way (e.g., Hibbard et al., 1994; Koribalski & Dickey, 2004). Whilst such studies can reveal the detailed structure of HI in mergers, they cannot be used for statistical determination of changes in gas fractions. Scudder et al. (2015) presented a larger sample of close pairs observed with the Very Large Array (VLA), in order to investigate how star formation rate enhancement correlates with gas fraction. However, the lack of a suitable control sample meant that Scudder et al. (2015) were not able to quantify how enhanced/depleted these galaxies were compared to non-interacting galaxies.

• •

**Cold-dense gas (connection to observations):

**Our simulations predict that close interactions elevate cold-dense gas mass content in galaxies during the galaxy-pair period (median $\sim$18 % increase – Figure 9, lower-left panel). If we follow other published works in the literature (e.g., Bournaud et al., 2015; Orr et al., 2018) and treat our cold-dense regime is a proxy for molecular gas, our results are qualitatively in line with observations. Using a survey of CO(1-0) and CO(2-1) observations with the IRAM 30 meter telescope, Braine & Combes (1993) finds that CO luminosity is enhanced in tidally disturbed galaxies. Combes et al. (1994) find similar results using a sample of IRAS-detected galaxies in binary systems. These authors also find a strong correlation between normalised CO and FIR luminosities, suggesting that enhancement in molecular gas content is responsible for the triggering of star formation in interacting galaxies. More recently, Violino et al. (2018) used IRAM 30-m CO(1-0) observations of 11 galaxies with close companions to show that interactions boost the H2 gas content in these systems. They find an increase of 0.4 dex, slightly higher than our excess of $\sim$18% (Figure 9, lower-left panel). Kaneko et al. (2013), and Ueda et al. (2014) find similar results.

Although there is a broadly consistent qualitative agreement between our simulations, and observations, that the cold-dense gas component (traced by the molecular phase) is enhanced by the interaction, several caveats are necessary on both sides. For example, having defined the cold-dense component (Figure 5), it is trivial for us to recover the mass in this regime. However, observations must convert an observed CO luminosity to a molecular gas mass, via the adoption of a conversion factor ($\alpha_{CO}$, Narayanan & Hopkins, 2013; Bolatto et al., 2013). The value of $\alpha_{CO}$ is notoriously different in ‘normal’ disk galaxies and ULIRGs, and the choice of an appropriate value for mergers is germane (Narayanan et al., 2011, 2012). Secondly, our experimental design, which involves evolving the interacting galaxies also in isolation, makes a quantitative assessment of the merger straightforward. Observations need to carefully construct their control samples in order to avoid biases. Several recent papers (e.g., Violino et al., 2018, Pan et al., submitted) have made important progress in recognizing these vital issues. On the simulation side, a more direct comparison would require full radiative-transfer calculations of the $\alpha_{CO}$ factor. Work in this direction is currently in progress by our group (Bueno et al., in prep). Our simulations also show that excess in the cold-dense gas mass content is mildly correlated with enhanced SFR (Pearson coefficient $=$ +0.31, Figure 10). Our fiducial run also show connection between these two quantities: interactions simultaneously boost star formation and form a reservoir of cold-dense gas (Figures 6 and 7). However, even when star formation enhancement peters out, this reservoir remains. In general, this is accompanied by a slow recovery of a cool gas reservoir (Figure 9). In other words, as in the case of cool gas (and HI-gas in observations, Ellison et al., 2018), our simulations suggest that the cold-dense (or cool) gas content is not exhausted (or expelled) after episodes of elevated star formation in interacting galaxies.

Lastly, we also break down our cold-dense reservoir by density and find that it is dominated by gas that is too diffuse to form stars – i.e., only $\sim$0.15% of the cold-dense gas content achieves densities exceeding 1000 cm*-3* (Figures 16 and Figure 17). Emission line ratios, such as $L_{\rm CO(3-2)}/L_{\rm CO(1-0)}$ and $L_{\rm HCN}/L_{\rm CO(1-0)}$ (Gao & Solomon, 2004; Lada et al., 2012; Kauffmann et al., 2017; Onus et al., 2018), are good tracers of the fraction of cold-dense gas in the ultra-dense regime. Hopkins et al. (2013c) shows that, when applied to isolated galaxies, our model – 100% star formation efficiency in self-gravitating (at resolution scale), self-shielded gas at densities exceeding 1000 cm*-3* – is in reasonable agreement with observations. Additional measurements of these emission-line tracers have the potential to constrain the relative fraction of cold ultra-dense to cold-dense gas in interacting galaxies. It is plausible that the turbulent nature of the ISM, especially during and soon after intense periods of star formation, may make the ISM stable against collapse, even in the presence of an abundant gas supply (Alatalo et al., 2015; Smercina et al., 2018). This can potentially explain why molecular gas reservoirs in galaxy pairs are enhanced.

• •

**Hot gas (connection to observations):

Very few observations focused on the hot gas content in interacting galaxies have been performed to date. Henriksen & Cousineau (1999) analyse a sample of 52 galaxy pairs from the Catalog of Paired Galaxies (Karachentsev, 1972), 25 of which are spiral-spiral systems. These authors find that, at a fixed luminosity in the B-band, X-ray luminosity in spiral-spiral pairs is enhanced relative to normal spirals. Similarly, using published publicly-available catalogues (Vorontsov-Velyaminov, 1959; Arp, 1966; Arp & Madore, 1987), Casasola et al. (2004) find that X-ray luminosities from diffuse gas is higher in interacting galaxies than in normal ones. Our simulations predict that close galaxy interactions elevate the hot gas mass component in galaxies significantly (median $\sim$400% – Figure 9, lower-right panel), in line with the aforementioned observations. Smith et al. (2018) use archived Chandra data to measure the diffuse X-ray luminosity, $L_{\rm X}$(gas), from a sample of 49 equal-mass galaxy pairs at various stages of interaction. They find that, for galaxies SFR $>$ 1 M⊙yr-1, $L_{\rm X}$(gas)$/$SFR is not correlated with SFR or interacting stage. These authors do not report ratios of $L_{\rm X}$(gas)$/$SFR versus isolated controls, so a quantitative comparison is not possible. These authors claim that their results suggest that $\sim$2% of the total energy output from supernovae and stellar winds converts into X-ray flux. We remind the reader that our initial conditions do not include a hot gas atmosphere, which might alter our predictions.

4.4 Connection to Other Simulations

Many previous works investigate SFR enhancements and (total) gas depletion due to galaxy mergers using simulations, whereas fewer studies consider the effects of mergers on specific ISM regimes, usually because of limitations of the sub-grid models employed. We mention a few of the most relevant works here. Using two different numerical methods, Di Matteo et al. (2008) find that simulated low-redshift major mergers exhibit only modest SFR enhancements (rarely more than a factor of 5 during the coalescence-phase starbursts, and even less during the pre-coalescence phases). Cox et al. (2008) presents the results from a suite of simulations of low-redshift galaxy mergers to investigate the dependence of starburst strength on mass ratio. They find that the amplitude of the SFR enhancement at coalescence decreases sharply with increasing mass ratio, and even for the equal-mass mergers, the elevation in SFR during the pre-coalescence phase is less than a factor of a few. Various authors have also find that the strength of the SFR enhancement depends on the progenitor galaxies’ gas fractions. For example, Fensch et al. (2017) argue that mergers are less efficient at high redshift, owing to galaxies typically having higher gas fractions than at low redshift.

The conclusions about how mergers enhance the SFR and alter the phase balance of the ISM are unfortunately sensitive to both resolution and the implementation of stellar feedback, thus complicating the interpretation of simulation results. Teyssier et al. (2010, see also ) demonstrate that how well the clumpy structure of the ISM is resolved can have important consequences for how galaxy mergers drive starbursts. They argue that at sufficiently high resolution, interactions lead to enhanced fragmentation into cold clouds and that this enhanced fragmentation (rather than nuclear inflows, as in previous works that employed lower resolution and effective equation of state ISM models, Torrey et al., 2012; Patton et al., 2013; Moreno et al., 2015) is responsible for the SFR enhancement in mergers. Teyssier et al. (2010) find that galaxy-galaxy interactions cause the gas density probability density function (PDF) to shift to higher densities, more so in their higher-resolution simulation. This also occurs in our simulations, as evidenced by the enhancement in cold-dense gas – and especially the even stronger enhancement in the ultra-dense cold gas content (bottom panel of Figure 16). Since specific conclusions may be sensitive to the details of the stellar feedback model used in the simulations (e.g. Hopkins et al., 2013b; Perret et al., 2014); it would be consequently worthwhile to revisit the analysis of Teyssier et al. (2010) using our FIRE-2 simulations in future work.

The above caveats aside, most simulations – both low- and high-resolution and with more or less sophisticated models for stellar feedback – find that mergers have relatively modest effects on the SFR and gas phase structure during most of the pre-coalescence phase, consistent with the results of our work.

5 Conclusions

We present results from a large suite of idealised (non-cosmological) high-resolution (parsec-scale) galaxy merger simulations constructed with the “Feedback In Realistic Environments” (Version 2) model (FIRE-2, Hopkins et al., 2018b), a feedback-regulated galaxy formation framework capable of capturing the multi-phase structure of the interstellar medium (ISM). Our suite consists of 24 simulations, including sets of (nearly) prograde, polar and retrograde mergers, covering a broad range of spatial extents and durations (Figures 3 and 4). Our goal is to investigate how various temperature-density regimes of the ISM fuel and drain one another as galaxy-galaxy interactions induce episodes of enhanced star formation. We consider the following ISM regimes: warm, cool, cold-dense and hot, motivated by the observed ionised, atomic, molecular and hot gas phases in galaxies. See Table 4 and Figure 5 for details.

This papers focuses on the ‘galaxy-pair period’ of merging, defined as the period of time between first and second pericentric passage. Our main results are as follows:

• •

**Interactions enhance star formation in galaxies

**(average enhancement $=$ 30%).

• •

**Interactions do not alter the cool gas content

**(average mass excess $=$ 4%).

• •

**Interactions build up a cold-dense gas reservoir

**(average mass excess $=$ 18%).

• •

Interactions significantly elevate the ultra-dense gas content (average mass excess $=$ 240%).

• •

**Most cold-dense gas is too diffuse to create stars

(only $\sim$0.15% has densities $>$ 1000 cm-3).

The first point is qualitatively in line with observations (Woods et al., 2006; Ellison et al., 2008; Lambas et al., 2012; Ellison et al., 2013; Patton et al., 2013; Scott & Kaviraj, 2014) and predictions by earlier simulations (Mihos & Hernquist, 1996; Di Matteo et al., 2007, 2008; Teyssier et al., 2010; Renaud & Gieles, 2013; Hayward et al., 2014; Moreno et al., 2015). Regarding the second and third points: one might naïvely expect intense star-forming episodes to either consume the cold-dense gas content (and deplete the cool gas as this transforms into cold-dense gas), or expel this gas via feedback produced by new stars. Instead, our feedback-regulated model allows for the recovery of cool gas after brief depletion – and keeps most of the cold-dense gas in a diffuse state (fifth bullet point), thus maintaining the presence of a cold-dense gas reservoir as star formation unfolds, in line with observations of interacting galaxies with molecular gas mass excess (Violino et al., 2018).

Regarding the other two regimes in the ISM, we find:

• •

**Interactions enhance the depletion of warm gas

**(average suppression $=$ 11%).

• •

**Interactions strongly elevate the hot gas content

**(average mass excess $=$ 390%).

In particular, X-ray observations suggest that galaxy encounters augment the hot gas mass budget (Henriksen & Cousineau, 1999; Casasola et al., 2004). However, because our simulations do not include initial hot gas atmospheres, our results are only indicative of partial contributions to the hot gas budget by heating of the galactic ISM.

Future work includes: (1) investigations of how our results depend on the details of orbital configurations and mass ratios; (2) assessment of the spatial of distribution star formation and the structure of the ISM; (3) perform radiative-transfer calculations, coupled with chemical network solvers, to construct mock surveys for better comparisons with existing observations; (4) deeper investigations of the physical mechanisms driving the fate of cold-gas gas, which might either stay moderately-dense, become ultra-dense and be consumed by star formation, or be cycled back into the warm or hot gas budget by stellar feedback; (5) and conduct similar studies at and after coalescence.

Acknowledgements

The computations in this paper were run on the Odyssey cluster supported by the FAS Division of Science, Research Computing Group at Harvard University. Support for JM is provided by the NSF (AST Award Number 1516374), and by the Harvard Institute for Theory and Computation, through their Visiting Scholars Program. The Flatiron Institute is supported by the Simons Foundation. DK is supported by the NSF (AST Award Number 1715101), and by a Cottrell Scholar Award from the Research Corporation for Science Advancement. We thank Jillian Scudder and George Privon for illuminating discussions, Jason Brown for helping us create our galaxy mergers website, and the two anonymous reviewers, whose insightful suggestions improved the quality of this paper. JM thanks Dra. Nicole Cabrera Salazar (Movement Consulting) for co-mentoring a large group of outstanding undergraduate students of colour, including MB. We honour the invaluable labour of the maintenance and clerical staff at our institutions, whose contributions make our scientific discoveries a reality. This research was conducted on Tongva-Gabrielino Indigenous land.