The Interplay between the Initial Mass Function and Star Formation Efficiency through Radiative Feedback at High Stellar Surface Densities

The observed rest-UV luminosity function at cosmic dawn (z ∼ 8–14) measured by JWST revealed an excess of UV-luminous galaxies relative to many prelaunch theoretical predictions. A high star formation efficiency (SFE) and a top-heavy initial mass function (IMF) are among the mechanisms proposed for explaining this excess. Although a top-heavy IMF has been proposed for its ability to increase the light-to-mass ratio (ΨUV), the resulting enhanced radiative pressure from young stars could decrease the SFE, potentially driving galaxy luminosities back down. In this Letter, we use idealized radiation hydrodynamic simulations of star cluster formation to explore the effects of a top-heavy IMF on the SFE of clouds typical of the high-pressure conditions found at these redshifts. We find that the SFE in star clusters with solar-neighborhood-like dust abundance decreases with increasingly top-heavy IMFs—by ∼20% for an increase of a factor of 4 in ΨUV and by 50% for a factor of ∼10 in ΨUV. However, we find that an expected decrease in the dust-to-gas ratio (∼0.01 × solar) at these redshifts can completely compensate for the enhanced light output. This leads to a (cloud-scale; ∼10 pc) SFE that is ≳70% even for a factor of 10 increase in ΨUV, implying that highly efficient star formation is unavoidable for high surface density and low-metallicity conditions. Our results suggest that a top-heavy IMF, if present, likely coexists with efficient star formation in these galaxies.


