VARIATIONS OF THE ISM COMPACTNESS ACROSS THE MAIN SEQUENCE OF STAR FORMING GALAXIES: OBSERVATIONS AND SIMULATIONS

We use the same suit of hydrodynamical simulations described in Lanz et al. (2014, L14 hereafter), with some additions. Here we will briefly describe those aspects of the simulations that are relevant for our discussion. A hydrodynamical code is combined with a radiative transfer code to obtain mock SEDs at different times for four isolated progenitor galaxies (M0, M1, M2, M3) and the 10 possible binary mergers arising from interactions between these progenitors. These four progenitor galaxies are similar to typical SDSS galaxies, with masses ranging from $6\times {10}^{8}\;{M}_{\odot }$ to $4\times {10}^{10}\;{M}_{\odot }$. Here we add a fifth progenitor galaxy (M4) with a starting stellar mass of $1.25\times {10}^{11}\;{M}_{\odot }$, which is three times more massive than M3. This progenitor is intended to be an analog of a massive spiral that has consumed most of its initial gas mass by forming stars. Its properties are those of the late-type spirals observed in the local universe. For this M4 galaxy, we only simulate a binary interaction with itself. We also include a more massive, gas-rich merger intended to represent a typical submillimeter galaxy (SMG) without a well-defined bulge. Such galaxies are exclusively found at high redshifts ($z\gt 1$). We list the basic properties for our simulated mergers in Table 2.

Table 2.  Galaxy Models for the Simulations

	M0	M1	M2	M3	M4	SMG
M* (${10}^{10}\;{M}_{\odot }$)	0.061	0.38	1.128	4.22	12.5	16.0
${M}_{{\rm{tot}}}$ (${10}^{10}\;{M}_{\odot }$)	5.0	20.0	51.0	116.0	260.0	940.0
${M}_{{\rm{gas}}}$ (${10}^{10}\;{M}_{\odot }$)	0.035	0.14	0.33	0.80	1.8	24.0
${N}_{{\rm{DM}}}$	30000	50000	80000	12000	26400	60000
${N}_{{\rm{gas}}}$	10000	20000	30000	50000	32000	48000

Download table as:  ASCIITypeset image

The hydrodynamical calculations are performed using the Gadget-3 code (Springel 2005), which solves the gas dynamics using smoothed particle hydrodynamics (SPH).⁹ The radiative transfer is calculated using the 3D Monte Carlo dust radiative transfer code Sunrise (Jonsson 2006; Jonsson et al. 2010). The methods for combining the output from the hydrodynamical code with the dust physics for particular types of galaxies are described in Narayanan et al. (2010a, 2010b) and Hayward et al. (2011, 2012).

Unresolved star formation is accounted for by assuming that gas particles with densities above $n\sim 0.1\;{{\rm{cm}}}^{-3}$ form stars according to the volume-density-dependent Schmidt–Kennicutt (S–K) law $\dot{{\rho }_{*}}\propto {n}_{{\rm{gas}}}^{\alpha }$, with $\alpha =1.5$ (Schmidt 1959; Kennicutt 1998). Due to the limited mass resolution of the simulations, individual stars are not created. Instead, equal mass star particles are generated stochastically such that the SFR obtained in this way agrees with the S–K law. For all snapshots in each simulation, we register the following quantities: the instantaneous SFR (${{\rm{SFR}}}_{{\rm{inst}}}$), defined as the sum of the SFRs of the individual gas particles, calculated based on their gas densities and the assumed sub-resolution star formation prescription; the separation between the nuclei of the interacting pair (${d}_{{\rm{BH}}}$); and the mass of all stars formed from the start of the simulation (${m}_{{\rm{N}}}$).

For the radiative transfer, the Milky Way R = 3.1 dust model from Weingartner & Draine (2001) is used, with the Draine & Li (2007) update. Sunrise allows two different treatments of the sub-resolution dust structure. In the "default ISM" (DISM) approach, the dust associated with cold clouds in the Springel & Hernquist (2003) sub-resolution models is ignored, whereas in the "alternate ISM" (AISM) approach, the total dust mass is used. The net effect of this in the radiative transfer is that in the case of the AISM models, photons are propagated through more dust, and the actual amount of dust varies from cell to cell, because the fraction of ISM contained in cold gas depends on the local conditions. The variations in the SED shape due to a particular choice of dust sub-structure are thoroughly discussed in L14. It is important to note that the choice of the sub-structure dust model is perhaps the most significant uncertainty in the radiative transfer calculations. In this paper we use the DISM subresolution model only. The reason for this is that we are mostly interested in merger stages where star formation ins enhanced, and emission is dominated by the emission from H ii regions and PDRs, with little contribution from cold dust.

The emerging SEDs of these systems are computed for ∼100 snapshots spanning a time range of $\sim 6\;{\rm{Gyr}}$. The time between snapshots is 100 Myr, although a finer time resolution (10 or 20 Myr) is used near the peaks of star formation, for example during the coalescence phase of the interactions. SEDs are obtained for seven different viewing angles. In the case of interacting systems, the galaxy pairs were put in a specific parabolic orbit with initial separations increasing with the mass of the larger galaxy. We obtain mock photometry of the simulated systems by convolving (assuming z = 0) the resulting SEDs with the filter response of the following bands: GALEX-FUV, GALEX-NUV, U, B, V (Johnson), J, H, ${K}_{{\rm{S}}}$ (2MASS), Spitzer-IRAC 3.6 $\mu {\rm{m}}$, Spitzer-IRAC 4.5 $\mu {\rm{m}}$, Spitzer-IRAC 5.8 $\mu {\rm{m}}$, Spitzer-IRAC 8 $\mu {\rm{m}}$, IRAS 12 $\mu {\rm{m}}$, Spitzer-MIPS 24 $\mu {\rm{m}}$, IRAS 25 $\mu {\rm{m}}$, IRAS 60 $\mu {\rm{m}}$, Herschel-PACS 70 $\mu {\rm{m}}$, IRAS 100 $\mu {\rm{m}}$, Herschel-PACS 100 $\mu {\rm{m}}$, Spitzer-MIPS 160 $\mu {\rm{m}}$, Herschel-PACS 160 $\mu {\rm{m}}$, Herschel-SPIRE 250 $\mu {\rm{m}}$, Herschel-SPIRE 350 $\mu {\rm{m}}$, and Herschel-SPIRE 500 $\mu {\rm{m}}$.

We describe the statistical details of our fitting algorithm in the Appendix. Here we briefly describe the astrophysical aspects of the models on which the fitter is based. Chiburst is based on the original star-forming galaxy SED models described in Dopita et al. (2005, 2006a, 2006b) and Groves et al. (2008). These models compute the SED of a star-forming galaxy as the combination of three main components (each of them accounting for a normalized continuous SFR over different timescales) that can be scaled to adjust the star formation history (SFH) of each galaxy: (i) a starburst (ionizing stars + H ii regions) population of young stars with a continuous SFR averaged over a period of 10 Myr (SFR₁₀); (ii) a population of stars formed at a constant rate between 10 Myr and 100 Myr ago (SFR₁₀₀); and (iii) a component of very recent ($\lt 1$ Myr) star formation represented by Ultra-Compact H ii regions (UCHIIRs) responsible for dust heating at temperatures of ∼300 K (SFR₁). In order to model more realistic galaxies, we have used Starburst99 (Leitherer et al. 1999) to compute a fourth component (additional to the prescription in Groves et al. 2008) and included it in our models to account for the population of even older (up to 5 Gyr) field stars, which we parametrize according to their total mass (M_*). This addition extends the range of galaxy types that we can study with the models to include systems where older populations significantly contribute to optical and near-infrared (near-IR) wavelengths.

