© The Authors 2024

1. Introduction

Observations of star-forming regions in the Milky Way and nearby disk galaxies have revealed both close-to-universal scaling relations of the cloud properties (e.g., Larson 1981; Rosolowsky et al. 2003; Rosolowsky 2007; Heyer et al. 2009; Miville-Deschênes et al. 2017), and a clear influence of the galactic dynamics (Schinnerer et al. 2013; Hughes et al. 2013; Sun et al. 2020, also supported theoretically, see e.g., Renaud et al. 2013; Meidt et al. 2013; Fujimoto et al. 2014), yet with a subtle impact on the star formation efficiency itself (Querejeta et al. 2021). This suggests complex and multiple coupling mechanisms between the scales of the galaxy, the clouds, and star (cluster)-forming knots.

Thanks to pioneer studies with the Hubble Space Telescope (HST, e.g., Elmegreen et al. 2007; Guo et al. 2012), and more recently with ALMA, JWST, and the power of strong gravitational lensing, a wide window is now opened on the gaseous and stellar inner structures of more distant galaxies, in particular at cosmic noon where the cosmic star formation history peaks (z ≈ 1 − 3, Madau & Dickinson 2014). For instance, with strong lensing, physical resolutions of a few 10 pc can be reached by ALMA at z = 1 (e.g., Dessauges-Zavadsky et al. 2023), and JWST complements such studies by providing integrated properties of the stellar objects (e.g., Claeyssens et al. 2023). Disk galaxies at these redshift host more extreme physical conditions than in the local Universe (e.g., Swinbank et al. 2011), with most of their morphologies being dominated by clumpy structures (Cowie et al. 1995; Elmegreen et al. 2013; Guo et al. 2015; Shibuya et al. 2016; Fisher et al. 2017; Sattari et al. 2023), particularly bright in restframe UV (Förster Schreiber et al. 2009; Zanella et al. 2015; Cava et al. 2018). In the vast majority of cases and up to very high redshift (z ≈ 7, Treu et al. 2023), such clumpy appearance has an internal (“in situ”) origin (e.g., Förster Schreiber et al. 2009; Tacconi et al. 2013; Girard et al. 2018), as opposed to the multi-component structure that galaxy mergers would yield (e.g., Puech 2010; Calabrò et al. 2019).

The magnification from lensing, and/or the gain in angular resolution and sensitivity from JWST revealed the hierarchical nature of (at least some of) the young stellar clumps that host sub-structures below the scales of ∼1 kpc and ∼10^8 − 9 M_⊙ (Meštrić et al. 2022; Messa et al. 2022), contrarily to previously thought (Elmegreen & Elmegreen 2005; Förster Schreiber et al. 2011; Guo et al. 2012). Similar conclusions have also been proposed by adopting observational limitations to detect clumps in simulated galaxies, and finding overestimations by up to a factor of 10 for the clump’s stellar masses (Huertas-Company et al. 2020, see also Faure et al. 2021). These results ask the question of the hierarchical organization of the dense gas and of its star-forming phase in gas-rich clumpy galaxies at z ≈ 1 − 3.

Although clumpy galaxies are found in the main sequence of galaxy formation (Schreiber et al. 2015; Fisher et al. 2019), the peculiar organization of their interstellar medium (ISM), clustering of star formation, and elevated star formation rate (SFR) influence the distribution and luminosity of the traditional tracers of dense gas and young stars (see e.g., Daddi et al. 2015; Renaud et al. 2019a on the α_CO conversion factor in clumpy galaxies), further complicating the interpretation of not-always resolved observational data. Despite the significant recent improvements in precision, completeness, depth, and redshift range, observation programs are still not conclusive on the similarities or differences of star-forming gas structures between nearby galaxies and gas-rich clumpy disks. In particular, it is still not clear whether (and/or at which scale) the structure of the star-forming clouds or clumps and the star formation process are universal across cosmic time, or if new regimes must be invoked in the extreme physical conditions of gas-rich galaxies.

The clumpy nature of gas-rich disks has been well-retrieved in many simulations of galaxies in isolated and cosmological context (among many others, see Noguchi 1999; Agertz et al. 2009; Dekel et al. 2009; Bournaud 2010; Ceverino et al. 2012; Hopkins et al. 2012; Bournaud et al. 2014; Tamburello et al. 2015; Behrendt et al. 2016; Mandelker et al. 2017; Fensch & Bournaud 2021; Renaud et al. 2021a). This peculiar morphology has been attributed to violent disk instabilities triggered by gas inflows onto an existing disk galaxy (Agertz et al. 2009; Dekel et al. 2009; Bournaud 2010), or to hierarchical growth by mergers of smaller clumps (Behrendt et al. 2016). In the first paper of this series (Renaud et al. 2021a, hereafter Paper I), we show that it is the gas-rich nature of these disks which allows for the development of new regime of instability, replacing the classical Toomre’s found at low gas fractions (below 20% of the baryonic mass of the galaxy). Such instabilities are still driven by the stellar component, but with a significant contribution of the molecular phase in the inter-clump regions. This explains the formation of local clumpy structures, instead of the galactic-wide spiral pattern at lower gas fractions. The highly turbulent nature of these clumps imprints a large velocity dispersion to the young stars, to levels potentially sufficient to explain the kinematics of the thick Galactic disk (van Donkelaar et al. 2022, hereafter Paper II).

While qualitatively retrieved in many simulations of gas-rich disks, the quantitative results on the clump properties vary a lot across the literature. Different numerical techniques, in particular for stellar feedback, lead to different conclusions on fundamental properties of the clumps such as their mass (Mandelker et al. 2016; Mayer et al. 2016) and lifetimes (Genel et al. 2012; Bournaud et al. 2014; Oklopčić et al. 2017; Fensch & Bournaud 2021). Conversely, Ejdetjärn et al. (2022, hereafter Paper III) showed that the same level of turbulence could be reached in disk galaxies with and without stellar feedback, thus favoring gravitational instabilities as the main driver of turbulence at galactic scale.

It remains to determine whether the star-forming clumps in gas-rich galaxies are scaled-up versions of the giant molecular clouds in nearby galaxies, requiring a mere extrapolation of the scaling relations established in the local Universe, or if they fundamentally differ in nature. This question of paramount importance is the key to either allow for a transposition to gas-rich disks of our detailed understanding of the physics of the ISM and star formation obtained in multi-wavelength, deep, and highly resolved studies of local galaxies, or on the contrary, to call for new physical regimes specific to the conditions found at high redshift.

In this fourth opus of the From giant clumps to clouds series, we use numerical models of isolated disk galaxies, varying only the gas fraction. These simulations lead to spiral galaxies and to clumpy morphologies at low and high gas fractions respectively. The subgrid models adopted here have been used in previous works focusing on the diversity of clouds and the role of feedback in the organization of the ISM, which have concluded to remarkable matches with observations of the ISM of Milky Way-like galaxies (Grisdale et al. 2017, 2018, 2019). Here, we examine the properties of gas clumps at various gas fractions, with the aim of testing the hypothesis of universal scaling relations across redshift.

2. Method

2.1. Simulations

The simulations and numerical method are introduced in Paper I, and briefly summarized here. We use the RAMSES code (Teyssier 2002) to model a 7 × 10¹¹ M_⊙ galaxy comprising a dark matter, a central bulge, a gaseous and a stellar disk initially setup with the MAGI code (Miki & Umemura 2018). For the sake of control on the parameters, the cosmological context is ignored. Using these initial conditions, we vary the relative mass of the stellar and gaseous components to generate three models, every other parameter being left unchanged. After running the simulations for an initial relaxation phase of ∼100 − 150 Myr, the gas fraction f_gas = M_gas/(M_⋆ + M_gas) in the galaxies is 10%, 25% and 40%. These models are then called ℱ10, ℱ25, and ℱ40 respectively.