Introduction
Presupernova feedback via radiation, jets, and winds emitted by young stars has been recognized to play a pivotal role in regulating star formation and dictating the life cycle of giant molecular clouds (GMCs) in galaxies (Chevance et al. 2023;Burkhart et al. 2024;Jeffreson et al. 2024).This feedback disrupts GMCs in order ∼unity dynamical timescales via the energy and momentum they impart (Krumholz & Matzner 2009;Fall et al. 2010;Thompson & Krumholz 2016) and drives turbulent motions that could further provide support against collapse (e.g., Mac Low & Klessen 2004;Krumholz et al. 2006;Federrath et al. 2010b;Gallegos-Garcia et al. 2020;Menon et al. 2020Menon et al. , 2021;;Appel et al. 2022).Numerical simulations have demonstrated that this limits the integrated star formation efficiency (SFE; ò * = M * /M gas )-defined as ratio of the final stellar mass M * formed to the available gas mass in the parent molecular cloud M gas -to values of 10% in environments typical of the local Universe (Geen et al. 2016;Kim et al. 2016Kim et al. , 2018Kim et al. , 2021;;Raskutti et al. 2016;Burkhart 2018;Grudić et al. 2018Grudić et al. , 2022;;He et al. 2019;Fukushima & Yajima 2021;Lancaster et al. 2021a).
However, it has become increasingly evident that this is not the case for GMCs typical of high interstellar medium (ISM) pressure environments (P/k B  10 8 K cm −3 ), for which both models (Fall et al. 2010;Thompson & Krumholz 2016) and numerical simulations (Grudić et al. 2018;Fukushima & Yajima 2021;Lancaster et al. 2021a;Menon et al. 2022aMenon et al. , 2023;;Polak et al. 2023) suggest efficiencies ò *  80% because the energy/momentum deposition rate of feedback in this regime is unable to counteract gravity.Such pressures correspond to GMCs with surface densities (Σ  Σ crit = 10 3 M ☉ pc −2 )-2-3 orders of magnitude higher than typical of GMCs in the local Universe-which are the likely sites of so-called super star cluster formation (e.g., McCrady et al. 2005;Portegies Zwart et al. 2010;Turner et al. 2015;Smith et al. 2020); observational estimates seem to be consistent with a high value of ò * for these conditions (Turner et al. 2017;Emig et al. 2020;Rico-Villas et al. 2020;Smith et al. 2020;Costa et al. 2021;He et al. 2022;McKinney et al. 2023;Sun et al. 2024).The environments that host these conditions are relatively rare in the local Universelimited to scenarios such as nuclear starbursts (e.g., Leroy et al. 2018;Emig et al. 2020;Levy et al. 2021), merging luminous infrared (IR)-bright galaxies (e.g., Johnson et al. 2015;Finn et al. 2019;Inami et al. 2022), and localized starbursts in dwarf galaxies (e.g., Ochsendorf et al. 2017;Oey et al. 2017; Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence.Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.Turner et al. 2017).On the other hand, the higher densities, gas fractions, merger rates, and accretion rates of galaxies at higher redshift suggest that high-pressure conditions are more commonly realized at these epochs; indeed, conditions observed in dusty starburst galaxies (Casey et al. 2014), prequiescent massive compact galaxies (Diamond-Stanic et al. 2012;Rupke et al. 2019), and proto-globular-cluster candidates resolved via gravitational lensing (Vanzella et al. 2022a(Vanzella et al. , 2022b;;Pascale et al. 2023) reflect these conditions.
It is therefore timely in this context that JWST has revealed that this dense, clumpy, and compact mode of star formation may well be ubiquitous in the reionization era through the discovery of extremely blue, ultraviolet (UV)-luminous, compact galaxies at redshifts z  10 ( Casey et al. 2024;Finkelstein et al. 2023a;Harikane et al. 2023;Morishita et al. 2024;Robertson et al. 2023;McLeod et al. 2024).The observed sizes (0.5 kpc) of these objects indicate stellar surface densities (Σ  10 4 -10 5 M ☉ pc −2 ) that are comparable to or possibly somewhat higher than those in local super star clusters (see, e.g., Figure 6 of Casey et al. 2024).Highly magnified regions of lensed fields reveal systems at z ∼ 8-10 that are comprised of multiple dense, intensely star-forming clusters, possibly representing the formation sites of presentday globular clusters (Adamo et al. 2024;Mowla et al. 2024).The observed numbers of these bright z  10 galaxies are in excess of the predictions of nearly all prelaunch models of galaxy formation, including both semianalytic models (SAMs) and numerical hydrodynamic simulations (e.g., Dayal et al. 2017;Kannan et al. 2023Kannan et al. , 2022;;Wilkins et al. 2023;Hassan et al. 2023;Yung et al. 2024).Moreover, models almost uniformly predict a much more rapid evolution of the number density of bright galaxies with redshift at these early epochs than the observations indicate (Finkelstein et al. 2023b).Some of the proposed solutions to this tension allude to possibly distinct conditions in star-forming clouds in and around these galaxies, resulting in higher SFE and/or weaker stellar feedback (e.g., Williams et al. 2024;Yung et al. 2024), or to the possibility of a top-heavy initial mass function (IMF), which could lead to higher light-to-mass ratios (Inayoshi et al. 2022;Harikane et al. 2023;Yung et al. 2024).
For example, the Feedback-free Burst Model (FFB; Dekel et al. 2023) posits that when both the gas density9 and surface density in star-forming clouds are high enough and the metallicity is low but not negligible (Z ∼ 0.01-0.1 Z ☉ ), star formation occurs in a burst over a freefall time of ∼1 Myr prior to the onset of supernova feedback and with only weak effects from stellar winds and radiative feedback.This leads to globally efficient star formation in z  10 galaxies, many of which are expected to satisfy these conditions.Li et al. (2023) show that this model produces predictions that are consistent with the JWST observations.
On the other hand, several studies have suggested that the IMF could be top-heavy at these redshifts due to a higher cosmic microwave background temperature (Chon et al. 2022), low metallicities (Sharda & Krumholz 2022;Chon et al. 2024;Sharda et al. 2023), and/or a contribution of Population III stellar populations for which there is general agreement on the possibility of top-heaviness (Larson 1998;Omukai et al. 2005;Harikane et al. 2023;Klessen & Glover 2023).The associated higher UV luminosity per unit mass from a top-heavy IMF could help reconcile the UV luminosity functions without requiring a high SFE (e.g., Inayoshi et al. 2022).Indeed, Yung et al. (2024) show that their fiducial (without changing ò * ) SAM can reproduce the observed UV luminosity function at z ∼ 11 when they increase the UV luminosity-to-mass ratio by a factor of ∼4.
However, ò * and the IMF are not necessarily independent of each other; a top-heavy IMF and the associated increased level of radiative and wind feedback due to a higher fraction of massive stars (that dominate these modes of feedback) is very likely to affect ò * .The metallicity could also affect ò * through its impacts on the dust abundance and cooling physics.Quantifying the interdependence of ò * with the IMF/ metallicity is crucial to shedding light on potential solutions to these surprising findings.This is also relevant in the context of star formation in extreme environments at lower redshifts, where regions of high surface density seem to show possible evidence of top-heavy IMFs (Schneider et al. 2018;Zhang et al. 2018;Upadhyaya et al. 2024).While several previous authors have studied the impact of a top-heavy IMF on star cluster formation (Chon et al. 2024;Fukushima & Yajima 2023), they focused on clouds with mass surface densities and escape speeds that are lower (10 3 M ☉ pc −2 ; v esc  20 km s −1 ) than the extreme cases being found with JWST.These studies have also been done use using numerical methods with less accurate radiation moment closures (Wünsch 2024) and a reduced speed-of-light approach that becomes increasingly computationally expensive at the high optical depths achieved in this regime of surface densities (Skinner & Ostriker 2013).Our goal in this paper is to make use of the more accurate radiative transfer methods developed by Menon et al. (2022b) to explore the effects of a top-heavy IMF and varying dust-to-gas ratio in precisely the conditions that JWST is now probing.
In this paper, we run idealized radiation hydrodynamic numerical simulations of star cluster formation and their radiative feedback with varying levels of UV luminosity-tomass (to emulate differing levels of top-heaviness) to quantify how ò * changes.The paper is organized as follows.In Section 2, we describe the numerical prescriptions we use and the initial conditions of our clouds and the parameter space we explore.In Section 3, we present the evolution of our model clouds and the ò * values we find over our parameter space and discuss the feedback physics driving the trends we find.In Section 4, we provide a discussion on the implications of our results in the context of the JWST results and enumerate the missing physics in our simulations and their possible effects on our outcomes.In Section 5, we conclude with a brief summary of our results.

Simulation Setup
Our simulation setup is very similar to that described in Menon et al. (2023); we briefly summarize the salient features of the setup below and refer the reader that paper for further details.Our simulations represent an isolated cloud of mass M cloud and radius R cloud , which correspond to a mass density of and a mass surface density of . We place our clouds in an ambient medium of density ρ = ρ cloud /100 in pressure equilibrium in a computational domain of size L = 4R cloud .We initialize the fluid with turbulent velocities that follow a power spectrum E(k) ∝ k −2 with a natural mixture of solenoidal and compressive modes for 2 ,6 4, generated with the methods described in Federrath et al. (2010b), and using the implementation of these methods provided in Federrath et al. (2022).We scale the velocity dispersion of the cloud σ v such that our clouds are marginally bound, i.e., α vir = 2, where α vir is given by cloud .We use diode boundary conditions for the gas quantities wherein we permit gas to escape the boundaries but no inflows through them.
We model radiation feedback in two wavelength bands, the UV and the IR; the former is technically a combination of the Lyman continuum (hν 13.6 eV) and far-UV (6.6 eV hν < 13.6 eV) bands, which we collectively refer to as "UV" for simplicity.The only sources of UV radiation are the sink particles that form in our simulations, which represent stellar populations.We adopt a constant UV luminosity-tomass ratio (Ψ UV ) for a given simulation such that the radiative output from a sink of mass M sink is L UV = M sink Ψ UV .On the other hand, the IR emission can come from dust grains that are heated due to the absorption of these UV photons.In addition, we also account for the dust-reprocessed IR field and the associated heating of grains and radiation pressure.We assume that the dust temperature (T d ) is instantaneously equal to the radiation temperature set by the equilibrium between dust emission and UV + IR photon absorption (see Menon et al. 2023 for a justification of this assumption).This assumption might cease to hold true in optically thin conditions; in Section 4.2, we discuss this caveat and argue that it should not affect our results.We use T d to estimate the Planck-Rosseland emission and absorption opacities in the IR using the temperature-dependent Semenov et al. (2003) model; Menon et al. (2022a) show that ignoring this temperature dependence can strongly overestimate the effectiveness of the IR radiation pressure.For the UV, we assume a fixed opacity (identical Planck and Rosseland opacities) of κ UV = 1000 cm 2 g −1 for all our runs,10 consistent with typical estimates of the gray radiation pressure cross section per H atom to blackbody radiation peaking at UV wavelengths (blackbody temperatures of ∼few × 10 4 K; Draine 2011; Kim et al. 2023).These opacities are for Z = Z ☉ ; for other metallicities, we scale our opacities linearly with Z with the underlying assumption that the dust-to-gas ratio scales with metallicity in a linear fashion, which is consistent with observations to zeroth order (e.g., De Vis et al. 2019).It is possible that this assumption overestimates the dust-to-gas ratio at low Z due to the lack of efficient gas-phase accretion (see, e.g., Feldmann 2015; Choban et al. 2022); however, it shall become clear that a more accurate treatment of the metallicity dependence of the dust-to-gas ratio would only reinforce the conclusions we reach below.We initialize our clouds with zero radiation energy/flux in the UV and an IR radiation field corresponding to a dust temperature of T d = 40 K; this is consistent with dust temperatures in observed high-z starburst galaxies (Sommovigo et al. 2022).We adopt Marshak boundary conditions (Marshak 1958) for the radiation with the background value set to match the initial conditions.
We note that we do not include photoionization, stellar winds, protostellar outflows, and magnetic fields in these simulations.In Section 4.2, we discuss (and show in Appendix B) that the former omission would not affect our results, and we discuss the implications of the other missing physics.

Parameter Space
All our clouds have M cloud = 10 6 M ☉ with R cloud = 10 or 3.2 pc to achieve a target Σ cloud = 3.2 × 10 3 M ☉ pc −2 and 3.2 × 10 4 M ☉ pc −2 , respectively.We adopt these values to mimic the high ISM pressure conditions 22 ) expected and now being observed at high redshifts (see Section 1).The lower (higher) Σ cloud value is (approximately) equal (above) the critical surface density beyond which early stellar feedback is unable to regulate the SFE (Fall et al. 2010;Grudić et al. 2018;Lancaster et al. 2021a;Menon et al. 2023), which is a key input in models predicting efficient star formation in galaxies such as the FFB model (Dekel et al. 2023).
For these two Σ cloud cases, we explore variations in Z and Ψ UV , respectively.To mimic increasingly top-heavy IMFs, we explore values of Ψ UV = 1, 4, and 10 times the value for a standard Chabrier IMF (Chabrier 2005) ). Parameterizing the top-heaviness of the IMF with Ψ UV allows us to be agnostic to the degenerate ways in which one can achieve an IMF with an excess of massive stars; regardless, in Appendix A, we present outputs from the Stochastically Lighting Up Galaxies (SLUG) stellar population synthesis code for how these values map to the slope of the high-mass end of the IMF and/or the maximum mass of the star in the stellar population, for the sake of providing intuition.We note that the Ψ UV = 4 case is additionally motivated by empirical estimates of the factor by which the UV luminosities need to be enhanced to reasonably reproduce the bright end of the UV luminosity functions at z ∼ 10 with JWST (Finkelstein et al. 2023b;Yung et al. 2024).We run each of these cases at Z = 10 −2 Z ☉ in addition to solar metallicity to test the effects of the lower metallicities (and implied lower dust-to-gas ratios) expected at high redshifts.We also run a case with Z = 4 Z ☉ and Ψ UV = Ψ Fiducial motivated by possible evidence of supersolar metallicities found in super star clusters in our local Universe (Turner et al. 2015).We summarize our suite of simulations and their parameters in Table 1.