The starburst component comes in two flavors that can be added linearly: a naked H ii region (stars + atomic gas) and an H ii region fully covered by a photon-dominated region (PDR) layer, which is the interface between the ionized gas in the H ii region and the molecular gas birth cloud from which the cluster was formed. We interpret their linear combination as a covering fraction of PDR surrounding the H ii region (${f}_{{\rm{PDR}}}$). As thoroughly described in Groves et al. (2008), the exact shape of the SED for these two components is controlled by two main parameters: the ISM metallicity and a dimensionless quantity related to the dust heating flux: the compactness parameter (${ \mathcal C }$).

The motivation to introduce the concept of compactness arises from the fact that, when dealing with highly complex systems such as starbursts, where different factors play a role in the heating of the dust, it seems rather inaccurate to describe the ISM conditions in terms of a single dust temperature (${T}_{{\rm{dust}}}$), or even in terms of a linear combination of a few values of ${T}_{{\rm{dust}}}$. Instead, we can think in terms of a distribution of dust temperatures $\langle {T}_{{\rm{dust}}}(t)\rangle$ that changes in time as the starburst progresses. We can parametrize the time evolution of $\langle {T}_{{\rm{dust}}}(t)\rangle$ directly rather than assigning single temperature values. Since the latter is determined at any time by the intensity of the stellar radiation field produced by a cluster with luminosity L_*, at radius R in a spherical nebula $\langle {T}_{{\rm{dust}}}(t)\rangle$ is a function of the time-dependent dust heating flux $\langle {L}_{*}(t)\rangle /\langle 4\ast \pi \ast R{(t)}^{2}\rangle$. Hence, models that preserve the run of ${T}_{{\rm{dust}}}$ with time should also preserve the quantity ${L}_{}/{R}^{2}$ averaged over time:

$$\begin{eqnarray}&&\mathrm{log}{ \mathcal C }\propto \displaystyle \frac{\langle {L}_{*}(t)\rangle }{\langle R{(t)}^{2}\rangle }.\end{eqnarray} \tag{ 1 }$$

The compactness parameter is proportional to this ratio, and hence sets the evolution of the ${T}_{{\rm{dust}}}$ distribution with time. By assuming that $\langle L*(t)\rangle$ scales with cluster mass, and by considering the time evolution of the radius and pressure of a mass-loss bubble as in the Castor et al. (1975) approximation, the parameter ${ \mathcal C }$ can be written in terms of the average star cluster mass ${M}_{{\rm{cl}}}$, and the average pressure of the ISM normalized by the Boltzmann' constant(${P}_{0}/k$) in a given system:

$$\begin{eqnarray}&&\mathrm{log}{ \mathcal C }=\displaystyle \frac{3}{5}\mathrm{log}\left(\displaystyle \frac{{M}_{{\rm{cl}}}}{{M}_{\odot }}\right)+\displaystyle \frac{2}{5}\mathrm{log}\left(\displaystyle \frac{{P}_{0}/k}{{{\rm{cm}}}^{-3}\ {\rm{K}}}\right).\end{eqnarray} \tag{ 2 }$$

In Section 6.2 we will discuss the usefulness of the compactness parameter and how measuring this parameter is different from measuring dust temperatures from the SED.

We have fitted the multi-wavelength photometry of the simulated interactions described in Section 3 using Chiburst. This also includes the M4–M4 binary interaction and the gas-rich SMG simulation. As noted, the simulations cover a broad range of galactic masses, interactions stages, luminosities, and morphologies; we therefore expect them to be a fair sample of the conditions found in real galaxies. In total, we have obtained fits for over a thousand mock SEDs, that correspond to different stages along the binary interactions of the galaxies in Table 2. Specifically, the snapshots fitted account for the entire interaction sequence, starting with the initial approach and ending when the merger reaches a post-coalescence passively evolving stage (total times are between 2.5 and 6 Gyr, depending on the total mass of the merger). A higher time resolution is used near the coalescence phase, and in the case of M4M4, for the entire duration of the simulation. In this paper we limit our study of the mock SEDs to a single viewing angle, but further down in this section we briefly explore the effect of viewing angle on the derived parameters.

To illustrate the evolution of the SED along the interaction sequence for a particular simulation, Figure 1 shows an overview of SED fits at different stages. Plotted in the lower left corner is SFR${}_{{\rm{inst}}}$ for the M2–M3 interaction as a function of time, as calculated in the hydrodynamical simulation. SFR${}_{{\rm{inst}}}$ increases during the first passage at around 1 Gyr after the start of the simulation and then increases abruptly during the coalescence phase, at around 3 Gyr. We have indicated three particular times along the interaction sequence that serve as representative stages of the interaction: right after the first passage, at the peak of SFR${}_{{\rm{inst}}}$, and during the phase after coalescence.

Zoom In Zoom Out Reset image size

For the selected snapshots, we show the mock SEDs and the corresponding best Chiburst fits in the insets surrounding the time evolution plot. The most remarkable change in the SED as the interaction evolves is in the relative contribution of far-IR emission (coming from the H ii regions) to the bolometric luminosity of the system, which peaks during coalescence. The contribution of optical and near-IR emission from the oldest field stars remains rather constant along the interaction, although it increases as more stars are formed during the burst. The UV emission (to which both the H ii regions and the stars formed during the 100 Myr prior to a given snapshot contribute) follows the far-IR emission, although it significantly decreases in the relaxation phase after coalescence. A similar behavior is observed in all remaining simulations, with the SEDs of less massive interactions and isolated galaxies evolving much less dramatically than in the strong mergers.

Also shown in the small insets are the probability density functions (PDFs) derived for three of the model parameters: $\mathrm{log}\;{ \mathcal C }$, SFR₁₀, and M_*. The shaded regions in these PDFs correspond to the 1σ (blue) and 90% (red) confidence levels and represent the most general description of the uncertainties involved in the fitting when all possible points of the parameter space are considered. In a few cases, the best-fit values are outside the 90% confidence region, which suggests that they can differ very significantly from the median-likelihood values. Model degeneracies, as well as the marginalization of a multi-dimensional joint probability distribution, are likely responsible for this behavior.

The 1σ confidence regions estimated from the PDFs are small compared to the parameter allowed ranges (typically a third or a fourth of the range), and therefore we are able to infer real differences in the parameters between simulations. One example is the SFR averaged over the last ten million years (${{\rm{SFR}}}_{10}$), estimated from the contribution of the H ii+PDR region to the bolometric SED. Near coalescence, Chiburst gives a 1σ range for this parameter between ∼0.50 and ∼0.80 ${M}_{\odot }\;{{\rm{yr}}}^{-1}$ for the M1–M1, whereas for the M1–M2 interaction the obtained solutions range between 0.30 and 0.65 ${M}_{\odot }\;{{\rm{yr}}}^{-1}$. This relative difference in the instantaneous SFR agrees with the values predicted by the hydrodynamical simulation. Further down, in Section 5.2, we will investigate how the derived PDFs relate to the actual values from the simulations.