The simulations reach a resolution of 12 pc and include heating from a cosmological background calibrated at z = 0 (Haardt & Madau 1996), atomic and molecular cooling, star formation at an efficiency (SFE) of 10% per free-fall time in the gas denser than 100 cm⁻³. The stellar feedback model follows that of Agertz et al. (2021), that is winds, radiation pressure, type-II and Ia supernovae (SNe) and the production of oxygen and iron. In post-processing, the fractions of atomic and molecular gas are computed from cooling tables assuming collisional ionization equilibrium in every cell of the simulation, and following Krumholz et al. (2009).

As detailed in Paper II, the gas fractions and resulting star formation activity in these runs can be related to that of disk galaxies at different redshifts. An approximate correspondence with the SFRs of Milky Way-like galaxies and supposed progenitors observed in van Dokkum et al. (2013) matches our ℱ10, ℱ25, and ℱ40 cases with the redshifts ≈0, ≈1.5 − 2, and ≈2 − 2.5 respectively. We stress however that our controlled experiments are designed to keep the same masses and sizes for all galaxies as we only vary the gas fraction, which implies that our three cases are not representative of different stages along the evolution of the same galaxy.

2.2. Clump detection and statistics

Our working resolution of 12 pc is sufficient to detect gas overdensities (clouds and clumps), but not to resolve their inner morphologies and sub-structures. While such fine details are irrelevant for the present paper in the case of giant molecular clouds found at low gas fractions, the situation could be different in much larger and more massive clumps in gas-rich disks, which could host distinct star-forming sub-structures (e.g., Behrendt et al. 2016). By keeping this limitation in mind, we focus on a specific level of the hierarchical ISM, that is the cloud-like gas structures with a high density contrast from their surroundings. Even if our simulations do resolve some of them, sub-clumps are not examined in this work.

We identify gas clumps as over-densities in face-on surface density maps (> 300 M_⊙ pc⁻²). Their size is then determined by iteratively adding concentric shells centered on the center of mass of the over-density, until the density of gas of a shell is less than 10% of that of the peak. By doing so, we focus our analysis on the densest structures, and avoid spiral-like elongated shapes. Once this identification is done, the clump is cataloged as a sphere centered on the maximum density¹. In the case where several density peaks are not separated by a diffuse enough medium and thus overlap, the associated clumps are considered as sub-structures of a single object. As mentioned above, our working resolution is not sufficient to systematically properly capture the formation and evolution of sub-structures. Therefore, the overlapping clumps are merged such that only the largest (containing) one is retained. As a consequence, the center of mass of a clump can be offset from its maximum density. We come back to this point when examining the density profiles of clumps in Sect. 3.3.1.

We have checked that this method provides an identification similar to what a visual examination would lead to, as illustrated in Figs. 1 and 2. In accordance with the objectives of this work, our identification method relies on a structural argument rather than on the gas phase. The main motivation for this choice is to adopt an identification process which works equally well at all gas fractions, without having to adjust parameters between the runs as it could be required by the different morphologies and properties of the gas over-densities. This implies that clumps in our definition do not necessarily coincide with the molecular phase of clouds, but rather with gas overdensities, including their atomic envelopes. Therefore, our clumps are significantly more massive than GMCs in the Milky Way, even in the low gas fraction case. Including the atomic envelope in our selection also provides consistency with the resolution limitations discussed above: molecular sub-clumps of the same structures are merged when they share a (dense enough) atomic envelope. We note that Grisdale et al. (2018) used the same sub-grid recipes as ours and initial conditions comparable to our ℱ10 case to conduct a direct comparison with observed molecular clouds. The good match between their simulations (using mock CO projected maps) and the observations such as the power spectrum of the ISM and the clump mass function validates the numerical we also use. This allows us to focus now on physical mechanisms involved in the clumps formation and evolution, in cases where the match with real objects is more difficult to assess. Therefore, we do not adopt a selection method motivated by observations, but rather by physical considerations.

Fig. 1. Surface density of gas in the three simulations, at the instant considered in the entire analysis. Red circles indicate the position and size of the gas clumps detected (see text).

Fig. 2. Example of a massive clump hosting sub-structures, as a zoom-in on Fig. 1 (right panel, fgas = 40%). The contours indicate the surface densities of molecular gas at 10, 100 and 1000 M⊙ pc−2. The red circles shows our clump identification. As visible in this example, some of the largest clumps host clear sub-structures, which are inexistent in smaller clouds (either physically or because of limited numerical resolution).

As the bottom of the potential well, the central region of disk galaxies is the convergence point of many gaseous and stellar structures, such that the properties of these structures rather fall under the influence of complex evolutionary processes such as accretion and mergers. This makes the interpretation of the clump properties very involved. For this reason, we ignore the gas clumps in this volume, by taking an exclusion zone within a conservative galactic radius of 2 kpc.

This methodology gives approximately 70, 110 and 150 clumps, with average radii of ≈120, 160 and 180 pc, for the 10%, 25% and 40% cases respectively. To improve the statistics, this detection is repeated on 5 consecutive snapshots (spanning about 20 Myr) for each simulation.

Throughout the paper, we use robust statistics to limit the weight of outliers in our spare datasets. We use the median instead of the mean, we compute the robust standard deviation as the median absolute deviation (MAD) divided by 0.6745, and we perform median-based fits rather then least-squares fits, as in Romeo et al. (2023, see their Sect. 2.3 for more information). Furthermore, to robustly identify the outliers of the relation between the measurements {x₁,  …,x_N} and {y₁,  …,y_N}, we compute the best-fitting relation y = f(x), and the threshold T defined as:

$$\begin{array}{c}\hfill T=\sqrt{2ln(N)}\frac{\mathrm{median}\left({\left\{|y_i-f(x_i)|\right\}}_{i=1}^N\right)}{0.6745},\end{array}$$(1)

where N is the number of points, and the fraction is the robust standard deviation (see Eqs. (9), (10), and (12) of Romeo et al. 2004, for details). A data point (x_i, y_i) is robustly considered as an outlier if |y_i − f(x_i)| > T.

3. Results

3.1. Disruptive mechanisms

Different regimes of instability between the low and high gas fractions in Milky Way-mass disk galaxies have been shown to imprint the morphology of dense gas clumps, with a transition between the two regimes occurring at gas fractions of 20% (Paper I, see also Paper II on kinematics). We examine here whether these differences are amplified or mitigated by environmental effects within and around the clumps.

Figure 3 shows, in logarithmic bins of gas density, the intensity of turbulence, tides, and shear, measured at a scale of 120 pc (i.e., over 10 resolution elements, and of the order of the average scale of clumps). The similarities in the gas density probability distribution functions (PDFs, Fig. 3 in Paper I) justify using the same density binning for all three gas fractions. The adopted bins approximately encompass the very diffuse ISM found in supernovae bubbles and galactic outskirts (∼10⁻¹ cm⁻³), the inter-arm or inter-clump regions (∼10¹ cm⁻³), the structures connected to the clumps, such as spiral arms (∼10² cm⁻³), and the star-forming clumps themselves (≳10³ cm⁻³). The measurements are done on a regular Cartesian grid, without centering on the clumps of any other structure. The tidal force is computed as the maximum eigenvalue of the tidal tensor (i.e., minus the Hessian of the gravitational potential) via first-order finite differences of the gravitational acceleration. Therefore, it is not the azimuthally averaged tidal strength, but rather the force along the direction where local tides are the strongest (see Renaud et al. 2011, 2017 for details). The shear force is evaluated as the quadratic sum of the off-diagonal terms of the Jacobian matrix of the velocity field, also computed with first-order finite differences (see Renaud et al. 2015). This quantity can be interpreted as the squared rotational frequency term in the expression of the local centrifugal force. These raw values are rather difficult to interpret and thus we present them in normalized forms, using arbitrarily the case at f_gas = 40% as reference.