Numerical Methods
We solve the equations of self-gravitating radiation hydrodynamics for all our simulations.We use the FLASH magnetohydrodynamics code (Fryxell et al. 2000;Dubey et al. 2008) for our simulations, with the explicit Godunov method in the split and the five-wave HLL5R (approximate) Riemann solver (Waagan et al. 2011) for the hydrodynamics.The Poisson equation for the self-gravity is solved using a multigrid algorithm implemented in FLASH (Ricker 2008).Sink particles are used to follow the evolution of gas at unresolved scales, the formation of which is triggered when gas properties satisfy a series of conditions to test for collapse and star formation (Federrath et al. 2010a).Gravitational interactions of sink particles with gas and other sinks are considered, and a second-order leapfrog integrator is used to advance the sink particles (Federrath et al. 2010a(Federrath et al. , 2011)).To model the radiative transfer and the associated energy and momentum transfer to gas, we use the Variable Eddington Tensor-closed Transport on Adaptive Meshes (VETTAM) method (Menon et al. 2022b).VETTAM solves the nonrelativistic, angleaveraged, moment equations of radiative transfer in the mixed-frame formulation (Mihalas & Klein 1982), retaining terms that are of leading order in all limiting regimes of radiation hydrodynamics (see, e.g., Krumholz et al. 2007).It uses the Variable Eddington Tensor closure obtained with a time-independent ray-trace solution (Buntemeyer et al. 2016) to close the moment equations; this approach yields much more accurate solutions for problems with multiple radiation sources than any purely local approximation for the Eddington tensor (e.g., the M1 approximation).VETTAM uses an implicit global temporal update for the radiation moment equations for each band, accounting for the coupling between the bands due to dust reprocessing.The radiative output from sink particles is included as a smoothed source term in the moment equations, where we have tested convergence in the smoothing parameters (see Menon et al. 2022a).We used a fixed uniform grid resolution of 256 3 for all our simulations; although VETTAM fully supports adaptive mesh refinement, we chose a fixed modest resolution for simplicity and the computational feasibility required to explore our broad parameter space; we demonstrate convergence in the SFE (within 5%) at these resolutions in similar numerical setups in Menon et al. (2023).
We pause to comment that our numerical model has been used to study the competition between star formation and feedback set by radiation pressure in our previous work (Menon et al. 2022a(Menon et al. , 2023)).These works showed that the integrated SFE increases with the gas surface density of the cloud, producing efficiencies approaching unity for Σ cloud  10 4 M ☉ pc −2 .However, for gas surface densities representative of the local Universe (Σ cloud ∼ 100 M ☉ pc −2 ), we found efficiencies of ∼30%, which is substantially higher than other works in the literature (e.g., Kim et al. 2017).We can confirm that this was simply because this work did not include photoionization, which becomes important in that parameter regime; indeed, in our more recent work with photoionization (S.H. Menon et al. 2024, in preparation), we find efficiencies of ∼10% for clouds in this parameter regime, in strong agreement with other numerical simulations and observed estimates (Chevance et al. 2023).We make this point here to clarify for the reader that our numerical model produces consistent results in regions of parameter space that have been studied widely in previous work.