The SED fit plots also show how different model parameters control the fluxes in different wavelength bands: the UV is dominated by emission from stars younger than 100 Myr, either those associated with the unobscured, youngest and most massive systems in the H ii regions and parametrized as SFR₁₀ (dashed cyan line), or field A-type stars lasting 10 times longer and parametrized by SFR₁₀₀ (single-dot dashed red line); the optical and near-IR are dominated by the photospheres of the oldest field stars parametrized by M_* (dashed purple line), whereas the mid-IR and far-IR emission comes almost entirely from the H ii regions and the PDRs with bright polycyclic aromatic hydrocarbon (PAH) emission (because no heating of diffuse dust from older stellar populations is included in the models), again parametrized by SFR₁₀ (the dashed light blue line again). In certain systems, there might be a non-negligible contribution to the mid-IR from UCHIIRs, parametrized by SFR₁ (triple-dot dashed yellow line).¹⁰

To illustrate how the method works for mergers of different masses, in Figures 2(a)–(c) we show fits to the mock SEDs of the simulations corresponding to the interaction pairs M0–M1, M1–M1, and M2–M3. In all three cases, the SEDs correspond to the coalescence phase. Typically, more massive interactions have larger bolometric luminosities, and the relative contribution of far-IR (8–1000 $\mu {\rm{m}}$) emission to the total luminosity tends to increase with the total mass of the interacting galaxies, an early indication that more massive interactions have larger sSFRs near the coalescence phase.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The latter happens at first order only, because the mass ratio of the interacting pair affects the relative contribution of dust emission to the SED. The M1–M1 interaction, for example, has a larger far-IR contribution at coalescence than the M1–M2 system, which is more massive. In fact, our results show that the largest fractional far-IR luminosities happen preferentially in interacting pairs where the mass ratio is close to 1. This agrees with results from the numerical simulations of Cox et al. (2008), who show that the strength of the starburst induced by a merger decreases as the mass ratio between the progenitors increases. There is also a trend for the far-IR bump peaking at shorter wavelengths in those systems where the far-IR emission is larger (see also Hayward et al. 2012). We will later interpret this peak shift in terms of the compactness parameter (Equation (2)), and will argue that this parameter and its correlation with the sSFR reveal insightful properties of the ISM physics in interacting systems.

As discussed in L14 and in Hayward & Smith (2015), the viewing angle mostly affects the SED in the UV bands, where the obscuring effect of dust is larger. At these wavelengths, the typical variations in $\lambda {F}_{\lambda }$ are of the order of 0.25–0.5 dex, and reach a maximum immediately after coalescence, when a large amount of UV photons from young massive stars are absorbed by thick layers of dust. To study the effect of viewing angle on our derived parameters, we have fitted a post-coalescence snapshot SED from the M3–M3 simulation viewed from two different, orthogonal angles. The derived parameters PDFs do not change significantly, except perhaps for a difference of 0.1 dex in the derived value of $\mathrm{log}\;{A}_{V}$. More relevant for the present study is the fact that differences in viewing angle do not significantly affect the estimation of the compactness parameter or the SFR, because these parameters are better constrained using the mid-IR and far-IR emission, where dust obscuration is negligible or nonexistent. Moreover, we know that in the post-coalescence phase, when the effect of viewing angle is more important, the infrared indicators significantly overestimate the SFR (Calzetti et al. 2010; Hayward et al. 2014a), and so we interpret the post-coalescence parameters with caution.

In this section we compare the parameter values derived from the fitting of mock SEDs for the original L14 simulations with the true values from the hydrodynamical simulations. As mentioned, in the simulations unresolved star formation is accounted for by assuming that gas particles with densities above certain threshold form stars according to the volume-density-dependent S–K law. The resulting star particles are assigned Starburst99 SEDs with a Kroupa initial mass function (IMF), which are the input for the radiative transfer code. Since Chiburst assumes the same IMF, we assume the true SFR to be the instantaneous SFR from the simulations. In Figure 3(a) we show the correlation between SFR${}_{{\rm{inst}}}$ and the Chiburst-derived value of SFR₁₀.

Chiburst overestimates the instantaneous SFRs, but not dramatically: for the most massive interactions (M2–M3, M3–M3), we obtain SFRs that are an average factor of 2 above the true values, whereas for the less massive systems SFR${}_{{\rm{inst}}}$ and SFR${}_{{\rm{10}}}$ are in good agreement within the 1σ uncertainties derived from the PDFs. The outliers of this correlation, significantly above the majority of points (up to a factor of 100), correspond to the post-coalescence phase of the most massive mergers, where star formation is abruptly quenched. One possible reason for the post-coalescence overestimation of the SFR is heating from stellar populations older than 10 Myr (Groves et al. 2012; Hayward et al. 2014a), which is not accounted for in the models. When dust heating from older stellar populations is included in the SED modeling, the SFR can be recovered much more accurately in the post-coalescence phase (Hayward et al. 2014a). After coalescence, the fraction of recently formed stars decreases dramatically, and this effect becomes more important. Note, however that the correlation spans over a range of 3 orders of magnitude in ${{\rm{SFR}}}_{{\rm{inst}}}$. The SMG simulation is also in good agreement with this correlation, which indicates that Chiburst can determine the true values of SFR within a factor of 2 over a dynamic range covering almost five orders of magnitude.

We also compare the stellar mass in each simulation with our estimation of the M_* parameter from Chiburst. In Figure 3(b) we show this comparison for the original L14 simulations. We slightly underestimate the total stellar mass, obtained by summing up the initial disk and bulge stellar masses, and the mass of stars formed during the simulation. We obtain stellar masses that are a factor of 1.25 or less below the true stellar masses, with better agreement near the coalescence phase. The combined effect of simultaneously overestimating the SFR and underestimating the stellar mass implies that our derived values for the sSFR (${\rm{sSFR}}\equiv {{\rm{SFR}}}_{{\rm{inst}}}/{M}_{*}$) are a factor of 2.5 (or 0.4 dex) above the true sSFR values.

These is a systematic effect and therefore does not affect the conclusions we make below regarding relative correlations between model parameters. In general, the uncertainty in the determination of the sSFR is method-dependent. By collecting published values of the sSFR as a function of stellar mass and redshift using various methods, Behroozi et al. (2013) have estimated the uncertainty in sSFR due to the use of different techniques. They show that such uncertainty varies from 0.3 dex to 0.4 dex for stellar masses between ${10}^{10.5}\;{M}_{\odot }$ and ${10}^{9.5}\;{M}_{\odot }$. Our derived values for sSFR are therefore within the uncertainty associated with the use of a particular method.

Using Chiburst, we have fitted the multi-wavelength SEDs of the galaxies listed in Table 1. For most of these systems, the wavelength coverage is similar to that of the mock photometry presented in the last section, although some of the bands are missing for specific galaxies. The intermediate-z LIRGs appear heavily obscured, and therefore no GALEX detections are available. For these obscured systems, however, the rest-frame UV emission that for local systems falls within the GALEX bands is redshifted into the optical, making the lack of GALEX data a less severe problem for our purposes.