Fig. 3. Average Mach number (star symbols), and intensity of tidal (circles) and shear (triangles) forces in gas density bins, normalized to that in the fgas = 40% case. Empty symbols denote the forces normalized to the local self-gravity. All quantities are measured at a scale of 120 pc. Symbols have been slightly offset horizontally for clarity. The approximate distinctions between three different media are identified at the bottom, for reference. The two most gas-rich cases (red) yield comparable values at all densities for all these quantities, significantly higher than in the low gas fraction galaxy (blue).

The turbulence, as traced by the Mach number ℳ, is very similar between the two gas-rich cases with less than a 10% difference at any density (as shown by the red star symbols in Fig. 3). It evolves from the transonic regime (ℳ ≈ 0.6) at the lowest densities, to highly supersonic (ℳ ≈ 25) inside the star-forming clumps. However, the low gas fraction stands apart with a turbulence about 1.5 times weaker in ℱ10 (blue star symbols), with this ratio showing no strong nor clear dependence on the gas density. The fact that gas rich disk galaxies have significantly higher levels of turbulence than those with more modest gas fractions has been well established observationally (see e.g., Förster Schreiber et al. 2009; Fisher et al. 2017; Girard et al. 2019). Possible explanations invoke the continuous injections of turbulence by internal mechanisms (such as feedback from a sustained star formation and kpc-scale disk dynamics, see e.g., Hoffmann & Romeo 2012; Agertz et al. 2015; Fisher et al. 2019; Brucy et al. 2020), and stirring of the ISM by external factors such as gas accretion and more frequent mergers at high redshift (Agertz et al. 2021; Renaud et al. 2021b; Forbes et al. 2023). Using simulations of Milky Way-like progenitor galaxies in isolation, several works reported however that external factors are not necessary to develop and maintain high levels of turbulence in gas-rich disks at redshifts ≈1 − 3 (Clarke et al. 2019; Khoperskov et al. 2021; Paper II; Paper III), and that galaxies evolve intrinsically from strong to moderate levels of turbulence as the gas fraction decreases. Our works show that, all other things being equal, this transition is sharp, with strong turbulence maintained down to gas fractions around ≈20%, before it drops (see e.g., the Fig. 2 of Paper II). Owing that turbulence support is of paramount importance in maintaining a gas disk geometrically thick and dynamically hot (the other main factor being the tidal excitation from satellites and mergers, see Renaud et al. 2021b), our results imply that, in the absence of mergers, disk galaxies in the mass range we consider undergo a rapid transition from thick to thin disk when their gas fraction reaches ≈20%.

Similar trends are also found for the tides and shear, although the differences between the two gas-rich cases are more important than for the Mach number. Yet, they remain within a factor of ≈1.3. Here again, larger differences are found between the ℱ10 and ℱ40 cases: tides are ≈1.5 − 12 times weaker and shear is ≈2 − 6 times weaker at low gas fraction. This further adds to the similarities between the two gas-rich cases, over the entire range of densities in the galaxies, and a sharp transition when moving towards lower gas fractions. This contrasts with the results of Fensch & Bournaud (2021) who reported a shear parameter increasing with decreasing gas fraction. We note however that their parameter does not measure the full shear, but only the differential rotation of the disk via the difference of tangential galactic velocities around clumps. The high velocity dispersions of clumpy disks however implies that the radial and vertical components of the velocity field can also carry non-negligible contributions to the total shear. In addition, the shear parameter of Fensch & Bournaud (2021) does not account for the scale at which it is measured (while ours does, through the transverse velocity gradients). As the characteristic scale of instabilities and structures differs between giant clumps in gas-rich galaxies and spiral arms at low gas fraction (Paper I), these two differences likely explain why our definition of shear leads to opposite conclusions to that of Fensch & Bournaud (2021).

Tides and shear result from a complex interplay of phase-space structures, with contributions from both large and small scales, including from the gas clumps themselves. To estimate the effective effects of these forces on the gas structures, we normalize them by the local self-gravity (open symbols in Fig. 3). This tends to reduce the differences between the galaxies, especially in the inter-clump medium (∼10¹ − 10³ cm⁻³), but the low gas fraction case still stands apart from the other two. The presence of massive clumps is thus an important driver of the phase-space organization in their vicinity (i.e., in the inter-clump medium), via their local gravitational influence (as shown by the tidal force) and other kinematic effects (possibly outflows from feedback interfering with the large-scale disk rotation, as shown by the shear). We come back to this in Sect. 3.3.1.

Overall, Fig. 3 indicates that the main agents capable of disrupting or destroying structures are significantly stronger at high f_gas, at all densities. This combines with a higher SFR, and thus a more energetic disruptive effect from stellar feedback. This could point toward a shorter lifetime of the clumps at high than low gas fraction. However, the details of the organization of the ISM, the size of the structures and their density profiles, that is their fragility to destructive factors, must be accounted for before concluding on the rapid dissolution of clumps. We discuss these aspects in the next sections.

3.2. Clump mass function

The gas density PDFs (Fig. 3 of Paper I, reproduced in Fig. A.1 below) show important similarities between the three cases, including in the densest phase. In particular, the two gas-rich cases do not host much denser gas than ℱ10. This naturally implies that the large clumps at high f_gas are not significantly denser, but rather more massive than the clouds at low gas fractions. This is confirmed by the clump mass function shown in Fig. 4. Recall that our clump identification method does not probe the low-mass end of this distribution, and we thus focus here on the high-mass end. The maximum clump mass increases with f_gas, leading to a dependence of the high-mass end slope of the mass function with the gas fraction.

Fig. 4. Clump mass functions for the three cases (using staking from 5 consecutive snapshots, to improve the statistics). The dashed lines indicate the best-fit of the high-mass end using a power-law. This gives indexes of −2.0, −1.5 and −1.6 for the ℱ10, ℱ25, and ℱ40 cases respectively. The dotted curves are fits to the data of fgas = 40% using a Schechter function (see text). This functional form does not provide a good match to the other cases.

The best-fit power-law relations to the mass functions yields indexes of −2.0, −1.5 and −1.6 for the three cases respectively, all within the range of values observed across galactic environments (≈ − 1.5 to −2.9, see e.g., Rosolowsky 2005). The mass function of the ℱ10 case is roughly similar to the −2 index well-established in initial star cluster mass functions (e.g., Whitmore et al. 1999; Larsen 2002; Bik et al. 2003; Chandar et al. 2014), and often presented as an outcome of gas fragmentation driven by gravo-turbulence (Elmegreen & Efremov 1997). However, higher gas fractions yield slightly shallower and top-heavy clump mass functions with respect to the low gas fraction, but also in comparison with the observed slopes of −1.7 and −2.0 (depending on the method adopted) in gas-rich disks at z ≈ 1 − 3.5 (Dessauges-Zavadsky & Adamo 2018). This difference means that gas-rich disks do not simply host massive clumps in addition to those found at lower f_gas, but that they also count more numerous clumps at intermediate masses. The shapes of the mass functions of the two gas-rich cases are very similar, but they differ by the overall number of clumps.