Competition between Star Formation and Feedback
The initial evolution of all our models is relatively similardensity enhancements due to the turbulent fluctuations undergo self-gravitational collapse and go on to form sink particles (that represent subclusters), which then accrete and continue to increase the total stellar mass (and therefore ò * = M * /M cloud ) in the cloud.Turbulent fluctuations also introduce some nonnegligible mass loss through the computational boundaries at early times (∼10%-20%).This occurs because local gas patches can become unbound and escape through our isolated boundary conditions in the initial phases, although the cloud is globally marginally stable.The stellar mass continues to grow rapidly for t  t ff , after which the evolution between the clouds starts to differ due to the regulating effects of radiative  , where m H is the mass of atomic hydrogen.σ v : turbulent velocity dispersion of the cloud.v esc : escape velocity of the cloud.t ff : freefall time of the cloud.Ψ UV : UV luminosity per unit mass scaled by the value for a Chabrier IMF (Ψ Fiducial ); see Figure 4 for the mapping between the degree of top-heaviness of the IMF and this quantity.Z: metallicity in units of solar metallicity (Z e ); this is also used to scale the dust abundance assuming a linear trend with Z.All our simulations use a resolution of 256 3 .feedback.We can see this in Figure 1, which shows the time evolution of ò * for all our model clouds.The rate of star formation-interpreted from the rate of change of ò * with time -for clouds with progressively higher Ψ UV slows down earlier and more dramatically.This is due to the stronger feedback around the radiating sources that reverses the accretion flow in its vicinity and starts to drive this gas locally outward.We can see this visually in Figure 2, which shows the projected gas density and velocity fields at t = 3t ff for the higher Σ cloud runs.We can see that the sinks are still accreting for Ψ UV = Ψ Fiducial , whereas increasing amounts of gas are outflowing for higher Ψ UV , and the resulting ò * is lower.This can be understood due to the stronger levels of UV radiative feedback for these cases at any given stellar mass.However, consistent with Figure 1, the Z ∼ 0.01 Z ☉ cases seem to show much more modest effects from the feedback even for the higher Ψ UV cases.The stellar masses accumulated in the same time are also higher.This suggests that the lower dust-to-gas ratio in these runs skews the feedback-star formation competition in favor of the latter.We can see that the effects of the dust-to-gas ratio are much less pronounced in the lower Σ cloud case.We explain the reason for this behavior with Z and Σ cloud in Section 3.3.