In Figure 4 we show example fits for representative galaxies in each of our three groups: local, average luminosity interactions, local LIRGs, and the intermediate-z LIRGs from Magdis et al. (2014). In most cases, residuals of the fit are within ±0.2 dex over the entire wavelength range. One exception is BOOTES1 which has a flat, featureless mid-IR spectrum. Our hydrodynamical simulations, as well observational evidence from IRS observations of embedded active nuclei (Imanishi et al. 2007), suggest that such flat spectrum is due to a significant contribution from mid-IR AGN thermal emission to the bolometric luminosity (Snyder et al. 2013). None of the other systems appears to have a mid-IR spectrum compatible with a significant contribution from an AGN. As we have mentioned before, minor AGN contributions can slightly affect the SFR estimates, but not dramatically.

Zoom In Zoom Out Reset image size

Table 3 lists the Chiburst-derived parameters for the observed systems. The range of stellar masses and SFRs that we obtain for the local interactions fall within the ranges spanned by these same parameters in the simulated interactions (see Section 5.5 below), which allows a comparison between the observed systems and corresponding simulated stages of the interactions. Not surprisingly, observed systems that are in later stages of the interaction (i.e., strongly interacting, morphologically disturbed systems such as NGC 3034 and NGC 3690) have SFRs and stellar masses similar to those of the most massive simulated interactions near the coalescence phase, whereas observed systems in early phases of the interactions usually have lower SFRs for similar masses, and have similar derived parameters as the early stages of the simulated interactions. Such is the case for M51a and NGC 3031. To illustrate this association between observed and simulated systems, in Figure 5 we show a comparison of the derived probability distributions for observed and simulated systems.

Zoom In Zoom Out Reset image size

Table 3.  Fitting Results

Galaxy	Type	$\mathrm{log}$ SFR10	$\mathrm{log}$ SFR100	$\mathrm{log}$ M*	$\mathrm{log}$ SFR1	${F}_{{\rm{PDR}}}$	$\mathrm{log}\;{ \mathcal C }$	$\mathrm{log}{A}_{V}$
		(${M}_{\odot }\;{{\rm{yr}}}^{-1}$)	(${M}_{\odot }\;{{\rm{yr}}}^{-1}$)	(${M}_{\odot }$)	(${M}_{\odot }\;{{\rm{yr}}}^{-1}$)			(mag)
NGC 2976	Local merger	−0.610	−1.143	9.188	−1.714	0.67	4.88	−0.357
NGC 3031	Local merger	0.008	0.113	10.611	−0.936	0.75	4.36	−0.004
NGC 3034	Local merger	1.007	−1.124	10.154	1.819	0.84	6.04	0.665
NGC 3185	Local merger	−0.120	−0.761	10.073	−0.733	0.81	5.02	−0.057
NGC 3187	Local merger	−0.005	−0.160	9.593	−0.921	0.71	4.72	−0.183
NGC 3190	Local merger	0.263	−1.400	10.652	−0.617	0.92	4.82	0.367
NGC 3395/3396	Local merger	0.848	0.962	10.144	0.154	0.62	5.55	−0.463
NGC 3424	Local merger	0.650	−0.835	10.289	−0.000	0.96	5.54	0.275
NGC 3430	Local merger	0.615	0.355	10.285	−0.749	0.70	5.27	−0.123
NGC 3448	Local merger	0.458	0.232	9.935	−1.035	0.79	5.32	−0.028
UGC 6016	Local merger	−1.414	−0.438	8.140	−2.375	0.82	4.79	−0.722
NGC 3690/IC 694	Local merger	2.213	0.622	10.718	1.953	0.95	6.36	0.374
NGC 3786	Local merger	0.382	−0.456	10.331	0.274	0.81	5.17	−0.032
NGC 3788	Local merger	0.271	−0.254	10.239	−0.401	0.81	4.72	−0.263
NGC 4038/4039	Local merger	1.405	0.939	10.888	0.476	0.81	5.49	−0.390
NGC 4618	Local merger	−0.556	0.150	9.300	−1.487	0.83	4.90	−0.768
NGC 4625	Local merger	−0.985	−0.798	8.891	−2.538	0.78	4.98	−0.305
NGC 4647	Local merger	0.320	0.323	10.175	−0.793	0.87	5.28	0.383
M51a	Local merger	0.869	0.429	10.571	−0.010	0.76	5.14	−0.213
M51b	Local merger	−0.081	−1.031	10.213	−0.619	0.85	5.10	0.263
NGC 5394	Local merger	1.140	0.346	10.428	1.206	0.90	5.84	0.472
NGC 5395	Local merger	1.204	0.513	11.089	0.298	0.84	4.94	−0.058
M101	Local merger	0.741	0.847	10.428	−0.107	0.53	4.94	−0.853
NGC 5474	Local merger	−1.129	−0.300	8.976	−1.996	0.81	4.49	−0.843
NGC 2623	Local LIRG	1.681	0.385	10.352	0.990	0.89	5.70	0.576
UGC 4881	Local LIRG	1.957	0.361	11.017	0.979	0.95	5.55	0.463
VV 283	Local LIRG	1.813	0.305	10.702	0.820	0.96	5.77	0.533
Mrk 273	Local LIRG	2.236	0.559	10.697	1.991	0.96	6.48	0.534
VV 705	Local LIRG	2.075	0.516	10.861	2.297	0.91	6.13	0.428
NGC 6090	Local LIRG	1.827	0.240	10.557	1.720	0.86	6.08	0.324
ELAISS	Interm. z LIRG	2.085	2.319	11.081	1.328	0.80	4.70	1.016
CDFS2	Interm. z LIRG	2.149	2.278	10.876	1.221	0.85	4.88	1.123
CDFS1	Interm. z LIRG	2.355	2.201	11.160	1.719	0.82	4.81	0.931
SWIRE4	Interm. z LIRG	2.250	1.846	10.777	1.360	0.84	5.53	0.756
SWIRE5	Interm. z LIRG	2.609	1.720	11.590	1.709	0.85	5.21	0.608
SWIRE2	Interm. z LIRG	2.362	1.703	11.417	1.633	0.82	4.36	1.018
SWIRE7	Interm. z LIRG	2.879	1.845	11.142	2.146	0.91	5.49	0.464
BOOTES2	Interm. z LIRG	2.237	1.599	10.893	1.335	0.84	5.28	0.626
BOOTES1	Interm. z LIRG	3.187	1.671	11.262	3.755	0.84	6.09	0.868

Download table as:  ASCIITypeset image

All of the local and intermediate-z LIRGs are outside the range of sSFRs covered by the binary interactions for the same stellar masses (Figure 7) (i.e., they are outliers of the MS), but they appear to be somewhere in between the region of the parameter space covered by these progenitors and the region covered by the gas-rich SMG simulation. Their compactness values are also low relative to those of local interactions with similar sSFRs. This is relevant in the context of the current discussion regarding whether the enhanced sSFRs of intermediate- and high-z galaxies is due to higher gas fractions in galaxies earlier in cosmic history (e.g., Daddi et al. 2010a; Tacconi et al. 2013), or to higher efficiencies in the conversion of gas into stars, likely driven by mergers (e.g., Frayer et al. 2008; Daddi et al. 2009). Here we will quantify the effect of both interactions and gas content in rising a particular galaxy above the MS.