We note that the mass function from the ℱ40 run shows a remarkable match with a Schechter (1976) function of the form $\propto M_{\mathrm{gas}}^{-0.73}\exp [-M_{\mathrm{gas}}/(1.3\times {10}^8{M}_{\odot })]$ (dotted line in Fig. 4), while this is not the case for the other galaxies. Some studies proposed that the initial mass function of star clusters follows a Schechter form (Larsen 2009; Adamo et al. 2017), which could suggest a direct mapping between our clumps and observed clusters. However, the presence of an exponential cutoff differencing the Schechter functional form from the classical power-law could have a physical or a statistical origin. As discussed in Mok et al. (2020), small datasets could induce an apparent cutoff in the mass function, simply because of the statistical lack of objects to populate the high-mass end of the distribution. This would then differ from a cutoff imposed by physical mechanisms. Our analysis suffers from the same problem, and we are not able to conclude on the origin of the exponential cutoff. Yet, it is noteworthy that only the galaxy with the highest gas fraction (which is also the one hosting the most clumps and thus which is potentially the least sensitive to a statistical effect on the cutoff) is the only case yielding an agreement with the Schechter functional. If physically-driven, the reason for a Schechter mass function only in the gas-rich case is difficult to identify. It could result from a truncation of the most massive clumps, possibly by stellar feedback and/or environmental effects (Sect. 3.1). Hopkins et al. (2012) proposed an analytical model retrieving the exponential cutoff at the high-mass end of the mass function of GMCs. We note however that this is a direct consequence of the assumptions made on the shape of the power spectrum of the ISM, thus leaving open the question of the fundamental cause(s) of such a cutoff (if physical). From our experiments, it is not clear whether the truncation is intrinsic to massive clumps (regardless of other properties), or rather to the high gas-fraction. To investigate these hypotheses, the next section focuses on the structure of the most massive clumps in our galaxies.

To examine whether our conclusions on the influence of the gas fraction on the clump mass function is robust to a change of the subgrid recipes employed, we show in Appendix A examples of the same analysis conducted with different parameters and ingredients for the star formation and feedback models. We also show the density PDFs, for reference. The conclusion is that, within the relatively wide range of parameters explored, changing the subgrid models does not alter the clump mass function as strongly as the gas fraction does.

3.3. Clump properties

3.3.1. Clump density profile

Figure 5 shows the gas density radial profiles of the five most massive clumps in our three galaxies (excluding the central-most, bulge-like, overdensities). Two sets of profiles are computed, using the maximum density and the center of mass as centers. The difference originates from the finite resolution of the simulations², and the hierarchical structure of clumps, as discussed below (also recall Fig. 2). In the following, we ignore the central parts of the profiles (up to ≈30 − 40 pc), affected by the resolution limits. Potential sub-structures in these clumps do not necessarily exist at the same distance from the center of each clump, and thus their contributions are rapidly averaged out in the median radial profiles shown here. In this section, we ignore the size of the clumps and show the radial profiles from the clump centers up to large distances.

Fig. 5. Median density profile of the five most massive gas clumps in each simulation (excluding the central 2 kpc of the galaxies, where radial inflows at galactic scale play a non-negligible role in the accumulation of gas). The top panel shows the radial profiles with the maximum density as the center, and the bottom panel those centered on the center of mass. The shaded areas show the robust standard deviation over these 5 clumps. The dashed line represents a power-law of index −2, for reference.

Once again, the cases at f_gas = 25% and 40% are remarkably similar, while the ℱ10 case stands apart. Clumps in the two gas-rich cases are slightly but significantly denser than those in the ℱ10 galaxy by factors ranging from 1.5 to 3.1, up to several 100 pc, that is where reaching the inter-clump medium. However, the density of clumps from the two gas-rich disks are remarkably similar at all radii. Therefore, it is the organization of the ISM, initiated by the instabilities discussed in Paper I³, and not (directly) the increase of the gas mass of the disk which is responsible for the increased density of the most massive clumps in gas-rich galaxies.

When using the maximum density as the center of the clumps, the clumps profiles exhibit power-laws slightly steeper in the gas-rich disks (≈ − 2.7 to −2.2, depending on radius) than at low gas fraction (≈ − 1.8 to −1.5), but always close to the −2 index found for an isothermal sphere.

When centering on the center of mass, the presence of sub-structures near but not exactly at the center creates density peaks at various radii and azimuths, which average out as a flattening of the density profiles as illustrated in Appendix B. The central cores in Fig. 5 are thus a statistical effect, and are generically not representative of the actual profile of a given clump taken alone.

By selecting the most massive clumps of each simulation, we introduce a potential bias because the high-mass end of the clumps mass functions of the three simulations do not overlap (recall Fig. 4). Yet, it is not obvious whether this difference in the clump mass is the sole explanation for the differences in the density profiles. Figure 6 confirms this hypothesis by showing that the density profiles of clumps selected in a common mass range populated in all three galaxies (2 × 10⁷ − 1 × 10⁸ M_⊙) are remarkably similar. In fact, the differences between the three cases are even smaller than the clump-to-clump scatter within a given simulation (shown as shaded areas). Therefore, at a given mass, a gas clump has a density profile independent of the gas fraction of its host disk. The different profiles only arise because of the presence of more massive clumps in gas-rich galaxies. As explained above, these more massive cases are denser, and possibly yield a central core. (This conclusion could be affected by our resolution limitations.) Therefore, such clumps ought to be more resistant to external disruptions (Sect. 3.1), and internal effects such as feedback. This suggests that massive clumps in gas-rich galaxies can be longer-lived than local GMCs, which supports, at least qualitatively, the conclusions of previous works on clump survival over several orbital times (≳100 Myr, see e.g., Bournaud et al. 2015; Dekel et al. 2023, but see Genel et al. 2012; Oklopčić et al. 2017 for opposite findings). Our results illustrate the need for capturing the inner structure of the clumps and the spatial (and time) scatter of sources of feedback simultaneously with galactic dynamics to definitely conclude this long-lasting debate.

Fig. 6. Median density profile of gas clumps in each simulation (as in Fig. 5), but only considering clumps in the same mass interval. Clumps of comparable masses have similar density profiles regardless of the gas-fraction of their host galaxy.

3.3.2. Larson scaling relations

In the previous section, we found that some properties of gas clumps such as their density profiles are independent of the gas fraction of their host galaxy, at a given clump mass, but that a different regime of instability in gas-rich disks allow for the formation of more massive and denser clumps (Paper I). This duality questions the universality of scaling relations for gas clouds, which we examine here.

Figure 7 shows the distributions of all the clumps detected in the simulation, over five consecutive snapshots to improve the statistics, in mass M, radius r, and velocity dispersion σ. Observational census of molecular clouds in the Milky Way (e.g., Falgarone et al. 2009; Heyer et al. 2009; Miville-Deschênes et al. 2017), nearby galaxies (e.g., Hughes et al. 2013), and some examples are high redshift (Dessauges-Zavadsky et al. 2023) report close to universal relations between these quantities, as first noticed by Larson (1981). Recent estimates of these so-called Larson relations in the Milky Way are M ∝ r^{2.2 ± 0.2} and σ ∝ r^{0.63 ± 0.30} (Miville-Deschênes et al. 2017). As noted above, our method does not only detect the molecular contents of clumps, but also their atomic envelopes around the density peaks. Indeed, the diversity of clump sizes and morphologies between low and high gas fractions, and the presence of sub-structures in some of the most massive clumps (Fig. 2) imply that the molecular phase alone does not necessarily cover the entirety of a clump structure, specially in gas-rich disks. The consequence of this selection criteria is that, compared to molecular clouds, our clumps encompass a larger amount of gas at intermediate densities. Therefore, the mass–size relation of our clumps is necessarily shallower than that of molecular clouds, as confirmed by the fits in Fig. 7. These fits are computed in log-log space as least absolute deviation fittings, that is linear regressions similar to the least squares technique, but using absolute values instead of square values to reduce the influence of outliers. We perform this computation following the medfit algorithm of Press et al. (1992).

Fig. 7. Mass–size (top) and velocity dispersion-size (bottom) relations for the clumps in our simulations. The color of the points encodes the logarithm of the third quantity (velocity dispersion on the top panels, and mass in the bottom row), increasing from dark to light. The solid lines are the least absolute deviation fittings of the linear relations in log-log space. They are translated into a power-law in linear-linear space, with the indexes indicated at the top of each panel. The dashed lines show constant densities, for reference, and the dotted lines mark the threshold for the robust identification of outliers (see text, Eq. (1)).