Integrated SFEs
The aforementioned trends are also reflected in the final saturated level of ò * set by the star formation/feedback balance in our simulations-our key quantity of interest.We calculate this as the value at the point when there is less than 5% of the gas mass remaining in the computational domain, and it is indicated by the termination of the curves in Figure 1.When reporting this quantity, we normalize by its corresponding value for a control run without any feedback to account for the initial gas mass loss due to our isolated turbulent cloud numerical setup, since this gas does not participate in the feedback-star formation competition; we refer to this normalized final SFE as ò *,f .We also do this to place less emphasis on the exact value of ò * in the simulation-since this is expected to vary depending on the turbulent initial conditions-and more on the relative effect of the feedback for a given cloud.We show the values obtained for ò *,f across our simulation suite in Figure 3 as a function of the input value of Ψ UV we use for the stellar populations.As expected, we can see the general trend that ò *,f decreases with increasing Ψ UV .However, ò *,f increases with decreasing Z for a given Ψ UV .The dependence on Z is weak for Σ cloud ∼ 10 3 M ☉ pc −2 (10%).However, for the higher Σ cloud case, the lower dust content can more or less completely counteract the effects of the higher Ψ UV .We can also see that ò *,f for the run with Z ∼ 4 Z ☉ is almost identical to the solar metallicity run; for the higher Σ cloud case, the difference is slightly more evident.This, along with the corresponding trends for lower metallicities, suggests that the dependence on dust content is stronger for the higher Σ cloud case.We also overplot approximate trends with Ψ UV to guide the eye for Z = Z ☉ and Z = 0.01 Z ☉ .We can see that the trend is largely linear (albeit with different slopes) for the lower Σ cloud run, whereas it is clearly nonlinear for the higher Σ cloud at Z = Z ☉ but essentially flat for Z = 0.01 Z ☉ .This implies that for this case, the lower dust-to-gas ratio completely compensates for the (high) increase in the UV luminosity.This has important implications for star cluster formation at high redshifts, as we will discuss below.

Physics Driving Trends
In this section, we briefly explain the feedback physics that drives the trends with Ψ UV and Z in our simulations.The primary feedback mechanism that drives the dynamics in our clouds is the radiation pressure on dust grains-both the single-scattering UV force and the multiple-scattering force due to reemitted IR radiation by warm dust.The former applies a constant force over the absorbing shell of ∼L * /c as long as it is optically thick in the UV; this requires ( )  , which is satisfied across our parameter space even for the Z = 0.01 Z ☉ runs.On the other hand, to be  1).All our simulations undergo rapid collapse (over t  2t ff ) and star formation followed by a saturation in ò * at the point when feedback is able to counteract the collapse; this saturation point is clearly different across our runs. assuming an average κ IR = 5 cm 2 g −1 .Of the simulations we present in this paper, only the two with Σ cloud  10 4 M ☉ pc −2 with Z = Z ☉ and 4 Z ☉ satisfy this condition.These points imply that the IR radiation pressure is only an important contributor for this subset of runs-this is consistent with the findings reported in Menon et al. (2023).In these conditions, the trapped IR radiation field can impart a force ∼f trap L * /c, where f trap describes the trapping factor that quantifies the momentum gained by the multiple scattering of the IR radiation.Now the competition between radiation pressure and gravity can be quantified with the Eddington ratio.The Eddington ratio for a column of gas with surface density Σ exposed to a stellar population with UV luminosity per unit mass Ψ UV for singlescattering radiation pressure is (e.g., Thompson et al. 2015) Since there would be a distribution of Σ surrounding the stellar population as the cloud evolves (Thompson & Krumholz 2016), the above calculation suggests that the fraction of sight lines that become super-Eddington increases with Ψ UV .This could, to zeroth order, explain the trend we see in the Σ cloud = 3000 M ☉ pc −2 runs.Note that there is no dependence on Z for Γ UV , as long as the gas is optically thick in the UV.
On the other hand, for the multiple-scattering IR radiation force, the Eddington ratio is where we explicitly note the (assumed) linear dependence of κ IR -the IR opacity-on the metallicity. 11Note that this does not have a dependence on Σ, as long as the column is optically thick in the IR, i.e., ( ) We can now understand the stronger metallicity dependence for the higher Σ cloud case: it is in the regime where the dust-togas-ratio-dependent IR radiation pressure is the dominant feedback mechanism.This force only plays a relatively minor role 12 for Σ cloud  10 4 M ☉ pc −2 -even less so at lower Z-as it is optically thin in the IR.For the higher Σ cloud case, the IR radiation pressure is clearly the crucial force, as the UV radiation pressure has insufficient momentum to compete with 11 We stress that this expression is valid only if the gas is optically thick in the UV such that the full stellar luminosity gets reprocessed to the IR.In addition, this expression amounts to assuming f trap = κ IR Σ for the IR radiation force; in reality, f trap would depend on the dust temperatures through the column of gas (Menon et al. 2022a) and nonlinear radiation-matter interactions (Krumholz & Thompson 2012), both of which are captured in our simulations.The constant κ IR value we use in this expression is just a simplification we make to explain the qualitative trends we find in our simulations. 12However, this effect is not negligible; some sight lines can become optically thick in the IR due to the turbulent overdensities.This likely explains the ∼10% differences between the Z = Z ☉ and Z = 0.01 Z ☉ runs for this cloud.
gravity at this high Σ, even for Ψ UV = 10Ψ Fiducial .This is reflected in the value of ò *,f for Z = 0.01 Z ☉ , where only the UV acts; an inspection of Equation (2) clearly indicates sub-Eddington conditions for these parameters.The linear dependence on Ψ UV in Equation (2) also explains the linear trend seen with Ψ UV for the lower Σ cloud cases.Finally, we find that the nonlinear trend seen at Z = Z ☉ for the higher Σ cloud runs is due to the subtle effect of more efficient trapping of IR photonsi.e., higher f trap -for clouds with higher Ψ UV .This occurs because the (IR) radiation temperature is significantly higher for higher Ψ UV , which renders κ IR higher, thereby imparting more momentum per unit stellar mass.Connecting with Equation (3), it is the combination of increasing κ IR due to warmer dust along with the linear increase due to Ψ UV that leads to the nonlinear trend.