As originally defined in Groves et al. (2008) and Section 4.2, the compactness parameter (${ \mathcal C }$) is related to the internal conditions of the ISM, specifically to the ISM pressure and the radiation field strength of individual H ii regions. Despite the fact that the hydrodynamical simulations considered here do not resolve the small-scale structure of the ISM and simplify its physics, the self-similar structure of the ISM implies that an analogous interplay between radiation strength and pressure takes place at the scales resolved by our simulations. Also, the sub-resolution models that describe the ISM physics at smaller scales are chosen within the simulations according to the resolved properties of the ISM for a given particle. We therefore expect ${ \mathcal C }$ to be a meaningful description of the large scale geometry of the dust with respect to the stars, and hence a reasonable probe of the global galaxy compactness.

Here we characterize observed and simulated systems according to their ${ \mathcal C }$ parameter in order to study whether the geometry and heating conditions of dust in galaxies change as a function of particular stages along the interaction sequence. Figure 6(a) shows how the derived ${ \mathcal C }$ values for observed and simulated galaxies correlate with the derived sSFRs. The plotted values correspond to the maximum of the marginalized PDF for each parameter. Simulations of the isolated and interacting progenitor disk galaxies (M0 to M4) at different snapshots are marked with small circles, both before (filled circles) and after (empty circles) coalescence, whereas the SMG simulation is indicated by the large yellow circles. Following the classification of Table 1, we have separated the observed systems into local interactions (red triangles), local LIRGs (blue diamonds), and intermediate-z LIRGs (green squares).

Prior to coalescence¹¹ , all simulated systems, regardless of total mass or gas fraction, appear to fall along a single power law of the form $\mathrm{log}{ \mathcal C }\propto \mathrm{log}{({\rm{sSFR}})}^{\alpha }$. This correlation between compactness and sSFR provides a physical framework that can be used to explain the correlation found by Magnelli et al. (2014) for $z\lt 2$ galaxies between sSFR and ${T}_{{\rm{dust}}}$. In fact, Magnelli et al. (2014) argue that at a given redshift, the ${T}_{{\rm{dust}}}$–sSFR correlation implies that galaxies situated above the MS are dominated by star-forming regions exposed to high radiation fields. Elaborating along similar lines, Magdis et al. (2012) argue that the mean radiation field $\langle U\rangle$ is the only parameter controlling the shape of the SED of high redshift galaxies. The parametrization presented here implies that, in fact, at scales larger than a few hundred parsecs, dust has a more compact geometry in outliers relative to MS galaxies. The latter can be an effect of a higher radiation field due to more massive star clusters, an increased pressure field due to the merger, or both.

In Figure 6(a), observed local interactions follow a trend that is similar in slope and normalization (see Table 4) to that of the simulations, an indication that the latter constitute a reasonable control group on which we can test our assumptions on the star-forming properties of galaxies. This agreement between models and observations is even more remarkable when the morphology of the observed systems is considered. Isolated disks and early interactions such as NGC 3031 and NGC 5474 typically have low $\mathrm{log}\;{ \mathcal C }$ values, but compactness gradually increases as we go toward late stage mergers, with the highest $\mathrm{log}\;{ \mathcal C }$ values derived for strongly disturbed systems such as NGC 3034. That local interactions do not exactly overlap with the simulations in this diagram can be attributed to the fact that the latter, although designed to reproduce a broad range of galaxy types, are nevertheless imperfect constructions, simplify the ISM physics, and span a relatively small parameter space in terms of total gas masses and gas fractions.

More interesting is the fact that both local and intermediate-z LIRGs follow the same correlation, although the normalization constant in each case is different with respect to the local systems. Furthermore, both simulations and observations are consistent with late stage interactions having elevated values of compactness, which implies that the interaction stage of an unresolved galaxy can be inferred from the properties of its SED. This is consistent with L14, where the authors claim that the SED shape can be used to discriminate between coalescence and non-coalescence systems. The tight correlation of Figure 6 indicates that the SED shape alone is able to discriminate other stages during the interaction.

Zoom In Zoom Out Reset image size

The normalization of the MS declines from z = 2.5 to z = 0.0 (Daddi et al. 2007; Noeske et al. 2007; Magdis et al. 2010; Elbaz et al. 2011; Magnelli et al. 2014), which implies that the vertical dashed line in Figure 6(a) moves toward larger values of sSFR at higher redshifts. The relative offset of the intermediate-z LIRGs in the figure could be due to this effect. In order to compare all systems with respect to the location of the z = 0 MS, in Figure 6(b) we have shifted the intermediate-z LIRGS to match the sSFR that they would have if they were in the Local universe, according to the parametrization of Elbaz et al. (2011) (their Equation (13)). We then fit three different power laws ($\mathrm{log}{ \mathcal C }=a\times \mathrm{log}{\rm{sSFR}}+b$) to the resulting correlation, each corresponding to one of the following subsets: the simulations prior to coalescence (filled circles), the observed local interactions (red triangles)¹² , and all the LIRGs, both local and at intermediate z (blue diamonds and green squares, plus NGC 3690, the uppermost red triangle in the figure). In Table 4 we list the parameters obtained for the power laws.

In Figure 7 we plot the simulation snapshots and the observed systems in the SFR-M_* plane, the plane on which the MS correlation is usually defined. The same z = 0 MS as in Figure 6(a) is indicated by the dashed line and the shaded area, and we have added dotted–dashed lines to indicate a distance of 0.9 dex from the MS. We observe a gradient in compactness as a function of distance to the MS in the simulations, with compactness increasing from the bottom to the top envelope of the MS. Practically all simulations involving the M0 to M4 progenitors fall within the 0.9 dex limits of the MS, at all times¹³ , but only a fraction of them fall within the 0.3 dex limits. The gas-rich SMG simulation, on the other hand, is a clear outlier of the correlation, sitting clearly above the 0.9 dex limit. The local interactions overlap with the simulations on this plane, whereas both the local and intermediate-z LIRGS are significantly above this correlation.

Zoom In Zoom Out Reset image size

It is important to note here that stages close to coalescence are over-represented in the sample of snapshots shown here, except for M4–M4, because the time resolution used to generate the SEDs is larger near the peak of infrared luminosity. In other words, in the real world galaxies spend only a small fraction of their lifetimes as strongly disturbed mergers (at high compactness stages), and a similar diagram weighted with a realistic merger fraction would look more crowded in the low compactness end. This means that the scatter shown if Figure 7 is not representative of the scatter seen in the observed MS. Nevertheless, the trend shown here is relevant in the context of the evolution of dust geometry with respect to the stars along the merger sequence, and the effect that this evolution has in shaping the MS. Undisturbed disk galaxies are located somewhere near the lower envelope of the correlation, where the majority of real galaxies lie along the MS, and as they evolve toward coalescence they move away from the MS, and to higher values of compactness. In our hydrodynamical models, this increase in compactness translates into a smooth evolution of the large scale geometry of the dust with respect to the stars, driven by higher pressures and/or stronger radiation fields, as galaxy pairs approach coalescence.

For a binary merger to have a very high compactness ($\mathrm{log}\;{ \mathcal C }\gt 6.0$), its mass ratio has to be very close to 1. Only the most massive mergers reach a very compact stage, and this only during the coalescence phase. However, not even the massive M4–M4 merger can produce the observed outlier sequence of LIRGs. Looking at the correlation between sSFR and compactness for the simulations we observe that as mergers become more massive for a given gas fraction (${f}_{{\rm{gas}}}\quad \sim$ 0.2–0.3), they move upwards along the same diagonal correlation in Figure 7. The gas-rich SMG simulation (${f}_{{\rm{gas}}}\quad \sim$ 0.2–0.6), on the other hand, is the only one of all our simulations that appears significantly above the MS, indicating that only a significant increase in gas fraction can produce the population of LIRG outliers. This is true even if such gas-rich simulations are not yet in the coalescence phase.