We indicate in Fig. 7 the Spearman coefficient and its p-value for each relation⁴. For the three gas fractions, the mass and size of clumps yield relatively strong correlations of high significance. However, the correlations are much weaker for the velocity dispersion-size distributions, due to an important scatter which increases with the gas fraction. This scatter mainly originates from clumps having high velocity dispersions, interestingly at all masses and sizes. When reported to the escape velocity (estimated roughly by using the total mass and radius as $\sqrt{2GM/r}$), the high σ’s indicate that these clumps are the least bound, and on the verge of dissolving. We note that the outliers in the r − σ relation are not necessarily the same as those in the r − M relation. Therefore, the elevated velocity dispersion is not, or at least not always, a consequence of a mass higher than average.

A possible reason for an elevated velocity dispersion in a clump is the injection of energy and momentum by stellar feedback. When examining the indicators of recent star formation activities in these clumps (star formation rate, efficiency, and depletion time, see Sect. 3.3.3), we find that they cover a large range of values, not exclusive to the high-σ clumps. Thus, stellar feedback cannot be the sole factor of the high levels of velocity dispersion for all these clumps. As already noted in Sect. 3.1, (hydro)dynamical disruptive effects from the galactic environments are stronger at high gas fractions, which probably plays an important role, with feedback, in creating the scatter in the distributions of velocity dispersion. Telling apart the external (tides, shear, large-scale turbulence) and internal (feedback) effects on a clump-by-clump basis requires dedicated experiments, out of the scope of this paper. Nevertheless, we examine the star formation activity in the extreme clumps in the next section.

As hinted in Sect. 3.3.1, clumps in the gas-rich disks are not necessarily denser than their counterparts of the same mass at low gas fraction. The similarity of the best fit for the Larson’s relations between the three cases indicates in fact that the denser clumps, outliers of the mass–size relation which are more numerous at high gas fraction, are not necessarily the most extended nor the most massive.

The PHANGS-ALMA program performed a census of molecular clouds in nearby disk galaxies and examined their mass–size and velocity dispersion–size relations across a diversity of galactic environments (Rosolowsky et al. 2021). They reported similar average relations as for the Milky Way clouds, but with significantly elevated scatter, which they attribute to both technical limitations (e.g., resolution bias and source blending) and a diversity of physical conditions (e.g., external forces).

At higher redshift, Dessauges-Zavadsky et al. (2019, 2023) combined the high angular resolution of ALMA with strong lensing to detect massive clumps in gas-rich galaxies at redshift ≈1. When examining the properties of the massive molecular clumps, they found an offset with respect to the Larson’s relations from local GMCs, in a similar fashion as clouds in local starbursting galaxies: these clumps tend to be more massive, denser, and with a higher velocity dispersion than local GMCs. These results also show an increased scatter in the two scaling relations. This led the authors to question the universality of the gas clumps scaling relations, and to confirmed the prime importance of the galactic environment in setting up the properties of clouds, which had been previously noted observationally (e.g., Swinbank et al. 2011; Hughes et al. 2013; Sun et al. 2020) and theoretically (e.g., Renaud et al. 2013; Grisdale et al. 2017).

Our results are compatible with the findings of Dessauges-Zavadsky et al. (2019, 2023), but we nuance their interpretation: the census of gas structures in distant galaxies is biased toward dense and large objects which overcome the resolution and signal-to-noise limits imposed by the instrument, the gravitational lensing, and the methodology adopted. Therefore, it is likely that the smallest and least dense clumps are not detected, or do not satisfy the selection criteria. In our analysis, the densest and, to a lower extent, the largest clumps are precisely the ones the furthest from the median mass–size relation (top row of Fig. 7), and hence the main contributors to the increasing scatter of this relation with the gas fraction. Considering these extreme clumps only would then, indeed, lead to the conclusion that local scaling relations break in gas-rich galaxies.

From our results, the emerging picture is instead that of a rather universality in the average properties and scaling relations of gas clumps (mass, size, velocity dispersion, density profile), together with a scatter in these properties which significantly increases with the gas fraction of the host galaxy, and allows for more frequent extreme outliers. Important information on the variations of the star-forming regions with the galactic environment and also their evolution with cosmic time could then be revealed from the study of these outliers, as we show next.

3.3.3. Extreme clumps and star formation activity

Figure 7 shows that the three galaxies yield comparable scaling relations, including in terms of intercept or normalization. The main difference resides in the scatter which considerably increases with the gas fraction. Gas-rich disks thus provide the physical conditions allowing some of their gas clumps to depart from the median scaling relations and to reach extreme properties, at formation and/or after some dynamical evolution. To investigate the underlying reasons behind this scatter, we quantify the deviation of a clump from the best-fit relations plotted in Fig. 7, as the excess (or deficit) each clump has with respect to the relation. For instance, we compute the excess mass Δ_M of a clump of mass M and radius r as

$$\begin{array}{c}\hfill {\mathrm{\Delta }}_{\mathrm{M}}=log\left(\frac{M}{b{r}^a}\right),\end{array}$$(2)

where a and b are the best-fit parameters of the mass–size relation M = b r^a. In other words, Δ_M is the vertical distance between a point and the line in the top panels of Fig. 7. It is counted negatively for a deficit of mass with respect to the best-fit relation. We compute the excess mass Δ_M and excess size Δ_r using the mass–size relations, and the excess velocity dispersion Δ_σ using the velocity dispersion–size relations, but studying Δ_M provides the strongest conclusions.

Figure 8 shows the excess mass Δ_M with respect to the star formation rate (SFR), the depletion time (τ_dep), and the star formation efficiency per free-fall time (SFE, ϵ_ff) for each clump. To measure these quantities, we consider the timescale for detection of star clusters via free-free emission t_y chosen equal to 4 Myr, and the mass M_⋆, y of stars younger than t_y in the clump (see Grisdale et al. 2019 for details). This gives:

$$\begin{array}{c}\hfill \mathrm{SFR}=\frac{M_{\star ,y}}{t_y},{\tau }_{\mathrm{dep}}=M\frac{t_y}{M_{\star ,y}},{ϵ}_{\mathrm{ff}}=\frac{t_{\mathrm{ff}}}{t_y}\frac{M_{\star ,y}}{M_{\star ,y}+M},\end{array}$$(3)

Fig. 8. Relation between the star formation rate (top), depletion time (middle), and star formation efficiency (bottom), and the excess of mass of clumps with respect to the best-fit mass–size relations of Fig. 7. The vertical dashed lines mark ΔM = 0, i.e., a clump sitting on the best fit relation shown in Fig. 7. The vertical dotted lines show the thresholds for the robust definition of outliers (see text). Normalized marginal distributions are shown on the right column.

using the average free-fall time of the clump: $t_{\mathrm{ff}}=\sqrt{4{\pi }^2{r}^3/(32GM)}$.

For reference, we indicate in Figs. 7 and 8 the thresholds defining outliers of the relations (recall Eq. (1)). For instance, a clump with excess mass above T is statistically robustly identified as an outlier of the mass–size distribution. This leads to 1.7%, 4.8%, and 8.2% of the clumps to be classified as outliers in mass, for the ℱ10, ℱ25, and ℱ40 cases respectively, all having a positive excess of mass (Δ_M > 0).

The normalized distributions of the SFR (top-right panel in Fig. 8) reveal, once again, the similarities between the two gas-rich cases which share comparable distributions and ranges of values. At low gas fraction however, only the low-SFR end of the distribution matches that of the other galaxies. No clump in ℱ10 reaches SFRs as high as the extreme cases in the ℱ25 and ℱ40 cases. These high values are found in the clumps with the highest excess of mass, which are more numerous and further from average in the gas-rich disks. An analogy with the main sequence of star formation at galactic scale would suggest that the highest SFRs are found in the most massive clumps (e.g., Tacconi et al. 2018). We show in Appendix C that this is not the case: a much stronger correlation is found between the SFR and the excess of mass than with the clump mass itself. At first order, this can be interpreted as a result of elevated density at high Δ_M, but the complex structure of clumps and regulation effects from feedback likely complicate this picture. We thus present further analysis, using other indicators of the star formation activity.