Implications for Massive Galaxies at Cosmic Noon
Higher global SFEs than the local Universe and a top-heavy IMF are two of several proposed scenarios to reconcile the observed abundance of massive UV-bright galaxies at z ∼ 8-12 with prelaunch model predictions (Inayoshi et al. 2022;Finkelstein et al. 2023b;Harikane et al. 2023;Yung et al. 2024).The former is a key element of the FFB model (Dekel et al. 2023), which invokes high cloud-scale SFE and ineffective feedback by stellar and supernova-driven winds to achieve more efficient galaxy-scale star formation and hence boost the numbers of UV-luminous galaxies at early times.Li et al. (2023) showed that this model is consistent with observations when they adopted a cloud-scale SFE of 50%.On the other hand, a top-heavy IMF and the associated higher UV luminosities could also match the UV luminosity functions while still adopting the lower SFEs that seem to be typical of galaxy populations at z  8 (Tacchella et al. 2018); Yung et al. (2024) show that a boost of ∼4 in the UV luminosity-to-mass ratio can reproduce the UV luminosity function at z ∼ 11 without modifying the SFE or feedback strength.However, these studies explored the impact of the SFE and the IMF as if they were independent, which of course is not the case in reality.In this study, we have quantified how these two quantities depend on each other at the cloud scale.
Our results indicate that as long as clouds have surface densities Σ cloud  10 3 M ☉ pc −2 -a condition that seems to be commonly satisfied at z  10 based on observed galaxy sizes (e.g., Finkelstein et al. 2023b;Casey et al. 2024;Morishita et al. 2024;Adamo et al. 2024)-an SFE significantly higher than that typical of the local Universe (∼10%) is unavoidable even in the presence of a top-heavy IMF (we have investigated cases where the luminosity-to-mass ratio is up to 10 times the typical value; see Figure 4 and Appendix A).A top-heavy IMF results in only a moderate reduction in the SFE, and only if the dust abundance is similar to the solar neighborhood.For metallicities that seem to be typical at the highest redshifts where we have reliable estimates, z ∼ 8-10, i.e., Z ∼ 0.1-0.3Z ☉ (Curti et al. 2023;Nakajima et al. 2023), and assuming a linear relation between dust-to-gas ratio and metallicity, 13 the SFE is higher for a given Ψ UV and completely counteracts the effects of Ψ UV for highly compact clouds (10 4 M ☉ pc −2 ).Moreover, the nature of the IMF for conditions at z ∼ 10 is highly uncertain, with works suggesting that it could even be bottom-heavy (e.g., Conroy & van Dokkum 2012;Tanvir et al. 2022;Tanvir & Krumholz 2024); that being said, these studies probe the IMF in a mass range (1 M ☉ ) that does not contribute to Ψ UV .Even if this scenario were true for the highmass end of the IMF, it would only imply that even more efficient star formation would be required, as Ψ UV would be even lower than with a standard IMF.All of this suggests that highly efficient star formation at the cloud scales may be Figure 3.The final integrated cloud-scale SFE (ò *,f ) obtained in all of simulations-scaled by the value obtained for a run without feedback (ò *,NoFB )-shown as a function of Ψ UV We indicate lines and their slopes to guide the eye (no fitting).We can see that there is a general trend of decreasing ò *,f with Ψ UV ; however, a lower Z can at least partly compensate for this decrease-much more so for the higher Σ cloud cases (right panel).The trend for the higher Σ cloud case is also clearly nonlinear.These trends are likely driven by the differing levels of momentum imparted by radiation pressure on dust across our simulations (see Section 3.3). 13Dust-to-gas ratios are highly uncertain at the masses and redshifts of the z  10 galaxies.Observational constraints on dust-to-gas ratios at z ∼ 2 and metallicities of 12 + log(O/H) ∼ 8.5-8.8 are consistent with values in the nearby Universe.However, for local galaxies, Rémy-Ruyer et al. (2014) find that the relationship between dust-to-gas and metallicity is best fit by a broken double power law, while De Vis et al. (2019) find that it is well fit by a single power law.This could lead to a difference of up to an order of magnitude in dust-to-gas at the typical metallicities (12 + log(O/H) ∼ 7.5) of the JWST galaxies.
ubiquitous in high-redshift galaxies, irrespective of the properties of the stellar populations that populate them.
If we take this at face value, it is possible that the combination of a top-heavy IMF in addition to efficient star formation could overpredict the UV luminosity functions, since they both contribute to an excess at the bright end.However, there are two key subtleties to point out in this context.First, studies that found that a factor of ∼4 increase in Ψ UV is sufficient to reproduce the observations do not account for any possible dust extinction (e.g., Yung et al. 2024).Second, the SFE values and trends we quantify in this study are at the cloud scale (10 pc), whereas the quantity relevant for the luminosity functions is the baryon efficiency ratio defined over the whole galaxy (ò *,gal ).For instance, Li et al. (2023) find that the FFB model fits the tentative JWST data for ò *,gal ∼ 20%.This could either reflect the true SFE values within each star-forming cluster or, alternatively, reflect a duty cycle of star formation in the galaxy.The FFB scenario does predict a duty cycle,14 due to the need to accumulate enough accreted gas for triggering the fragmentation into star-forming clouds, which can lead to ò *,gal ∼ 20% in spite of the assumed SFE ∼ 100% at the scale of the individual clouds (Dekel et al. 2023;Li et al. 2023).The result obtained in the current work of higher SFE within the individual clusters is consistent with this duty-cycle interpretation of the lower SFE when averaged over time in the galaxy.
Alternatively, another possibility to reconcile our high cloudscale SFE values with relatively lower ò *,gal is that only a fraction of the gas in the galaxy participates in star formation in clouds.It is possible that the remaining gas is ejected in outflows by the feedback from older stellar populations, possibly explaining the dust-free nature of these galaxies (Ferrara 2024;Ferrara et al. 2023;Fiore et al. 2023), which is likely critical to simultaneously explain the observed UV luminosity functions and the blue UV continuum slopes (e.g., Cullen et al. 2024)-although see Li et al. (2023) for an alternative explanation of these findings.This possibility raises another effect of a top-heavy IMF that we cannot capture in our simulations: a top-heavy IMF15 would lead to more energy and mass-loaded winds, potentially further decreasing ò *,gal , such that its combination with a top-heavy IMF could be consistent with the observations.There is scope for studying the interaction of these three key parameters-SFE at the cloud scale, IMF, and the feedback effects in driving galaxy-scale winds-taking into account their respective dependencies on each other.By combining these joint constraints into a galaxy-scale SAM (Somerville & Davé 2015), we may be able to constrain the regions of parameter space permitted by the observations.