We now estimate the parameter ${\rm{IR8}}={L}_{\mathrm{IR}}/{L}_{8\mu {\rm{m}}}$ defined in Elbaz et al. (2011) for our set of simulated and observed interactions. Elbaz et al. (2011) show evidence that the vast majority of IR8 values for galaxies in the GOODS-Herschel field follow a Gaussian distribution with median ${\rm{IR8}}\sim 4$. They claim that such distribution defines an infrared MS. The outliers of this infrared MS ($\lt 20$%) form a tail toward larger values of IR8 and typically have IR8 ∼ 10. They also show that these outliers are systems with compact projected star formation densities. Here we investigate whether these features of the infrared MS are also present in our simulated and observed systems.

We measure IR8 from the best-fitting SEDs obtained with Chiburst. We integrated the rest-frame SEDs between 8 and 1000 $\mu {\rm{m}}$ to estimate the value of L_IR for each system. For observed and simulated systems we then used the (measured or mock) flux density at 8 $\mu {\rm{m}}$ as an estimate of L₈. For those observed systems for which we did not have measured rest-frame 8 $\mu {\rm{m}}$ flux densities available, we integrated the SEDs convolved with the IRAC 8 $\mu {\rm{m}}$ filter response, and used the resulting value as L₈. We took care in performing the same process for the intermediate-z LIRGs, using their rest-frame SEDs.

Figure 8 shows IR8 as a function of the 8 $\mu {\rm{m}}$ luminosity for all simulated and observed systems, color-coded according to their compactness. For the simulated systems (except the SMG simulation), the mean of the distribution is IR8 $=\quad 2$ (solid line in Figure 8), and the standard deviation is 0.12 dex. The majority of simulated snapshots are within 1σ of the mean value (indicated by the dashed lines), and outliers include systems with both higher and lower values of IR8. Simulated outliers with higher values of IR8 correspond to low luminosity systems (${L}_{8}\lt {10}^{8}\;{L}_{\odot }$, isolated M0, isolated M1, and M0-M0 interaction). This resembles recent results in a survey of local galaxies by Cook et al. (2014), who find a departure from the infrared MS toward higher IR8 values in low luminosity systems and attribute this to the suppression or destruction of PAH molecules. Chiburst also finds smaller PDR covering fractions(${f}_{{\rm{PDR}}}\lt 0.4$) in these systems. In this case, however, rather than PAH destruction, the effect is due to a smaller filling factor of PDRs. No clear indication of a high-IR8 tail can be inferred from these simulated interactions. Nevertheless, observed outliers of the MS in Figure 7 also appear as outliers of the IR8 distribution.

Zoom In Zoom Out Reset image size

Objects with high compactness are located at the luminous end of the distribution, and their IR8 values are within the 1σ boundaries. About half of the observed local interactions (triangles in Figure 8) have higher values of the IR8 parameter (mean is 2.9 and standard deviation of 0.13 dex). The SMG simulation has IR8 values well above the distribution delimited by the dashed lines, with a mean value of IR8 = 5.5. Unlike the low luminosity systems, however, the SMG simulation snapshot SEDs do not show attenuated PAH emission, and therefore in this case we associate the increase in IR8 to an increase in ${L}_{{\rm{IR}}}$, which in turn might be associated with the relatively higher fraction of ISM mass in this gas-rich simulation.

The local and intermediate-z LIRGs also have IR8 values more than 1σ above the local interactions and the simulations, but still below the values for the SMG simulation, with a mean IR8 value of 4.1 and a small standard deviation of 0.05. With the limited information we have on the gas content of these systems, it is difficult to assess whether the effect is due to an increases ISM mass, or to suppression/destruction of PAHs. These systems, which are clear outliers of the MS, are also far off the IR8 distribution. Our simulations show that an elevated gas fraction can account for the oddities observed in disks, both in the regular MS and the infrared MS.

The mean value of IR8 for the simulated systems in Figure 8 is a factor of 2 lower than the value reported in Elbaz et al. (2011), for both their local sample of galaxies and the high redshift sample. A number of reasons could explain this disagreement. First, because of the heterogeneous nature of their local sample, both their L_IR and L₈ were estimated based on the available photometry, and no homogeneous method was used. In particular, for high-z galaxies, they were relying on SEDs that were usually incomplete in the long-wavelength end, since this part of the spectrum had been redshifted beyond the Herschel bands. Furthermore, they use a different set of SED templates (namely the Chary & Elbaz 2001, templates) to perform the fits to the galaxy photometry. These templates have a similar treatment of the dust particles as our models, but differ in the treatment of the emission from PAHs. Whereas Chiburst uses an empirical PAH spectrum with optical properties similar to those of coronene, and fitted to two particular interacting galaxies, the Chary & Elbaz (2001) templates use a mixture of six different planar PAH, and their emissivities are explicitly calculated.

Convergence of the mock SEDs with respect to the SPH resolution of the Gadget-3 simulations is necessary in our analysis. In order to demonstrate that our results do not strongly depend on resolution effects, we have performed an additional version of the simulated M3–M3 merger with eight times more particles, which implies that the particle mass is eight times smaller. We have also used softening lengths that are half those of the standard runs. In Figure 9(a) we show the resulting difference in the SFH due to the increase in resolution. In general, the behavior is very similar for both resolutions, and we only observe minor differences along the duration of the simulation. The largest deviations occur between the first and the second encounter (between 1 and 2 Gyr from the start of the simulation), and at the peak of star formation during coalescence (at ∼2.5 Gyr). At these stages, the SFRs can differ by a factor of about 2 between the low and high resolution simulations. These differences do not significantly affect our result in Figure 6, because changes in sSFR due to resolution effects are also accompanied by consistent changes in compactness, as can be seen in Figure 9(b), where we show that snapshots from the high resolution simulation overlap with those from the standard run in the $\mathrm{log}$ sSFR-$\mathrm{log}{ \mathcal C }$ space. This implies that the correlation derived between sSFR and compactness holds at different SPH resolutions, within the uncertainties indicated. In support of this result, Hayward et al. (2014b) show that in the case of mergers up to coalescence (although possibly not after coalescence), the physics derived using SPH methods are robust to the inaccuracies inherent to SPH itself.

Zoom In Zoom Out Reset image size

Figure 9(b) is also useful to visualize how the $\mathrm{log}$ sSFR-$\mathrm{log}{ \mathcal C }$ correlation relates to the time evolution of the merger: at initial approach, the simulation starts off with low sSFR and low compactness. As the merger approaches coalescence, there is a boost of the star formation activity, the radiation field intensity in a given volume of gas increases, and gas compression produces an increase in pressure. Therefore, both compactness and sSFR increase until they reach a maximum at coalescence, and in the process, a particular galaxy moves along the correlation in Figure 6. For the gas fractions considered in the simulations, that range from ${f}_{{\rm{gas}}}\sim 0.16$ to ${f}_{{\rm{gas}}}\sim 0.36$ (see Table 2), this evolution occurs along the same power law and the same normalization. Only the SMG simulation, which has a relatively larger gas fraction (${f}_{{\rm{gas}}}\sim 0.6$) appears shifted toward larger values of sSFR in this diagram.