The middle row of Fig. 8 shows that all the clumps at high Δ_M yield short τ_dep, that is fast star formation, regardless of the galactic gas fraction. However, the small scatter of the mass–size relation at low gas fraction implies that only a few such clumps exist in the ℱ10 case. Figure 9 complements this analysis by showing that the most extreme clumps (Δ_M ≳ 0.8 and log(τ_dep/1 Myr) ≲ 1.2) are also those with the highest average density ⟨ρ⟩=M/(4πr³/3). However, clumps with moderate excess mass (Δ_M ≈ 0.3 − 0.8 dex) still have short depletion times (a few 10 Myr), but yield a wide range of densities (⟨ρ⟩∼0.1 − 10 M_⊙ pc⁻³).

Fig. 9. Distributions of average density (total mass divided by total volume) and excess mass of the clumps, with color indicating the depletion time. Correlations between these three quantities are only found for the extreme clumps at high ΔM.

Examining the density radial profiles of individual clumps suggests that the contrast between the density in the center and the outskirts increases with Δ_M. Doing this analysis systematically is difficult because of the dual effect of Δ_M and M itself on the density profile. To illustrate this, we select clumps of a fixed mass, and study the influence of their excess mass on their density profiles. We choose to select all the clumps in the mass interval 2 × 10⁷ M_⊙ < M < 3 × 10⁷ M_⊙. This interval is large enough to contain several clumps in all the three simulations and thus limit the statistical noise, and narrow enough to restrict the analysis to a precise mass. We then group clumps in three bins of Δ_M: close to average (−0.25 < Δ_M < 0.25), intermediate (0.25 < Δ_M < 0.75), and extreme (0.75 < Δ_M), and compute the median of the density profiles in each bin. The results are displayed in Fig. 10. Note that clumps are less numerous in bins of high Δ_M, and thus the profiles are more subject to statistical fluctuations. Nevertheless, clumps with an excess mass (Δ_M > 0.25) yield steeper density profiles than those closer to the scaling relations (Δ_M ≈ 0), meaning that such extreme clumps tend to be more centrally concentrated, rather than only simply denser on average (as the sole information of their excess mass for their respective size could indicate). The high peak density in extreme clumps can be connected to their short depletion times. In fact, previous works have shown that reduced depletion times results from an excess of dense gas at galactic (Renaud et al. 2014; Segovia Otero et al. 2022) and cloud scales (Renaud et al. 2019b), triggered by tidal and turbulent compression in galaxy mergers. For the present study of disks in isolation, we find that elevated Δ_M’s correspond to steep density profiles, and short depletion times, in line with the conclusions in galaxy mergers. The source of the compression of clumps in our isolated galaxies is not clear, but is likely an outcome of the instability regime highlighted in Paper I. The non-linearity of the growth of these instabilities, and their collapse leading to star formation make it difficult to push this investigation further without dedicated and controlled experiments at high resolution.

Fig. 10. Median density radial profiles of clumps of radius between 2 × 107 M⊙ and 3 × 107 M⊙, grouped in three bins of excess mass ΔM. The numbers in square brackets indicate the slope of each curve measured in log–log space in the inner 100 pc. (The conclusions do not change when considering larger radial ranges.) Higher excess masses correspond to more centrally concentrated clumps (steeper profiles).

The distribution of depletion times (right column of Fig. 8) is slightly narrower at low gas fraction, due to the lack of clumps with extreme mass excesses (as discussed above), but also the rarity of average clumps with very long depletion times. The two cases at high gas fraction yield very similar distributions of τ_dep.

Finally, the bottom row of Fig. 8 shows the distribution of SFEs. Clumps with an excess of mass yield only high SFEs, while those closer to average display a large range of efficiencies, over several orders of magnitude. Contrary to the SFR and the depletion time, the distributions of SFEs are very similar for the three gas fractions. These results support the conclusions of Grisdale et al. (2019) of a wide diversity of efficiencies in clouds, with no clear correlation with the cloud properties. As found in Grisdale et al. (2019), the ad hoc value of the SFE per free-fall time implemented in the subgrid recipe for star formation (10% here) is not representative of the SFE at cloud-scale. Part of this discrepancy originates from the range of densities (and thus local free-fall times) found within a given clump. This range tends to be limited in the densest cases, which could explain the peak of the SFE distribution around 10% in Fig. 9, especially at high Δ_M. However, other factors must be taken into account to explain the scatter in SFE, such as the fact that clumps are found (and observed) at various stages of their star formation histories, and thus with evolving stellar-to-gas mass ratios. In addition, Grisdale et al. (2019) reported that the stellar feedback and galactic dynamics make the evolution of the SFE of clumps non-monotonic but rather complex, further increasing the diversity across a galactic disk. Some of the clumps in our gas rich disks undergo erosion and mergers with their neighbors, usually in a relatively dense background environment, which makes their tracking in time very difficult. We have not been able to systematically follow the evolution of our clumps, and thus leave these questions for future works.

For all three gas fractions, the excess velocity dispersion does not show strong nor significant correlations with the star formation indicators (SFR, τ_dep, SFE), as illustrated in Fig. C.2. Such a relation would be expected if stellar feedback would be the main source of turbulence within the clumps. However, galactic dynamics play an important role in the setting the turbulence cascade and thus the velocity dispersion in clumps (e.g., Hoffmann & Romeo 2012; Renaud et al. 2013; Ginzburg et al. 2022), which blurs the effects of feedback. Indeed, Fig. 7 already hinted that the excess mass does not correlate with the excess velocity dispersion (the darkness of the points does not correspond to their vertical distance to the solid line). This suggests that different factors set the mass and the star formation activity on the one hand, and the size and large-scale velocity dispersion of the clumps on the other.

4. Conclusion

Using controlled simulations of isolated disk galaxies, we have studied how the gas fraction influences the properties of the star-forming gas clumps, and whether universal scaling relations could be identified. Our main conclusions are as follows.

• Destructive factors such as turbulence, tides, and shear are significantly stronger in gas-rich galaxies than at low gas fractions, at all densities of the ISM. The sharp transition found in the instability regime (Paper I) and in the kinematics (Paper II) at f_gas ≈ 20% thus also occurs in the strength of the destructive effects.

• The mass function of gas clumps is slightly shallower at high gas fraction, due to the presence of more massive clumps in gas-rich galaxies. The presence of a Schechter-like exponential cutoff is possible at high f_gas, although subject to statistical effects.

• The radial density profile of a gas clump of a given mass is independent of the gas fraction of its parent galaxy. However, gas-rich disks (f_gas > 20%) host more massive clumps, which are more centrally concentrated (steeper profiles).

• The scaling relations between the mass, size and velocity dispersion of clumps yield the same slope and intercept at all gas fractions, in a similar way as Larson relations. However, the scatter of clumps around the median relations significantly increases with the gas fraction. In other words, on average and at a given mass, clumps have essentially the same properties regardless of the gas fraction, but gas-rich galaxies allow for the formation of extreme clumps not found at low gas fractions, which causes the increased scatter of the relations. Outliers are more massive, more turbulent and tend to be denser than average, which could bias the observational census of clump populations.

• The excess of mass of a clump above the median scaling relation correlates with a high SFR and a short τ_dep, that is abundant and fast star formation. The most extreme cases are found in the most centrally concentrated clumps, that is with steep density profiles. Such extreme star formation activity correlates much more strongly with the excess of mass for a given size, than with the mass itself.

• Clumps close to the scaling relation yield a diversity of SFEs, with no clear link with their properties. However, only high SFEs are found in outliers to the scaling relations, with large mass excesses.

• The variations of the gas fractions have a stronger impact on the properties of clumps and of the ISM than changing the sub-grid recipes for star formation and feedback.