Missing Physics and Possible Implications
It is important to note that several physical mechanisms are missing in our numerical simulations; we list them here and discuss how they might affect the outcomes.
We only model the radiative feedback on dust and do not include photoionization and therefore the momentum from the associated thermal pressure of ionized gas.However, we argue that this would make little difference to the outcome of our simulations, as the clouds we model have escape speeds v esc  ∼2-5c s,ion , where c s,ion ∼ 10 km s −1 is the ionized gas sound speed. 16These arguments are consistent with results presented in models and Figure 4.The UV luminosity per unit mass of stars: (i) for a given slope at the high-mass end of the IMF (α; left panel) and (ii) for a fixed Salpeter slope but different upper stellar mass limits (M max ), both normalized by its counterpart for a standard Chabrier IMF-i.e., a Salpeter slope (i.e., α = -2.35)and an upper mass limit of 120 M ☉ .These values have been obtained with the SLUG stellar population synthesis code.The peak mass of the IMF and its shape are kept identical in this calculationonly the high-mass slope and the upper mass limit are varied independently.We use this range of values (Ψ UV ∼ 1-10 Ψ Salpeter ) to motivate the parameter space we explore.
previous numerical simulations that show that radiation pressure on dust is the dominant radiative feedback mechanism in this regime for regulating star formation (Krumholz & Matzner 2009;Dale et al. 2012;Kim et al. 2016Kim et al. , 2018)).We also demonstrate that this is the case in Appendix B by rerunning one of our models with the effects of photoionization included.We can justify the omission of protostellar outflows along similar lines; v esc is much higher than the ∼1 km s −1 threshold suggested by Matzner & Jumper (2015) for effective gas ejection by jets.
We also do not model stellar winds.While this might at first seem like a major omission, we note that the effectiveness of stellar wind feedback has been shown to be reduced compared to analytic estimates due to efficient cooling at turbulent interfaces in the multiphase gas, rendering it momentumlimited (Lancaster et al. 2021b(Lancaster et al. , 2021c)); this has been shown to be especially true in the regime of high Σ cloud  10 3 M ☉ pc −2 (Lancaster et al. 2021a) that we focus on here.However, this still implies that there would be a force  p w -where is the wind momentum for a mass-loss rate  M w and wind velocity v w -acting on the gas. p w is expected to be ∼L * /c for a stellar population (see Figure 3 in Lancaster et al. 2021b), suggesting that this should induce an order unity correction to our obtained values of ò *,f .In other words, it is possible that the ò *,f obtained at Ψ UV = 10 would be obtained for Ψ UV = 5 with the additional effect of stellar winds.That being said, it is highly likely that the two feedback mechanisms do not interact in a simple additive fashion.Moreover, at metallicities Z  0.1 Z ☉ , winds from massive O stars are considerably weaker (Leitherer et al. 1992;Vink et al. 2001), meaning that significant stellar wind feedback is delayed until the onset of Wolf-Rayet winds (Lancaster et al. 2021b;Dekel et al. 2023).This time delay may be too long to have a significant impact on star formation in clouds of these densities.Numerical simulations that combine wind and radiative feedback would provide more formal quantification of the resultant SFEs in such conditions.
We assume perfect coupling between gas and dust temperatures and radiative equilibrium for the radiation and dust temperatures.These assumptions are quite reasonable when our clouds are optically thick in the IR; however, for our simulations with Z ∼ 0.01 Z ☉ , they start to break down.For instance, dust and gas temperatures likely decouple in these conditions except in very high-density regions (n  10 7 cm −3 ), which renders our estimates for the gas temperatures incorrect.However, the dynamical impact of this would be minor, since the thermal pressure is not a significant force in the systems we are investigating.The fragmentation properties in our simulations would be affected by this error, but we do not resolve individual stars anyway, and our scope is limited to studying the net competition between radiation forces and gravity in clouds, which is unlikely to be affected.
Our assumption that the dust and the radiation field are in LTE also starts to break down at low dust-to-gas ratios when the dust becomes optically thin in the IR.In this limit, the color temperature of the IR radiation field at any spatial location is not equal to the local dust temperature-an effect that can only be captured by a numerical method that models the evolution of the full spectral energy distribution through the cloud.The way in which this assumption directly affects our numerical model is that our estimated dust temperatures would be incorrect in the optically thin limit, directly affecting the IR dust absorption opacities, which then subsequently affects the IR component of the radiation force (which is ∝ opacity).However, we estimate that the impact of this would be minor, since this error only applies in the limit where the dust is optically thin in the IR, in which case we are in the single-scattering regime anyway, and the IR radiation force is negligible; the latter becomes important only when optically thick in the IR, in which case our assumption is valid.One might question if the very statement that the cloud is optically thick/thin in the IR might itself be affected by the (indirect) error we make in the dust opacity.We estimate that this is unlikely, since the range of (gray) IR opacities for dust warmer than 40 K varies by at most a factor of ∼a few (see Figure 1 in Menon et al. 2022a).Therefore, even if we assume significant error in the dust temperature (which is itself unlikely17 ), it results in the IR dust opacity being underestimated by a factor of ∼a few, which is insufficient to alter the regime of the problem from the singlescattering to the multiple-scattering regime for our Z = 0.01 Z ☉ clouds, which are optically thin by at least 1-2 orders of magnitude.Hence, the impact of this assumption on our results is unlikely to be significant.That being said, this is a subtle effect that could affect systems that are marginally optically thick; there is scope for future (frequency-dependent) calculations to quantify the impact of this assumption in such conditions.
We also do not include magnetic fields, which could provide additional support against gravitational collapse and therefore possibly render higher fractions of gas unbound (e.g., Burkhart 2018;Krumholz & Federrath 2019;Kim et al. 2021).In addition, we do not have the influence of an external larger-scale turbulent environment, which could provide additional stabilization (Kim et al. 2021;Orr et al. 2022;Forbes et al. 2023) through a turbulent cascade acting on the scales of our clouds but also possibly additional compressive modes (Appel et al. 2023).Both of these could slightly affect our obtained values of ò *,f .We therefore urge caution in interpreting the exact values of ò *,f we report.We emphasize that our main takeaway is the trends we find with the IMF and dust content (metallicity).