In intermediate and high redshift surveys of star-forming galaxies, the observed scatter across the MS as well as the existence of a population of outliers pose challenging questions regarding the processes controlling star formation inside galaxies. A paradigm has emerged according to which galaxies that belong to the MS have secular star formation histories maintained by the accretion of gas from the intergalactic medium, whereas outliers of the MS form stars in timescales significantly shorter than their dynamical timescales. Based on morphological studies of a large sample of galaxies with z up to 2.5, Wuyts et al. (2011) find that a majority of galaxies in the MS show a disk-like morphology, which supports this scenario. Several studies, however, have found a significant population of mergers (up to 20%) on the MS (Kartaltepe et al. 2012; Hung et al. 2013). In this scenario, the intrinsic scatter of the MS arises from the stochasticity in the gas accretion history¹⁴ (Kelson 2014; Whitaker et al. 2014).

In this context, some relevant questions are: Do the heating conditions of the dust change across this correlation? Are mergers in a stage of coalescence responsible for the population of outliers of the MS? Here we have shown that gas content is not the only possible source of scatter in the MS. Instead, our results (Figures 6 and 7) reveal that, for moderate gas fractions, the large-scale geometry of the dust evolves smoothly across the MS as interacting pairs approach coalescence, from less to more compact, and that this process also moves galaxies away from the MS, toward higher star formation efficiencies. This change in compactness is manifested in the observed sSFR-${T}_{{\rm{dust}}}$ correlation. Evidence for variations in the ISM physics across the MS for a given redshift has been collected by several authors (Daddi et al. 2010b; Saintonge et al. 2011; Magnelli et al. 2014), although some studies have suggested that the distance from a particular galaxy to the locus of the MS depends only on ${f}_{{\rm{gas}}}={M}_{{\rm{gas}}}/{M}_{*}$ (Magdis et al. 2012).

Mergers therefore result in mobility of galaxies across the MS, but because of the short duration of the coalescence phase, the sharp increase of star formation during coalescence has only a minor impact on the scatter of the observed MS. Also, extreme outliers of the MS such as LIRGs fall significantly above the sequence, regardless of the mass of the interacting galaxies. This implies that other mechanisms should be invoked in order to explain their extremely high sSFR. One possibility is that they have significantly higher gas fractions (${f}_{{\rm{gas}}}$), as our SMG simulation demonstrates. The fact that LIRGs are also outliers of the infrared MS suggests that the change in physical conditions associated with the jump to more efficient star formation also makes the ratio of total infrared luminosity to PAH emission weaker. Our results indicate the this is not always associated with suppression of the PAH emission.

A fundamental difference between compactness and ${T}_{{\rm{dust}}}$ relies on the fact that, whereas ${ \mathcal C }$ is a parametrization of the compactness of star formation—that is, a measurement of the geometry of the star-forming ISM with respect to the (ionizing) stars—${T}_{{\rm{dust}}}$ is a result of this geometry (and of other factors affecting the dust heating in galaxies, as discussed in Groves et al. 2012).

In the Groves et al. (2008) models, compactness parametrizes the time evolution of ${T}_{{\rm{dust}}}$ as a function of the ISM radiation field and pressure. In the context of this paper, changes in compactness are related to variations in the distribution of dust temperatures as a function of sSFR. Even though the two quantities ${ \mathcal C }$ and ${T}_{{\rm{dust}}}$ are closely related, they do not necessarily provide the same information. For example, ${T}_{{\rm{dust}}}$ can be influenced by heating sources other than the young massive stars in clusters (see, for example Groves et al. 2012) and by variation in the radiation field (i.e., young versus old bursts), whereas compactness parametrizes purely the geometry of the dust associated with star formation. As an example of the differences between ${ \mathcal C }$ and ${T}_{{\rm{dust}}}$, we note that in the high sSFR end ($\mathrm{log}{\rm{sSFR}}\gt -10\;{{\rm{yr}}}^{-1}$), the compactness parameter is significantly more sensitive to changes in the sSFR than the ${T}_{{\rm{dust}}}$, as shown in Figure 10. Stages with such large values of sSFR correspond in our simulations to the coalescence phase of the most massive interactions, when the star formation activity peaks (and presumably, when the ISM conditions are changing more rapidly and the bolometric luminosity has the highest contribution from star-forming ISM). By using compactness, we can therefore characterize the most intense epochs of star formation as a function of the ISM properties, something that it is not possible based on ${T}_{{\rm{dust}}}$ measurements only.

Zoom In Zoom Out Reset image size

The large scale ($\gt 200$ pc) geometrical structure of the star-forming ISM is reflected in the compactness values that we have estimated for the simulated interactions. At these scales, a larger value of ${ \mathcal C }$ implies a higher projected star formation density, because both the average radiation field and the gas density increase in areas of concentrated star formation. Therefore, the properties of the correlation found between ${ \mathcal C }$ and sSFR (Figure 6 and Table 4) contain information regarding the large scale geometry of star-forming gas as a function of the total star formation taking place in these systems. We use the slope of the correlation to quantify this relation.

We start with the fundamental definition of compactness in Equation (1), and assume that the proportionality holds for all times. The amount of luminosity being produced in a patch of ISM with radius R scales linearly with the sSFR within that particular patch, and hence:

$$\begin{eqnarray}&&\mathrm{log}{ \mathcal C }=\mathrm{log}{\rm{sSFR}}-2\mathrm{log}R+\mathrm{constant}.\end{eqnarray} \tag{ 3 }$$

Now, we have already estimated the slope of the correlation to be ∼1.6. Therefore we have:

$$\begin{eqnarray}&&\mathrm{log}{ \mathcal C }=1.6\;\mathrm{log}{\rm{sSFR}}+\mathrm{constant}.\end{eqnarray} \tag{ 4 }$$

If we compare Equations (3) and (4), then we can infer that if ${{\rm{sSFR}}}^{-0.3}\propto R$, then we recover the measured slope of 1.6 for the $\mathrm{log}$ sSFR-$\mathrm{log}{ \mathcal C }$ correlation. A possible interpretation for this is that systems with high sSFR have smaller physical sizes. For a fixed gas mass, the volumetric S–K relation implies a similar effect, i.e., that the sSFR should be inversely proportional to the volume of the region in consideration. Qualitatively, this is also consistent with observations of local LIRGs. For example, approximately half of the luminosity of Arp 220 is emitted from within its central 100 pc, which means that this object with very high sSFR is also very compact.

The observed segregation between normal galaxies and LIRGs in (Figure 7), resembles the bimodality found by some authors in the S–K law, that relates the projected molecular gas density with the projected star formation density. For example, Daddi et al. (2010b) find that disk galaxies and starbursts occupy different regions in the gas mass versus SFR plane. They argue that the existence of these two regimes is related to different fractions of molecular gas in the ISM, and that only accounting for the dynamical timescales of galaxies recovers a single, universal star formation law. A different approach is adopted in Krumholz et al. (2012), where the authors claim that the observed variations are due to geometrical effects, and that all objects, from small molecular clouds to LIRGs, follow a universal volumetric star formation controlled by the local free fall times.