For virtually all the quantities and processes studied in this series of papers, we found remarkable similarities between all galaxies with gas fractions > 20%. This suggests an independence of these quantities and processes with redshift for disk galaxies above z ≳ 1. We suspect however that, if confirmed, this statement might break at higher redshifts (z ≳ 5 − 6), when the galaxy is not yet rotation-supported (Segovia Otero et al. 2022) and thus when other instabilities are probably active, and when the gas contents of small structures and their star formation activity might be affected by re-ionization.

Our results suggest that the properties of star formation (rate, efficiency, depletion time) span large ranges of values with weak or no correlations with the intrinsic properties of the clumps. However, stronger correlations appear for the clumps more massive than average for their respective size. Such clumps are more frequently found in gas-rich galaxies due to a special instability regime favoring their formation (Paper I). This regime influences the environmental factors acting on the clumps (tides, shear, large-scale turbulence), independently of the clump’s intrinsic properties.

The rapid dynamical evolution of clumps which repeatedly merge and split makes it very difficult to track them in a self-consistent manner. This is why we have refrained from tackling the long-lasting question of survivability and lifetime of the massive gas clumps. Nevertheless, we propose here some considerations, to be taken with caution. Dessauges-Zavadsky et al. (2023) found separations of a few 100 pc between their detected gas clumps (identified as above a signal-to-noise threshold in CO (4−3) emission) and young star clusters (10 − 100 Myr, with large uncertainties). They then argued that such separations would require an important drift of the stars if they would have formed in the clumps, and that the stellar ages would rather suggest that their natal gaseous nursery would have already been dispersed by feedback (owing a short lifetime for GMCs in nearby galaxies, ≈10 Myr). Keeping in mind that our definition of clumps largely differs from that of Dessauges-Zavadsky et al. (2023), we note nevertheless that the massive clumps in our gas-rich cases yield molecular sub-clumps separated by a few 100 pc (recall Fig. 2, see also Faure et al. 2021), and with velocity dispersions⁵ up to several 100 km s⁻¹ at these scales (Fig. 7). Therefore, it is possible that sub-structures within a massive clump could be separated by several 100 pc, either intrinsically, and/or after a dynamical evolution of a few Myr. This raises a cautionary note on using physical separations between gas and stellar structures to infer the lifetimes of clumps in gas-rich galaxies.

The denser and possibly cored structure of the extremely massive clumps found in our gas-rich galaxies (Fig. 5) make it possible for these clouds to be long-lived. Clumps closer to the average scaling relation, and thus with properties comparable to the clouds at low gas fractions, would be as short-lived as GMCs in the present-day Milky Way, or even shorter-lived due to stronger disruptive mechanisms (Sect. 3.1). Therefore, our conclusions support the idea of two regimes of cloud lifetimes proposed by Dekel et al. (2023), with the long-lived regime only existing for the most massive clumps in gas-rich galaxies. Reaching more affirmative conclusions on this matter is made difficult by the evolving nature of gaseous structures, an aspect that should be taken into account in future studies.

Acknowledgments

We thank the reviewer for comments that helped improve the clarity of the paper. We thank Marc Huertas-Company, Avishai Dekel, and Joel Primack for fruitful discussions. F.R. acknowledges support provided by the University of Strasbourg Institute for Advanced Study (USIAS), within the French national programme Investment for the Future (Excellence Initiative) IdEx-Unistra. O.A. acknowledges support from the Knut and Alice Wallenberg Foundation, the Swedish Research Council (grant 2019-04659) and the Swedish National Space Agency (SNSA Dnr 2023-00164).

References

Appendix A: Robustness to changes of subgrid models

A.1. Method

To test the influence of the subgrid models used for star formation and feedback on our conclusions, we restarted the simulations presented in the main text while varying one parameter or ingredient at a time:

• high SFE: the SFE used in the star formation model is multiplied by 3 (from 10% to 30% per free fall time),

• low SFE: the SFE used in the star formation model is divided by 10 (from 10% to 1% per free fall time),

• energetic SNe: the energy of individual supernova is multiplied by 10 (from 10⁵¹ erg to 10⁵² erg),

• no pre-SN: the stellar feedback only comprises the type-II and type-Ia SNe (winds and radiation pressure are turned off),

which leads to 3 gas fractions × 5 modifications = 15 runs in total.

To make the comparisons as direct as possible with the fiducial runs, these additional runs are started from an evolved snapshot of the fiducial set, and restarted with the new parameters for 30 Myr. This duration ensures that the modifications effectively impact the stellar feedback budget, at least at cloud scale. Capturing a full baryonic cycle would require waiting longer for the re-cooling, re-condensation, re-accretion of the feedback ejecta at coronal and galactic scales (see e.g., Semenov et al. 2017). This would then be done at the expense of severe modifications of the galaxies during this timespan, and would then void our approach of controlled experiments. Our restart duration of 30 Myr is thus a compromise between not allowing the modifications to induce more side-effects altering the models, and letting the modifications effectively change the star-forming clouds.

A.2. Gas density PDF

Figure A.1 shows the gas density PDFs of the 15 simulations, split for clarity into variations of the SFE (top panel), and of the feedback recipe (bottom panel). For any of the variations considered, gas rich cases (ℱ25 in red and ℱ40 in black) still have slightly larger fractions of dense gas than their ℱ10 equivalent. This difference between high and low f_gas amplifies when reducing the SFE, as the dense gas is consumed at a slower pace. This is however still significantly weaker than the excess of dense gas responsible for the starburst activity in galaxy mergers (Renaud et al. 2014, 2022; Segovia Otero et al. 2022). This increase of the fraction of dense gas is accompanied by a mild deficit of gas at ∼0.1 − 10 cm⁻³. Separating the atomic and molecular phases as in Paper I reveals that this deficit is found in the atomic outer envelopes of the clumps, which then fuels the molecular centers of the star-forming regions.

Fig. A.1. Gas density probability distribution function in all our 15 simulations. The color codes the gas fraction (10% in blue, 25% in red, 40% in black), and the line-style shows the modifications from the fiducial cases (solid lines, identical to Figure 3 of Paper I) as indicated in the legends.

Decreasing the SFE from 30% to 1% per free-fall time boosts the maximum density by about 2.5 dex. We note that this increase of the maximum density corresponds to a reduction of the free-fall time by a factor of ≈20, which is not enough to fully compensate for the division of the SFE by 30. This indicates limits in the self-regulation of the star formation activity by feedback caused by changes in the ISM density structure.

No strong differences are found when varying the feedback model, neither in the shape nor maximum density. Contrarily to the previous case, this hints towards efficient self-regulation of the star formation process at these time and spatial scales.

In summary, varying the parameters and recipes for star formation and feedback does influence the shape of the gas density PDFs, especially for a low SFE. However, the moderate but significant dependences on the gas fraction found in the fiducial runs (Paper I) are maintained in all cases: the gas-rich galaxies always yield a slightly higher fraction of dense gas than their ℱ10 counterparts, but not a significantly higher maximum density, regardless of the sub-grid recipe used. Gas-rich clumpy galaxies do not contain denser gas, but only a slightly higher fraction of dense gas than contemporary disks, and this trend is largely independent of the star formation efficiency and feedback recipes.

A.3. Clump mass function

Figure A.2 shows the clump mass function of the 15 simulations, again split into variations of the SFE (top panel), and of the feedback recipe (bottom panel). As for the PDF, our variations in the subgrid recipes, specially with the feedback models, have a weaker influence than the change of the gas fraction. In particular, the models at low-gas fraction remain well-separated than all the others, regardless of the variations of the subgrid models.

Fig. A.2. Clump mass function, as in Figure 4, in all our 15 simulations. The colors and line styles are as in Figure A.1. The shaded areas indicate the amplitude of the modifications induced by the change of subgrid models. At the high mass-end, the absence of significant overlap between these shaded areas shows that the gas fraction has a stronger imprint on the distribution of dense clumps than the SFE and the feedback model.