Conclusions
We study the efficiency of star formation set by radiative feedback for assumed IMFs that are (increasingly) top-heavy and at different dust-to-gas ratios (or metallicity Z, assuming a linear relation between the two).We focus on massive, dense, compact clouds with initial gas surface densities Σ cloud  10 3 M ☉ pc −2 , which are likely typical for galaxies that have been detected by JWST at z ∼ 10.Past theoretical studies have shown that clouds in this regime are expected to exhibit very high SFEs (e.g., Lancaster et al. 2021a;Menon et al. 2023;Polak et al. 2023) for a standard UV luminosity-tomass ratio assuming a Chabrier (2005) IMF (Ψ Fiducial ).We test the effects of increased feedback due to a top-heavy IMF on such clouds by assuming different values of the UV luminosity-to-mass ratios (Ψ UV )-up to 10Ψ Fiducial -for the stellar populations forming in our simulations.We also explore the effects of sub-and supersolar metallicities to mimic

Figure 1 .
Figure1.The SFE as a function of time for our simulations.Ψ UV = 1, 4, and 10 Ψ Fiducial are represented by dashed-dotted, dashed, and solid lines, respectively; colors indicate the different metallicities we explore, and the two panels represent the different Σ cloud values in our simulation suite (Table1).All our simulations undergo rapid collapse (over t  2t ff ) and star formation followed by a saturation in ò * at the point when feedback is able to counteract the collapse; this saturation point is clearly different across our runs.

Figure 2 .
Figure 2. Gas surface density distributions at a time t = 4t ff for the Σ cloud = 3.2 × 10 4 M ☉ pc −2 runs with increasing Ψ UV (left to right) at Z = Z e (top) and Z = 0.01 Z e (bottom).We can see the general trend of a stronger impact of feedback at higher Ψ UV and Z, as evidenced by the presence/absence of outflows in the velocity distribution, which are driven by radiation pressure on dust.The SFE (ò * ; annotated at the top left of each panel) is high in all these cases, except when both Ψ UV and Z are high.Comparison of the achieved ò * in the top right and bottom right panels demonstrates that a drop in Z can counteract the effects of a higher Ψ UV .

Table 1
Summary of Our Simulation Suite and Their Initial Condition Parameters Notes.The columns are described as follows.Σ cloud : mass surface density of the cloud given by .M cloud : mass of the cloud.R cloud : radius of the cloud.n cloud : number density of the cloud given by