In principle, we could explore the implications of our results (in particular the measured slopes in Table 4) under the assumption that the Krumholz et al. (2012) interpretation holds, and attempt to constrain the physics of the small scale star-forming ISM, such as pressure and average cluster mass. However, this would only be possible if the simulations used are fine enough to resolve such smaller scales ($\lt 100$ pc). Given our resolution limitations, we leave this endeavor for an upcoming paper.

The offset between LIRGs and local interactions in Figures 6 and 7 implies one or several of the following possibilities: they have different gas fractions (${f}_{{\rm{gas}}}$), they have different star formation efficiencies, or they have different conditions in pressure and radiation field that result in a different compactness for a given bolometric luminosity. Enhanced gas fractions can rise galaxies significantly above the MS, even in early quiescent, disk-like stages, as demonstrated by our results on the SMG simulation. We examine here whether an increased gas content in the observed LIRGs is responsible for their increased sSFR.

As we have mentioned, the normalization of the MS declines from z = 2.5 to z = 0.0, and the reasons for this decline are not entirely understood. Recent ALMA observations suggests that all galaxies at all cosmic times convert gas into stars with the same efficiency, but that galaxies at intermediate and high z have relatively larger gas reservoirs (Scoville et al. 2015). On the other hand, using a multi-wavelength stacking analysis of the SEDs of galaxies in the COSMOS field, Béthermin et al. (2015) show that even if the bulk of star formation up to $z\sim 4$ is dominated by secular processes, star formation efficiencies in strong starburst differ from those in MS galaxies.

Using the molecular gas mass estimates for the intermediate-z LIRGs in Magdis et al. (2014), which are based on SED fitting and an assumed gas-to-dust ratio, we find that these galaxies have typical ${f}_{{\rm{gas}}}$ of about 0.7–0.8 (and as high as 1.1 for SWIRE7, which brings it close to a typical SMG). This is about a factor of 3 higher than the value of ${f}_{{\rm{gas}}}$ in the simulations and local mergers, including the luminous stage 4 mergers (Sanders et al. 1991). This suggests that, in agreement with our SMG simulation, gas content plays an important role in creating the population of outlier LIRGs in the MS.

Finally, we want to emphasize that regardless of the normalization factor, the slope of the $\mathrm{log}$ sSFR–$\mathrm{log}\;{ \mathcal C }$ is very similar for simulations, normal galaxies, and LIRGs. This suggests that LIRGs might be gas-rich analogs of local galaxies, perhaps undergoing similar merging events. Confirming this would require morphological information that we currently lack for the intermediate-z LIRGs, but it is worth investigating it. Nevertheless, the large range of compactness values that these luminous systems span indicates that not all LIRGs have compact star formation, and that at least some of them should have disk-like morphologies. In fact, Magdis et al. (2014) show that the $z\sim 0.3$ LIRGs have ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{{\rm{IR}}}$ ratios similar to those of local disks. Regarding the morphology of LIRGs and ULIRGs, recent morphological and kinematical studies have revealed an increasing number of non-interacting, disk-dominated ULIRG-like systems with look-back time (e.g., Förster Schreiber et al. 2009; Kartaltepe et al. 2012; Wang et al. 2012). K. L. Larson et al. (2016, in preparation) also find evidence for non-interacting morphology well above the MS in a sample of nearby galaxies.

The hydrodynamical simulations used in this paper employ a simplified model for the ISM. First, with respect to other recent simulations, they assume a relatively low density threshold for the onset of star formation. Also, they employ a pressurized equation of state, which might have the effect of artificially smoothing the ISM, thus affecting the emerging far-IR emission.¹⁵ Other groups have used higher density thresholds and more sophisticated models for the ISM. Such simulations certainly feature more complicated small-scale ISM structure and might be able to more accurately reproduce its conditions at scales smaller that the resolution limit of the simulations presented here, which is of the order of a few hundred parsecs. This limitation in ISM physics notwithstanding, the usefulness of our approach lies on our ability to parametrize the large scale (>200–300 pc) structure of the ISM and to study how this relates to the star formation properties at even larger scales (e.g., the global compactness of galaxies as they evolve across the MS and the global SFR). Furthermore, we note that even though we do not resolve the small-scale structure of the ISM, the trends in star formation properties that we derive in Sections 5.4 and 5.5 for the simulations are in good agreement with what is observed in real interacting galaxies. This direct comparison of simulated and real galaxies through their parametrized SEDs therefore implies that the small scale physics that control star formation in galaxies are imprinted in the ISM structure at larger scales, and that it is currently possible to study such structure.

Suppose that you have obtained photometry of a galaxy in different bands, with certain observational uncertainties associated. Bayes Theorem states that, given those observations, the probability $P(M| D)$ of an SED model M being a true representation of the observed galaxy SED data D, is proportional to the product of the likelihood that your data-points can be obtained from your model ($P(D| M)$) times a prior distribution P(M) that contains independent evidence of certain model parameters having certain values. $P(M| D)$ is what we call the posterior PDF and is the solution we are after. The likelihood $P(D| M)$, or the probability of the data given the model, can be obtained from the distribution of reduced ${\chi }^{2}$ values if we assume that the observational errors are Gaussian, i.e., if we assume that multiple measurements of the flux at a particular band will distribute according to a Gaussian. The expression for the likelihood is then:

$$\begin{eqnarray}&&P(D| M)=\displaystyle \sum _{i}\mathrm{exp}(-1/2\;{\chi }_{{\rm{red}},i}^{2}),\end{eqnarray} \tag{ 5 }$$

where the sum is marginalized for each model parameter over all possible ${\chi }_{{\rm{red}},i}^{2}$ values for models with a given value of the parameter. The prior P(M) is a measure of any previous knowledge that you have on a particular parameter or set of parameters. For example, if previous independent measurements indicate that the stellar mass of a galaxy must be within certain values, then you can constrain the possible solutions to your fitting problem by constructing a prior on (${M}_{{\rm{cl}}}$) that is compatible with those values. In the present work we use uniform priors to specify ranges of reasonable values for our parameters, and to bias our posterior PDFs as little as possible. Finally, you need to apply a normalization factor to your posterior PDF to guarantee that the probability of at least one model being a representation for your galaxy equals one.

Given the infinite size of the parameter space (we interpolate from the original grid to allow any value of the parameters), it is not possible to calculate the posterior PDF for every single allowed value of the model parameters, especially as more data points are added and additional model parameters are considered. Instead, we use a Monte Carlo Markov Chain (MCMC) approach to step across the parameter space while properly sampling the posterior PDF. The idea is simple and is based on the Metropolis Algorithm: we start at a given location of the parameter space where we can calculate the value of $P{(M| D)}_{{\rm{old}}}$ and then randomly move to another location where we calculate the new value of $P{(M| D)}_{{\rm{new}}}$. We then calculate the ratio $P{(M| D)}_{{\rm{new}}}/P{(M| D)}_{{\rm{old}}}$ and compare it to a random number between 0 and 1. If the latter is lower than the calculated ratio, we accept the step and update the value of the model parameters to those of $P{(M| D)}_{{\rm{new}}}$. If, on the other hand, the ratio of probabilities is lower than the random number, we reject the step and the model remains unchanged. If the step size and the number of iterations are properly chosen, this process should converge to the posterior PDF. In other words, the histogram of models selected in this fashion should be a representation of the initially unknown posterior PDF, and tells us where in the parameter space are the most likely solutions located, given all the information at hand. This is the ultimate solution to our fitting problem.