Andersson et al. (2024) also varied the subgrid recipes for feedback in simulations of dwarf galaxies. They found that the shape of the cluster mass function strongly depends on the timing of stellar feedback, and its effect in stopping the star formation activity in individual clumps. In our analysis, we examine all gas clumps, regardless of their stellar contents. This necessarily implies that our clumps are at different stages of their evolution, including after the onset of star formation. Therefore, the conclusions of Andersson et al. (2024) are not directly translatable in our framework. The mapping of the clump into the cluster mass function is a complex process, subject to both internal alteration from feedback, and external disruptions. As discussed in Section 3.1, at least the latter vary with the gas fraction. Section 3.3.1 complement this analysis by examining the density structure of the clumps, which provides an insight on their resistivity to feedback. Yet, more detailed studies at higher resolution are needed to reach definitive conclusions on the variations with the gas fraction of the transformation of gas clumps into star clusters.

Appendix B: Examples of individual density profiles

Figure B.1 shows the density profiles centered on the center of mass for the three most massive clumps of each simulation. The presence of sub-clumps within the same structure shifts the center of mass away from a density peak, such that these density profiles are not monotonically decreasing with radius. When averaging over several clumps, the various radial shifts introduce a statistical flattening of the resulting profile, which is not representative of the structure of any clump.

Fig. B.1. Density profiles of the three most massive clumps in each simulation, centered on the center of mass. This illustrates how the diversity of position of density peaks within clumps affects the average profiles shown in Figure 5.

However, some of the clumps shown here do yield a flattening of their density profile in its innermost part. The size of this feature is close to our resolution limit (12 pc) and thus subject to caution. Yet, we note that such flat parts extend to different radii from clump to clump, some being captured by up to 6 resolution elements, and do not show rapid variations of the order of the cell size. This suggests that such an inflexion of the density profile could still be detected at higher resolution, and thus be a physical property of massive clumps. More detailed simulations are needed to reach a definitive conclusion.

Appendix C: Correlation between star formation, the clump mass, and the excess velocity dispersion

Figure C.1 is the equivalent of Figure 8, but replacing the excess mass with the mass of the clump. The correlations between the mass and the quantities related to star formation are much weaker (smaller absolute values of the Spearman coefficient) and less significant (higher p-values) than when using the excess mass. Therefore, the mass of gas clumps is a weaker driver of the star formation activity (rate, speed, and efficiency) than the excess of mass with respect to the average clumps of the same size.

Fig. C.1. Relation between the star formation rate, depletion time, star formation efficiency, and the gas mass of clumps. The number in brackets indicate the correlation coefficient and p-values of the same relations but for the excess mass, copied from Figure 8.

For completeness, Figure C.2 show the relation between the indicator of the star formation activity and the excess velocity dispersion. As hinted in Figure 7, the correlations are very weak (low Spearman coefficient) and not significant (high p-value).

Fig. C.2. Same as Figure 8, but with the excess velocity dispersion Δσ on the horizontal axis. In all panels, the correlations are much weaker and less significant than when considering the excess mass.

All Figures

	Fig. 1. Surface density of gas in the three simulations, at the instant considered in the entire analysis. Red circles indicate the position and size of the gas clumps detected (see text).
In the text

Fig. 2. Example of a massive clump hosting sub-structures, as a zoom-in on Fig. 1 (right panel, fgas = 40%). The contours indicate the surface densities of molecular gas at 10, 100 and 1000 M⊙ pc−2. The red circles shows our clump identification. As visible in this example, some of the largest clumps host clear sub-structures, which are inexistent in smaller clouds (either physically or because of limited numerical resolution).

In the text

Fig. 3. Average Mach number (star symbols), and intensity of tidal (circles) and shear (triangles) forces in gas density bins, normalized to that in the fgas = 40% case. Empty symbols denote the forces normalized to the local self-gravity. All quantities are measured at a scale of 120 pc. Symbols have been slightly offset horizontally for clarity. The approximate distinctions between three different media are identified at the bottom, for reference. The two most gas-rich cases (red) yield comparable values at all densities for all these quantities, significantly higher than in the low gas fraction galaxy (blue).

In the text

Fig. 4. Clump mass functions for the three cases (using staking from 5 consecutive snapshots, to improve the statistics). The dashed lines indicate the best-fit of the high-mass end using a power-law. This gives indexes of −2.0, −1.5 and −1.6 for the ℱ10, ℱ25, and ℱ40 cases respectively. The dotted curves are fits to the data of fgas = 40% using a Schechter function (see text). This functional form does not provide a good match to the other cases.

In the text

Fig. 5. Median density profile of the five most massive gas clumps in each simulation (excluding the central 2 kpc of the galaxies, where radial inflows at galactic scale play a non-negligible role in the accumulation of gas). The top panel shows the radial profiles with the maximum density as the center, and the bottom panel those centered on the center of mass. The shaded areas show the robust standard deviation over these 5 clumps. The dashed line represents a power-law of index −2, for reference.

In the text

	Fig. 6. Median density profile of gas clumps in each simulation (as in Fig. 5), but only considering clumps in the same mass interval. Clumps of comparable masses have similar density profiles regardless of the gas-fraction of their host galaxy.
In the text

Fig. 7. Mass–size (top) and velocity dispersion-size (bottom) relations for the clumps in our simulations. The color of the points encodes the logarithm of the third quantity (velocity dispersion on the top panels, and mass in the bottom row), increasing from dark to light. The solid lines are the least absolute deviation fittings of the linear relations in log-log space. They are translated into a power-law in linear-linear space, with the indexes indicated at the top of each panel. The dashed lines show constant densities, for reference, and the dotted lines mark the threshold for the robust identification of outliers (see text, Eq. (1)).

In the text

Fig. 8. Relation between the star formation rate (top), depletion time (middle), and star formation efficiency (bottom), and the excess of mass of clumps with respect to the best-fit mass–size relations of Fig. 7. The vertical dashed lines mark ΔM = 0, i.e., a clump sitting on the best fit relation shown in Fig. 7. The vertical dotted lines show the thresholds for the robust definition of outliers (see text). Normalized marginal distributions are shown on the right column.

In the text

	Fig. 9. Distributions of average density (total mass divided by total volume) and excess mass of the clumps, with color indicating the depletion time. Correlations between these three quantities are only found for the extreme clumps at high ΔM.
In the text

Fig. 10. Median density radial profiles of clumps of radius between 2 × 107 M⊙ and 3 × 107 M⊙, grouped in three bins of excess mass ΔM. The numbers in square brackets indicate the slope of each curve measured in log–log space in the inner 100 pc. (The conclusions do not change when considering larger radial ranges.) Higher excess masses correspond to more centrally concentrated clumps (steeper profiles).

In the text

	Fig. A.1. Gas density probability distribution function in all our 15 simulations. The color codes the gas fraction (10% in blue, 25% in red, 40% in black), and the line-style shows the modifications from the fiducial cases (solid lines, identical to Figure 3 of Paper I) as indicated in the legends.
In the text

Fig. A.2. Clump mass function, as in Figure 4, in all our 15 simulations. The colors and line styles are as in Figure A.1. The shaded areas indicate the amplitude of the modifications induced by the change of subgrid models. At the high mass-end, the absence of significant overlap between these shaded areas shows that the gas fraction has a stronger imprint on the distribution of dense clumps than the SFE and the feedback model.

In the text

	Fig. B.1. Density profiles of the three most massive clumps in each simulation, centered on the center of mass. This illustrates how the diversity of position of density peaks within clumps affects the average profiles shown in Figure 5.
In the text

	Fig. C.1. Relation between the star formation rate, depletion time, star formation efficiency, and the gas mass of clumps. The number in brackets indicate the correlation coefficient and p-values of the same relations but for the excess mass, copied from Figure 8.
In the text

	Fig. C.2. Same as Figure 8, but with the excess velocity dispersion Δσ on the horizontal axis. In all panels, the correlations are much weaker and less significant than when considering the excess mass.
In the text