starforge

While recent surveys have resolved hundreds of nearby protostellar discs, numerical simulations assuming ideal magnetohydrodynamics (MHD) have historically struggled to achieve disc formation due to efficient angular momentum removal by magnetic torques. Non-ideal MHD effects, relevant at the low ionization fractions typical of molecular clouds, have been shown to reduce the effectiveness of magnetic braking and promote disc formation. In this work, we present the results from a suite of calculations following the gravitational collapse of 50 $\mathrm{M}_{\odot }$ turbulent molecular cloud cores down to the formation and evolution of stellar systems and protostellar discs. We use the radiation-MHD code GIZMO including non-ideal MHD (Ohmic resistivity, ambipolar diffusion, and the Hall effect) and the STARFORGE (STAR Formation in Gaseous Environments) numerical framework for modelling star formation and stellar feedback. We compare the effects of assuming ideal versus non-ideal MHD and including sub-grid protostellar jet feedback on disc formation and evolution. Discs form in all of our models but are least massive in the model with ideal MHD and sub-grid jet feedback. Apart from the ideal MHD $+$ jets model, we do not observe any significant differences in disc properties between the ideal and non-ideal MHD models; however, ideal MHD discs are embedded in smaller rotating envelopes. Disc sizes are in general agreement with those of observed discs. Jet feedback increases core fragmentation and reduces final stellar masses. Our results suggest that magnetic braking does not efficiently suppress disc formation, regardless of whether ideal or non-ideal MHD is assumed, under the dynamical conditions in which multiple stellar systems form.

Star formation occurs within dense regions of giant molecular clouds (GMCs); however, exactly how gas collects and evolves to form individual stars and what role dense cores play remains unclear. We use the Lagrangian cell information in the STARFORGE simulation suite to track star-forming gas in three GMCs with varying magnetic field strengths. We find that once a protostar forms, the lifetime of the unaccreted gas correlates with the final stellar mass, where low-mass stars (M_* < 0.5M_⊙) accrete for 0.5─0.6 Myr from a relatively local reservoir of gas and high-mass stars (M_* > 2M_⊙) accrete over 3.3─4.7 Myr from a much larger volume. Although the protostellar accretion time increases weakly with magnetic field strength, the accreting gas radii, velocity dispersions, virial parameters, and magnetic energy ratios are largely insensitive to the global cloud properties. At the time of protostar formation, the unaccreted gas exhibits linewidth-size and mass-size relations characteristic of turbulently regulated, isothermal dense cores, following σ_v ∝ R^0.47−0.55 and M ∝ R^1.0−1.1, respectively. Low- and intermediate-mass stars undergo relatively continuous accretion, and their accretion histories are well-fit by isothermal sphere, turbulent core, or competitive accretion models, where no one model fits all masses. However, many high-mass stars experience intermittent accretion, and their accretion histories are not well-fit by any of these models. While the distribution of accreting gas is more extended than typically defined dense cores, the physical properties and structure of the star-forming gas resemble those of observed cores and are largely regulated by turbulence and feedback.

Cosmic rays (CRs) drive ionization and influence gas dynamics in molecular clouds (MCs), potentially impacting the resulting star formation outcomes. Although previous simulations of individual star formation have included methods for CR transport (CRT), none have been large enough to resolve the stellar initial mass function (IMF). We conduct numerical simulations following the collapse of a 20,000 M_⊙ MC and the subsequent star formation including CRT, both with and without CRs accelerated by winds from the young massive stars, and compare against a non-CRT simulation. We show that after the first massive stars form, the cavity produced by feedback is more pronounced in the CRT simulations because the external CRs are able to propagate inward and compress the gas into higher-density structures. This increases the subsequent star formation in the cloud; by the end of the simulation, the star formation efficiency (SFE) in the CRT simulation including stellar wind CRs is 43% higher than the non-CRT simulation. The IMF is also top-heavy in comparison, with a slope above 1 M_⊙ that is shallower by ∼20%. These effects are also present in the simulation without wind-accelerated CRs, but they are not as pronounced; the SFE is only 16% higher than the non-CRT simulation, and the IMF high-mass slope is shallower by ∼10%. These results may explain some of the observed top-heavy IMFs, which typically occur in high-CR environments such as the Galactic center. *This manuscript has been authored in part by UT-Battelle, LLC, under contract DE-AC05-00OR22725 with the US Department of Energy (DOE). The US government retains and the publisher, by accepting the article for publication, acknowledges that the US government retains a nonexclusive, paid-up, irrevocable, worldwide license to publish or reproduce the published form of this manuscript, or allow others to do so, for US government purposes. DOE will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan (http://energy.gov/downloads/doe-public-access-plan).

Star formation occurs within dusty molecular clouds that are then disrupted by stellar feedback. However, the timing and physical mechanisms that govern the transition from deeply embedded to exposed stars remain uncertain. Using the STARFORGE simulations, we analyze the evolution of "embeddedness," identifying what drives emergence. We find the transition from embedded to exposed is fast for individual stars, within 1.3 Myr after the star reaches its maximum mass. This rapid transition is dominated by massive stars, which accrete while remaining highly obscured until their feedback eventually balances, then overcomes, the local accretion. For these massive stars, their maximum mass is reached simultaneously with their emergence. Once these stars are revealed, their localized, pre-supernova feedback then impacts the cloud, driving gas clearance. Because massive stars dominate the luminosity, their fast, local evolution dominates the light emergence from the dust. We calculate the dependence of these processes on the mass of the cloud and find that emergence always depends on when massive stars form, which scales with the cloud's free-fall time. We also measure the evolution of dust and Hα luminosities, where for ∼2 Myr, these tracers outshine the emerging stellar continuum, reaching their peak when gas and dust remain tightly coupled to the massive stars. These results closely resemble observationally observed lifetimes, tying the observable dust and line emission directly to the same localized processes that drive stellar emergence, evidence that our simulated de-embedding physics is representative of real star-forming regions. Thus, because the initial embedding of the most luminous stars is highly local, the emergence of stars is a faster, earlier, more local event than the overall disruption of the cloud by gas expulsion.

We present a cosmological zoom-in radiation magneto-hydrodynamic (RMHD) simulation, using FORGE'd in FIRE, that follows the formation, growth, and evolution of a single metal-free Pop. III (proto)star at redshift $z \sim 14$. The simulation captures a rotationally supported circumstellar disk and protostellar jets, both resolved down to $<100$ au scales. We find the star grows to $\sim 27$ M$_{\odot}$ over $31,000$ years, with its final mass regulated by accretion and protostellar jets. Protostellar jets form because the magnetic mass-to-flux ratio lies within the regime that allows jet launching, and they are further enabled by a rotating circumstellar disk with sufficient gas-magnetic-field coupling, both present in this simulation. These jets regulate accretion onto the (proto)star and drive outflows that collide with infalling gas, slowing inflow at large radii due to the substantial momentum they carry. A circumstellar disk forms, extending out to $\sim 0.01$ pc, which remains gravitationally stable (Q $\gg 1$). The stability of the disk is maintained through both thermal support and turbulence. In this paper we focus on how jets play a critical role not only in shaping the final masses of Pop. III stars but also in directly influencing their surroundings by regulating accretion. These results will provide important insights into the initial mass function and feedback processes in the earliest star-forming regions of the Universe.

We systematically investigate how cloud-cloud collisions influence star formation, emphasizing the roles of collision velocity, magnetic field orientation, and radiative feedback. Using the first cloud-cloud collision simulations that model individual star formation and accretion with all stellar feedback mechanisms, we explore the morphological evolution, star formation efficiency (SFE), fragmentation, stellar mass distribution, and feedback-driven gas dispersal. Our results show that cloud collisions substantially enhance the rate and timing of star formation compared to isolated scenarios, though the final SFE remains broadly similar across all setups. Lower collision velocities facilitate prolonged gravitational interaction and accumulation of gas, promoting sustained star formation characterized by elongated filamentary structures. Conversely, high-velocity collisions induce rapid gas compression and turbulent motions, leading to intense but transient episodes of star formation, which are curtailed by feedback-driven dispersal. The orientation of the magnetic field markedly affects collision outcomes. Parallel fields allow gas to collapse efficiently along magnetic lines, forming fewer but more massive stars. In contrast, perpendicular fields generate significant magnetic pressure, which stabilizes the shock-compressed gas and delays gravitational collapse, resulting in more distributed and less massive stellar fragments. Radiative feedback from massive stars consistently regulates star formation, halting further gas accretion at moderate efficiencies (10-15%) and initiating feedback-driven dispersal. Although the cloud dynamics vary significantly, the stellar mass function remains robust across scenarios-shaped modestly by magnetic orientation but only weakly influenced by collision velocity.

Most main sequence stars, unlike our Sun, belong to multiple systems containing two or more stars. How and when these multiples come together and become bound is uncertain, as the earliest stages of star formation are difficult to resolve. Here we analyse simulations of star cluster formation in Milky Way-like conditions, including all key physics and stellar feedback mechanisms, to understand how multiple systems form. We show that ~70-80% of binaries are gravitationally bound from the moment the second star forms. Binaries evolve and accrete together, which will affect their planetary systems and chemical evolution. Half of the binaries are disrupted by the end of the star-formation epoch, such that ~40% of the final single stars belonged to a multiple at some point, with implications for the stellar initial mass function. Formation in multiples is the dominant mode of star formation, accounting for at least 57% of stars.

Molecular clouds (MCs) are active sites of star formation in galaxies, and their formation and evolution are largely affected by stellar feedback. This includes outflows and winds from newly formed stars, radiation from young clusters, and supernova explosions. High-resolution molecular line observations allow for the identification of individual star-forming regions and the study of their integrated properties. Moreover, state-of-the-art simulations are now capable of accurately replicating the evolution of MCs, including all key stellar feedback processes. We present ¹³CO(2─1) synthetic observations of the STARFORGE simulations produced using the radiative transfer code RADMC-3D, matching the observational setup of the SEDIGISM survey. From these synthetic observations, we identified the population of MCs using hierarchical clustering and analysed them to provide insights into the interpretation of observed MCs as they evolve. The flux distributions of the post-processed synthetic observations and the properties of the MCs, namely, radius, mass, velocity dispersion, virial parameter, and surface density, are consistent with those of SEDIGISM. Both samples of MCs occupy the same regions in the scaling relation plots; however, the average distributions of MCs at different evolutionary stages do not overlap on the plots. This highlights the reliability of our approach in modelling SEDIGISM and suggests that MCs at different evolutionary stages contribute to the scatter in observed scaling relations. We study the trends in MC properties, morphologies, and fragmentation over time to analyse their physical structure as they form, evolve, and are destroyed. MCs appear as small diffuse cloudlets in early stages, and this is followed by their evolution to filamentary structures before being shaped by stellar feedback into 3D bubbles and getting dispersed. These trends in the observable properties of MCs are consistent with other realisations of simulations and provide strong evidence that clouds exhibit distinct morphologies over the course of their evolution.

Recent studies suggest that most star-forming regions in our Galaxy form stellar associations rather than bound clusters. We analyse models from the STAR FORmation in Gaseous Environments (STARFORGE) simulation suite, a set of magnetohydrodynamical simulations that include all key stellar feedback and radiative processes following star formation through cloud dispersal. We create synthetic observations by introducing observational biases such as random spurious measurements, unresolved binaries, and photometric sensitivity. These biases affect the measurement of the group mass, size, and velocity dispersion, introducing uncertainties of up to 100 per cent, with accuracy improving as the number of system members increases. Furthermore, models favouring the formation of groups around massive stars were the most affected by observational biases, as massive stars contribute a larger fraction of the group mass and are often missing from astrometric surveys like Gaia. We compare the simulations to the Cepheus Far North (CFN) region, and show that CFN groups may have formed in a low-density environment similar to those modelled in STARFORGE but with massive stars not located preferentially in groups. We also question the effectiveness of the kinematic traceback method, showing that it is accurate within 20 per cent only for certain associations with actual virial parameters above 2. However, observational biases can artificially raise the virial parameter by up to a factor 10, making it difficult to evaluate the reliability of the traceback age. Additionally, since stars continue to form during the dispersal of the parent cloud, we find no relation between the stellar-dynamical age difference and the length of the embedded phase.

Among the remarkable strides made by JWST is the discovery of the earliest star clusters found to date. These have been proposed as early progenitors of globular clusters, which are known to come from the early stages of star formation in the Universe. This is an exciting development in modern astronomy, as it offers an opportunity to connect theoretical models of globular cluster formation to actual observations of these high-redshift structures. In this work, we aim to develop observational signatures of a star cluster formation route known as supersonically induced gas objects, which are dark matter─less gas clouds in the early Universe proposed as a potential origin of some globular clusters. For the first time, we follow the star formation process of these early Universe objects using high-resolution hydrodynamical simulations, including mechanical feedback. Our results suggest that the first dark matter─less star clusters are top heavy, meaning that they have a flatter initial mass function slope compared to very young low-metallicity star clusters in the local Universe, and they also have extremely high stellar mass surface densities compared to their local counterparts.

Constraining the physical and chemical evolution of molecular clouds is essential to our understanding of star formation. These investigations often necessitate knowledge of some local representative number density of the gas along the line of sight. However, constraining the number density is a difficult endeavor. Robust constraints on the number density often require line observations of specific molecules along with radiation transfer modeling, which provides densities traced by that specific molecule. Column density maps of molecular clouds are more readily available, with many high-fidelity maps calculated from dust emission and extinction, in particular from surveys conduction with the Herschel Space Observatory. We introduce a new probabilistic model which is based on the assumption that the total hydrogen nuclei column density along a line of sight can be decomposed into a turbulent component and a gravitationally dominated component. Therefore, for each pixel in a column density map, the line of sight was decomposed into characteristic diffuse (dubbed "turbulent") and dense (dubbed "gravitational") gas number densities from column density maps. The method thus exploits a physical model of turbulence to decouple the random turbulent column from gas in dense bound structures empirically using the observed column density maps. We find the model produces reasonable turbulent and gravitational densities in the Taurus L1495/B213 and Polaris Flare clouds. The model can also be used to infer an effective attenuating column density into the cloud, which is useful for astrochemical models of the clouds. We conclude by demonstrating an application of this method by predicting the emission of the [C II] 1900 GHz, [C I] 492 GHz, and CO (J = 1─0) 115 GHz lines across the Taurus L1495/B213 region at the native resolution of the column density map utilizing a grid of photodissociation-region models.

Dense gas in molecular clouds is an important signature of ongoing and future star formation. We identify and track dense cores in the STARFORGE simulations, following the core evolution from birth through dispersal by stellar feedback for typical Milky Way cloud conditions. Only ∼8% of cores host protostars, and most disperse before forming stars. The median starless and protostellar core lifetimes are ∼0.5─0.6 Myr and ∼0.8─1.1 Myr, respectively, where the protostellar phase lasts Myr. While core evolution is stochastic, we find that virial ratios and line widths decline in prestellar cores, coincident with turbulent decay. Collapse occurs over ∼0.1 Myr, once the central density exceeds ≳10⁶ cm⁻³. Starless cores, only, follow line-width─size and mass─size relations, σ ∝ R^0.3 and M ∝ R¹. The core median mass, radius, and velocity dispersion scale weakly with the cloud magnetic field strength. We cluster the core properties and find that protostellar cores have >80% likelihood of belonging to three particular groups that are characterized by high central densities, compact radii, and lower virial parameters. Overall, core evolution appears to be universally set by the interplay of gravity and magnetized turbulence, while stellar feedback dictates protostellar core properties and sets the protostellar phase lifetime.

The dust grain size distribution (GSD) likely varies significantly across star-forming environments in the Universe, but its impact on star formation remains unclear. This ambiguity arises because the GSD interacts nonlinearly with processes like heating, cooling, radiation, and chemistry, which have competing effects and varying environmental dependencies. Processes such as grain coagulation, expected to be efficient in dense star-forming regions, reduce the abundance of small grains and increase that of larger grains. Motivated by this, we investigate the effects of similar GSD variations on the thermochemistry and evolution of giant molecular clouds (GMCs) using magnetohydrodynamic simulations spanning a range of cloud masses and grain sizes, which explicitly incorporate the dynamics of dust grains within the full-physics framework of the STARFORGE project. We find that grain size variations significantly alter GMC thermochemistry: the leading-order effect is that larger grains, under fixed dust mass, GSD dynamic range, and dust-to-gas ratio, result in lower dust opacities. This reduced opacity permits interstellar radiation field and internal radiation photons to penetrate more deeply. This leads to rapid gas heating and inhibited star formation. Star formation efficiency is highly sensitive to grain size, with an order-of-magnitude reduction when grain size dynamic range increases from 10⁻³─0.1 μm to 0.1─10 μm. Additionally, warmer gas suppresses low-mass star formation, and decreased opacities result in a greater proportion of gas in diffuse ionized structures.

Stars form in dense cores within molecular clouds, and newly formed stars influence their natal environments. How stellar feedback impacts core properties and evolution has been the subject of extensive investigation. We performed a hierarchical clustering (dendrogram) analysis of a STARFORGE (STAR FORmation in Gaseous Environments) simulation, modelling a giant molecular cloud to identify gas overdensities (cores) and study changes in their radius, mass, velocity dispersion, and virial parameter with respect to stellar feedback. We binned these cores on the basis of the fraction of gas affected by protostellar outflows, stellar winds, and supernovae and analysed the property distributions for each feedback bin. We find that cores that experience more feedback influence are smaller. Feedback notably enhances the velocity dispersion and virial parameter of the cores, more so than it reduces their radius. This is also evident in the linewidth─size relation, according to which cores in higher-feedback bins exhibit higher velocities than their similarly sized pristine counterparts. We conclude that stellar feedback mechanisms, which impart momentum to the molecular cloud, simultaneously compress and disperse the dense molecular gas.

Stars form within dense cores composed of both gas and dust within molecular clouds. However, despite the crucial role that dust plays in the star formation process, its dynamics is frequently overlooked, with the common assumption being a constant, spatially uniform dust-to-gas ratio and grain size spectrum. In this study, we introduce a set of radiation-dust-magnetohydrodynamic simulations of star-forming molecular clouds from the STARFORGE project. These simulations expand upon the earlier radiation MHD models, which included cooling, individual star formation, and feedback. Notably, they explicitly address the dynamics of dust grains, considering radiation, drag, and Lorentz forces acting on a diverse size spectrum of live dust grains. We find that once stars exceed a certain mass threshold (∼2 M _⊙), their emitted radiation can evacuate dust grains from their vicinity, giving rise to a dust-suppressed zone of size ∼100 au. This removal of dust, which interacts with gas through cooling, chemistry, drag, and radiative transfer, alters the gas properties in the region. Commencing during the early accretion stages and preceding the main-sequence phase, this process results in a mass-dependent depletion in the accreted dust-to-gas (ADG) mass ratio within both the circumstellar disk and the star. We predict that massive stars (≳10 M _⊙) would exhibit ADG ratios that are approximately 1 order of magnitude lower than that of their parent clouds. Consequently, stars, their disks, and circumstellar environments would display notable deviations in the abundances of elements commonly associated with dust grains, such as carbon and oxygen.

Cosmic rays (CRs) are the primary driver of ionization in star-forming molecular clouds (MCs). Despite their potential impacts on gas dynamics and chemistry, no simulations of star cluster formation following the creation of individual stars have included explicit cosmic-ray transport (CRT) to date. We conduct the first numerical simulations following the collapse of a 2000M _⊙ MC and the subsequent star formation including CRT using the STAR FORmation in Gaseous Environments framework implemented in the GIZMO code. We show that when CRT is streaming-dominated, the CR energy in the cloud is strongly attenuated due to energy losses from the streaming instability. Consequently, in a Milky Way─like environment the median CR ionization rate in the cloud is low (ζ ≲ 2 × 10⁻¹⁹ s⁻¹) during the main star-forming epoch of the calculation and the impact of CRs on the star formation in the cloud is limited. However, in high-CR environments, the CR distribution in the cloud is elevated (ζ ≲ 6 × 10⁻¹⁸), and the relatively higher CR pressure outside the cloud causes slightly earlier cloud collapse and increases the star formation efficiency by 50% to ∼13%. The initial mass function is similar in all cases except with possible variations in a high-CR environment. Further studies are needed to explain the range of ionization rates observed in MCs and explore star formation in extreme CR environments.

Recently, we demonstrated self-consistent formation of strongly-magnetized quasar accretion disks (QADs) from cosmological radiation-magnetohydrodynamic-thermochemical galaxy-star formation simulations, including the full STARFORGE physics shown previously to produce a reasonable IMF under typical ISM conditions. Here we study star formation and the stellar IMF in QADs, on scales from 100 au to 10 pc from the SMBH. We show it is critical to include physics often previously neglected, including magnetic fields, radiation, and (proto)stellar feedback. Closer to the SMBH, star formation is suppressed, but the (rare) stars that do form exhibit top-heavy IMFs. Stars can form only in special locations (e.g. magnetic field switches) in the outer QAD. Protostars accrete their natal cores rapidly but then dynamically decouple from the gas and 'wander,' ceasing accretion on timescales ~100 yr. Their jets control initial core accretion, but the ejecta are 'swept up' into the larger-scale QAD flow without much dynamical effect. The strong tidal environment strongly suppresses common-core multiplicity. The IMF shape depends sensitively on un-resolved dynamics of protostellar disks (PSDs), as the global dynamical times can become incredibly short (< yr) and tidal fields are incredibly strong, so whether PSDs can efficiently transport angular momentum or fragment catastrophically at <10 au scales requires novel PSD simulations to properly address. Most analytic IMF models and analogies with planet formation in PSDs fail qualitatively to explain the simulation IMFs, though we discuss a couple of viable models.

It has recently become possible to zoom-in from cosmological to sub-pc scales in galaxy simulations to follow accretion onto supermassive black holes (SMBHs). However, at some point the approximations used on ISM scales (e.g. optically-thin cooling and stellar-population-integrated star formation [SF] and feedback [FB]) break down. We therefore present the first cosmological radiation-magnetohydrodynamic (RMHD) simulation which self-consistently combines the FIRE physics (relevant on galactic/ISM scales where SF/FB are ensemble-averaged) and STARFORGE physics (relevant on small scales where we track individual (proto)stellar formation and evolution), together with explicit RMHD (including non-ideal MHD and multi-band M1-RHD) which self-consistently treats both optically-thick and thin regimes. This allows us to span scales from ~100 Mpc down to <100 au (~300 Schwarzschild radii) around a SMBH at a time where it accretes as a bright quasar, in a single simulation. We show that accretion rates up to $\sim 10-100\,{\rm M_{\odot}\,yr^{-1}}$ can be sustained into the accretion disk at $\ll 10^{3}\,R_{\rm schw}$, with gravitational torques between stars and gas dominating on sub-kpc scales until star formation is shut down on sub-pc scales by a combination of optical depth to cooling and strong magnetic fields. There is an intermediate-scale, flux-frozen disk which is gravitoturbulent and stabilized by magnetic pressure sustaining strong turbulence and inflow with persistent spiral modes. In this paper we focus on how gas gets into the small-scale disk, and how star formation is efficiently suppressed.

In a companion paper, we reported the self-consistent formation of quasar accretion disks with inflow rates $\sim 10\,{\rm M_{\odot}\,yr^{-1}}$ down to <300 Schwarzschild radii from cosmological radiation-magneto-thermochemical-hydrodynamical galaxy and star formation simulations. We see the formation of a well-defined, steady-state accretion disk which is stable against star formation at sub-pc scales. The disks are optically thick, with radiative cooling balancing accretion, but with properties that are distinct from those assumed in most previous accretion disk models. The pressure is strongly dominated by (primarily toroidal) magnetic fields, with a plasma $β\sim 10^{-4}$ even in the disk midplane. They are qualitatively distinct from magnetically elevated or arrested disks. The disks are strongly turbulent, with trans-Alfvenic and highly super-sonic turbulence, and balance this via a cooling time that is short compared to the disk dynamical time, and can sustain highly super-Eddington accretion rates. Their surface and 3D densities at $\sim 10^{3}-10^{5}$ gravitational radii are much lower than in a Shakura-Sunyaev disk, with important implications for their thermo-chemistry and stability. We show how the magnetic field strengths and geometries arise from rapid advection of flux with the inflow from much weaker galaxy-scale fields in these 'flux-frozen' disks, and how this stabilizes the disk and gives rise to efficient torques. Re-simulating without magnetic fields produces catastrophic fragmentation with a vastly smaller, lower-$\dot{M}$ Shakura-Sunyaev-like disk.

Context. Massive stars play important roles throughout the universe; however, their formation remains poorly understood. Observations of jets and outflows in high-mass star-forming regions, as well as surveys of young stellar object (YSO) content, can help test theoretical models of massive star formation. Aims: We aim at characterizing the massive star-forming region AFGL 5180 in the near-infrared (NIR), identifying outflows and relating these to sub-mm/mm sources, as well as surveying the overall YSO surface number density to compare to massive star formation models. Methods: Broad- and narrow-band imaging of AFGL 5180 was made in the NIR with the Large Binocular Telescope, in both seeing-limited (~0.5″) and high angular resolution (~0.09″) Adaptive Optics (AO) modes, as well as with the Hubble Space Telescope. Archival continuum data from the Atacama Millimeter/Submillimeter Array (ALMA) was also utilized. Results: At least 40 jet knots were identified via NIR emission from H₂ and [FeII] tracing shocked gas. Bright jet knots outflowing from the central most massive protostar, S4 (estimated mass ~11 M_☉, via SED fitting), are detected towards the east of the source and are resolved in fine detail with the AO imaging. Additional knots are distributed throughout the field, likely indicating the presence of multiple driving sources. Sub-millimeter sources detected by ALMA are shown to be grouped in two main complexes, AFGL 5180 M and a small cluster ~15″ (0.15 pc in projection) to the south, AFGL 5180 S. From our NIR continuum images we identify YSO candidates down to masses of ~0.1 M_☉. Combined with the sub-mm sources, this yields a surface number density of such YSOs of N_* ~ 10³pc^‒2 within a projected radius of about 0.1 pc. Such a value is similar to those predicted by models of both core accretion from a turbulent clump environment and competitive accretion. The radial profile of N_* is relatively flat on scales out to 0.2 pc, with only modest enhancement around the massive protostar inside 0.05 pc, which provides additional constraints on these massive star formation models. Conclusions: This study demonstrates the utility of high-resolution NIR imaging, in particular with AO, for detecting outflow activity and YSOs in distant regions. The presented images reveal the complex morphology of outflow-shocked gas within the large-scale bipolar flow of a massive protostar, as well as clear evidence for several other outflow driving sources in the region. Finally, this work presents a novel approach to compare the observed YSO surface number density from our study against different models of massive star formation. The reduced images are available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr (ftp://130.79.128.5) or via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/682/A2

Most stars form in highly clustered environments within molecular clouds, but eventually disperse into the distributed stellar field population. Exactly how the stellar distribution evolves from the embedded stage into gas-free associations and (bound) clusters is poorly understood. We investigate the long-term evolution of stars formed in the STARFORGE simulation suite - a set of radiation-magnetohydrodynamic simulations of star-forming turbulent clouds that include all key stellar feedback processes inherent to star formation. We use NBODY6++GPU to follow the evolution of the young stellar systems after gas removal. We use HDBSCAN to define stellar groups and analyse the stellar kinematics to identify the true bound star clusters. The conditions modeled by the simulations, i.e. global cloud surface densities below 0.15 g cm^-2, star formation efficiencies below 15 per cent, and gas expulsion time-scales shorter than a free fall time, primarily produce expanding stellar associations and small clusters. The largest star clusters, which have ~1000 bound members, form in the densest and lowest velocity dispersion clouds, representing ~32 and 39 per cent of the stars in the simulations, respectively. The cloud's early dynamical state plays a significant role in setting the classical star formation efficiency versus bound fraction relation. All stellar groups follow a narrow mass-velocity dispersion power-law relation at 10 Myr with a power-law index of 0.21. This correlation result in a distinct mass-size relationship for bound clusters. We also provide valuable constraints on the gas dispersal time-scale during the star formation process and analyse the implications for the formation of bound systems.

When a detailed model of a stellar population is unavailable, it is most common to assume that stellar masses are independently and identically distributed according to some distribution: the universal initial mass function (IMF). However, stellar masses resulting from causal, long-ranged physics cannot be truly random and independent, and the IMF may vary with environment. To compare stochastic sampling with a physical model, we run a suite of 100 STARFORGE radiation magnetohydrodynamics simulations of low-mass star cluster formation in $2000M_\odot$ clouds that form $\sim 200$ stars each on average. The stacked IMF from the simulated clouds has a sharp truncation at $\sim 28 M_\odot$, well below the typically-assumed maximum stellar mass $M_{\rm up} \sim 100-150M_\odot$ and the total cluster mass. The sequence of star formation is not totally random: massive stars tend to start accreting sooner and finish later than the average star. However, final cluster properties such as maximum stellar mass and total luminosity have a similar amount of cloud-to-cloud scatter to random sampling. Therefore stochastic sampling does not generally model the stellar demographics of a star cluster as it is forming, but may describe the end result fairly well, if the correct IMF -- and its environment-dependent upper cutoff -- are known.

We explore the relation between stellar surface density and gas surface density (the star-gas, or S-G, correlation) in a 20,000 M _☉ simulation from the STAR FORmation in Gaseous Environments (STARFORGE) project. We create synthetic observations based on the Spitzer and Herschel telescopes by modeling contamination by active galactic nuclei, smoothing based on angular resolution, cropping the field of view, and removing close neighbors and low-mass sources. We extract S-G properties such as the dense gas-mass fraction, the Class II:I ratio, and the S-G correlation (Σ_YSO/Σ_gas) from the simulation and compare them to observations of giant molecular clouds, young clusters, and star-forming regions, as well as to analytical models. We find that the simulation reproduces trends in the counts of young stellar objects and the median slope of the S-G correlation. This implies that the S-G correlation is not simply the result of observational biases, but is in fact a real effect. However, other statistics, such as the Class II:I ratio and dense gas-mass fraction, do not always match observed equivalents in nearby clouds. This motivates further observations covering the full simulation age range and more realistic modeling of cloud formation.

Accurately quantifying the impact of radiation feedback in star formation is challenging. To address this complex problem, we employ deep-learning techniques known as denoising diffusion probabilistic models (DDPMs) to predict the interstellar radiation field (ISRF) strength based on three-band dust emission at 4.5, 24, and 250 μm. We adopt magnetohydrodynamic simulations from the STARFORGE project that model star formation and giant molecular cloud (GMC) evolution. We generate synthetic dust emission maps matching observed spectral energy distributions in the Monoceros R2 (MonR2) GMC. We train DDPMs to estimate the ISRF using synthetic three-band dust emission. The dispersion between the predictions and true values is within a factor of 0.1 for the test set. We extended our assessment of the diffusion model to include new simulations with varying physical parameters. While there is a consistent offset observed in these out-of-distribution simulations, the model effectively constrains the relative intensity to within a factor of 2. Meanwhile, our analysis reveals a weak correlation between the ISRF solely derived from dust temperature and the actual ISRF. We apply our trained model to predict the ISRF in MonR2, revealing a correspondence between intense ISRF, bright sources, and high dust emission, confirming the model's ability to capture ISRF variations. Our model robustly predicts radiation feedback distribution, even in complex, poorly constrained ISRF environments like those influenced by nearby star clusters. However, precise ISRF predictions require an accurate training data set mirroring the target molecular cloud's unique physical conditions.

Most observed stars are part of a multiple star system, but the formation of such systems and the role of environment and various physical processes is still poorly understood. We present a suite of radiation-magnetohydrodynamic simulations of star-forming molecular clouds from the STARFORGE project that include stellar feedback with varied initial surface density, magnetic fields, level of turbulence, metallicity, interstellar radiation field, simulation geometry and turbulent driving. In our fiducial cloud, the raw simulation data reproduces the observed multiplicity fractions for Solar-type and higher mass stars, similar to previous works. However, after correcting for observational incompleteness the simulation underpredicts these values. The discrepancy is likely due to the lack of disc fragmentation, as the simulation only resolves multiples that form either through capture or core fragmentation. The raw mass distribution of companions is consistent with randomly drawing from the initial mass function for the companions of $\gt 1\, \mathrm{M}_{\rm \odot }$ stars. However, accounting for observational incompleteness produces a flatter distribution similar to observations. We show that stellar multiplicity changes as the cloud evolves and anticorrelates with stellar density. This relationship also explains most multiplicity variations between runs, i.e. variations in the initial conditions that increase stellar density (increased surface density, reduced turbulence) also act to decrease multiplicity. While other parameters, such as metallicity, interstellar radiation, and geometry significantly affect the star formation history or the IMF, varying them produces no clear trend in stellar multiplicity properties.

One of the key mysteries of star formation is the origin of the stellar initial mass function (IMF). The IMF is observed to be nearly universal in the Milky Way and its satellites, and significant variations are only inferred in extreme environments, such as the cores of massive elliptical galaxies and the Central Molecular Zone. In this work, we present simulations from the STARFORGE project that are the first cloud-scale radiation-magnetohydrodynamic simulations that follow individual stars and include all relevant physical processes. The simulations include detailed gas thermodynamics, as well as stellar feedback in the form of protostellar jets, stellar radiation, winds, and supernovae. In this work, we focus on how stellar radiation, winds, and supernovae impact star-forming clouds. Radiative feedback plays a major role in quenching star formation and disrupting the cloud; however, the IMF peak is predominantly set by protostellar jet physics. We find that the effect of stellar winds is minor, and supernovae 'occur too late' to affect the IMF or quench star formation. We also investigate the effects of initial conditions on the IMF. We find that the IMF is insensitive to the initial turbulence, cloud mass, and cloud surface density, even though these parameters significantly shape the star formation history of the cloud, including the final star formation efficiency. Meanwhile, the characteristic stellar mass depends weakly on metallicity and the interstellar radiation field, which essentially set the average gas temperature. Finally, while turbulent driving and the level of magnetization strongly influence the star formation history, they only influence the high-mass slope of the IMF.

Stars form in dense, clustered environments, where feedback from newly formed stars eventually ejects the gas, terminating star formation and leaving behind one or more star clusters. Using the STARFORGE simulations, it is possible to simulate this process in its entirety within a molecular cloud, while explicitly evolving the gas radiation and magnetic fields and following the formation of individual, low-mass stars. We find that individual star-formation sites merge to form ever larger structures, while still accreting gas. Thus clusters are assembled through a series of mergers. During the cluster assembly process, a small fraction of stars are ejected from their clusters; we find no significant difference between the mass distribution of the ejected stellar population and that of stars inside clusters. The star-formation sites that are the building blocks of clusters start out mass segregated with one or a few massive stars at their centre. As they merge the newly formed clusters maintain this feature, causing them to have mass-segregated substructures without themselves being centrally condensed. The merged clusters relax to a centrally condensed mass-segregated configuration through dynamical interactions between their members, but this process does not finish before feedback expels the remaining gas from the cluster. In the simulated runs, the gas-free clusters then become unbound and breakup. We find that turbulent driving and a periodic cloud geometry can significantly reduce clustering and prevent gas expulsion. Meanwhile, the initial surface density and level of turbulence have little qualitative effect on cluster evolution, despite the significantly different star formation histories.

We analyse the first giant molecular cloud (GMC) simulation to follow the formation of individual stars and their feedback from jets, radiation, winds, and supernovae, using the STARFORGE framework in the GIZMO code. We evolve the GMC for $\sim 9 \rm Myr$, from initial turbulent collapse to dispersal by feedback. Protostellar jets dominate feedback momentum initially, but radiation and winds cause cloud disruption at $\sim 8{{\ \rm per\ cent}}$ star formation efficiency (SFE), and the first supernova at $8.3\, \rm Myr$ comes too late to influence star formation significantly. The per-free-fall SFE is dynamic, accelerating from 0 per cent to $\sim 18{{\ \rm per\ cent}}$ before dropping quickly to <1 per cent, but the estimate from YSO counts compresses it to a narrower range. The primary cluster forms hierarchically and condenses to a brief ($\sim 1\, \mathrm{Myr}$) compact ($\sim 1\, \rm pc$) phase, but does not virialize before the cloud disperses, and the stars end as an unbound expanding association. The initial mass function resembles the Chabrier (2005) form with a high-mass slope α = -2 and a maximum mass of 55 M_☉. Stellar accretion takes $\sim 400\, \rm kyr$ on average, but $\gtrsim 1\,\rm Myr$ for >10 M_☉ stars, so massive stars finish growing latest. The fraction of stars in multiples increase as a function of primary mass, as observed. Overall, the simulation much more closely resembles reality, compared to previous versions that neglected different feedback physics entirely. But more detailed comparison with synthetic observations will be needed to constrain the theoretical uncertainties.

Simulations of isolated giant molecular clouds (GMCs) are an important tool for studying the dynamics of star formation, but their turbulent initial conditions (ICs) are uncertain. Most simulations have either initialized a velocity field with a prescribed power spectrum on a smooth density field (failing to model the full structure of turbulence) or 'stirred' turbulence with periodic boundary conditions (which may not model real GMC boundary conditions). We develop and test a new GMC simulation setup (called TURBSPHERE) that combines advantages of both approaches: we continuously stir an isolated cloud to model the energy cascade from larger scales, and use a static potential to confine the gas. The resulting cloud and surrounding envelope achieve a quasi-equilibrium state with the desired hallmarks of supersonic ISM turbulence (e.g. density PDF and a ~k^-2 velocity power spectrum), whose bulk properties can be tuned as desired. We use the final stirred state as initial conditions for star formation simulations with self-gravity, both with and without continued driving and protostellar jet feedback, respectively. We then disentangle the respective effects of the turbulent cascade, simulation geometry, external driving, and gravity/MHD boundary conditions on the resulting star formation. Without external driving, the new setup obtains results similar to previous simple spherical cloud setups, but external driving can suppress star formation considerably in the new setup. Periodic box simulations with the same dimensions and turbulence parameters form stars significantly slower, highlighting the importance of boundary conditions and the presence or absence of a global collapse mode in the results of star formation calculations.

Star formation involves gravity, turbulence, magnetic fields, and feedback from new stars through jets, radiation and winds. The evolution of the density probability distribution function (ρ-PDF) is directly related to the star formation rate (SFR), forming the basis of several star formation models. We utilize two runs from the STARFORGE simulation suite that follow the evolution of molecular clouds, while resolving individual stars and including all gas and feedback physics. The two runs have different initial conditions, one is a periodic box with driven turbulence (Box), while the other is an isolated cloud without turbulent driving (Sphere). We find that the ρ-PDF for both runs is initially well-fit by a log-normal (LN) plus a power-law (PL) function. However, as the SFR peaks, the PDF for the Sphere run becomes well-fit by just a wide LN. Conversely, the Box run PDF remains well-fit by a LN+PL function for the entirety of the run.

Many astrophysical hydrodynamics simulations must account for gravity, and evaluating the gravitational field at the positions of all resolution elements can incur significant cost. Typical algorithms update the gravitational field at the position of each resolution element every time the element is updated hydrodynamically, but the actual required update frequencies for hydrodynamics and gravity can be different in general. We show that the gravity calculation in hydrodynamics simulations can be optimized by only updating gravity on a time-scale dictated by the already determined maximum time-step for accurate gravity integration ∆t_grav, while staying well within the typical error budget of hydro schemes and gravity solvers. Our implementation in the GIZMO code uses the time-scale derived from the tidal tensor $t_{\rm tidal} = \Vert \mathbf {T}\Vert ^{-1/2}$ to determine ∆t_grav and the force update frequency in turn, and uses the rate of change of acceleration evaluated by the gravity solver to construct a predictor of the acceleration for use between updates. We test the scheme on standard self-gravitating hydrodynamics test problems, finding solutions very close to the standard scheme while evaluating far fewer gravity forces, optimizing the simulations. We also demonstrate a $\sim 70{{\ \rm per\ cent}}$ speed-up in an example simulation of a giant molecular cloud. In general, this scheme introduces a new tunable parameter for obtaining an optimal compromise between accuracy and computational cost, in conjunction with, e.g. time-step tolerance, numerical resolution, and gravity solver tolerance.

We present STARFORGE (STAR FORmation in Gaseous Environments): a new numerical framework for 3D radiation magnetohydrodynamic (MHD) simulations of star formation that simultaneously follow the formation, accretion, evolution, and dynamics of individual stars in massive giant molecular clouds (GMCs), while accounting for stellar feedback, including jets, radiative heating and momentum, stellar winds, and supernovae. We use the GIZMO code with the MFM mesh-free Lagrangian MHD method, augmented with new algorithms for gravity, time-stepping, sink particle formation and accretion, stellar dynamics, and feedback coupling. We survey a wide range of numerical parameters/prescriptions for sink formation and accretion and find very small variations in star formation history and the IMF (except for intentionally unphysical variations). Modules for mass-injecting feedback (winds, SNe, and jets) inject new gas elements on the fly, eliminating the lack of resolution in diffuse feedback cavities otherwise inherent in Lagrangian methods. The treatment of radiation uses GIZMO's radiative transfer solver to track five frequency bands (IR, optical, NUV, FUV, ionizing), coupling direct stellar emission and dust emission with gas heating and radiation pressure terms. We demonstrate accurate solutions for SNe, winds, and radiation in problems with known similarity solutions, and show that our jet module is robust to resolution and numerical details, and agrees well with previous AMR simulations. STARFORGE can scale up to massive (>10⁵ M_☉) GMCs on current supercomputers while predicting the stellar (≳0.1 M_☉) range of the IMF, permitting simulations of both high- and low-mass cluster formation in a wide range of conditions.

The initial mass function (IMF) of stars is a key quantity affecting almost every field of astrophysics, yet it remains unclear what physical mechanisms determine it. We present the first runs of the STAR FORmation in Gaseous Environments project, using a new numerical framework to follow the formation of individual stars in giant molecular clouds (GMCs) using the GIZMO code. Our suite includes runs with increasingly complex physics, starting with isothermal ideal magnetohydrodynamics (MHD) and then adding non-isothermal thermodynamics and protostellar outflows. We show that without protostellar outflows the resulting stellar masses are an order of magnitude too high, similar to the result in the base isothermal MHD run. Outflows disrupt the accretion flow around the protostar, allowing gas to fragment and additional stars to form, thereby lowering the mean stellar mass to a value similar to that observed. The effect of jets upon global cloud evolution is most pronounced for lower mass GMCs and dense clumps, so while jets can disrupt low-mass clouds, they are unable to regulate star formation in massive GMCs, as they would turn an order unity fraction of the mass into stars before unbinding the cloud. Jets are also unable to stop the runaway accretion of massive stars, which could ultimately lead to the formation of stars with masses ${\gt}500\, \mathrm{M}_{\rm \odot }$. Although we find that the mass scale set by jets is insensitive to most cloud parameters (i.e. surface density, virial parameter), it is strongly dependent on the momentum loading of the jets (which is poorly constrained by observations) as well as the temperature of the parent cloud, which predicts slightly larger IMF variations than observed. We conclude that protostellar jets play a vital role in setting the mass scale of stars, but additional physics are necessary to reproduce the observed IMF.

Understanding the evolution of self-gravitating, isothermal, magnetized gas is crucial for star formation, as these physical processes have been postulated to set the initial mass function (IMF). We present a suite of isothermal magnetohydrodynamic (MHD) simulations using the GIZMO code that follow the formation of individual stars in giant molecular clouds (GMCs), spanning a range of Mach numbers found in observed GMCs ($\mathcal {M} \sim 10\!-\!50$). As in past works, the mean and median stellar masses are sensitive to numerical resolution, because they are sensitive to low-mass stars that contribute a vanishing fraction of the overall stellar mass. The mass-weighted median stellar mass M₅₀ becomes insensitive to resolution once turbulent fragmentation is well resolved. Without imposing Larson-like scaling laws, our simulations find $M_\mathrm{50} \,\, \buildrel\propto \over \sim \,\,M_\mathrm{0} \mathcal {M}^{-3} \alpha _\mathrm{turb}\, \mathrm{SFE}^{1/3}$ for GMC mass M₀, sonic Mach number $\mathcal {M}$, virial parameter α_turb, and star formation efficiency SFE = M_⋆/M₀. This fit agrees well with previous IMF results from the RAMSES, ORION2, and SPHNG codes. Although M₅₀ has no significant dependence on the magnetic field strength at the cloud scale, MHD is necessary to prevent a fragmentation cascade that results in non-convergent stellar masses. For initial conditions and SFE similar to star-forming GMCs in our Galaxy, we predict M₅₀ to be $\gt 20 \, \mathrm{M}_{\odot }$, an order of magnitude larger than observed ($\sim 2 \, \mathrm{M}_\odot$), together with an excess of brown dwarfs. Moreover, M₅₀ is sensitive to initial cloud properties and evolves strongly in time within a given cloud, predicting much larger IMF variations than are observationally allowed. We conclude that physics beyond MHD turbulence and gravity are necessary ingredients for the IMF.

We describe a new adaptive time-step criterion for integrating gravitational motion, which uses the tidal tensor to estimate the local dynamical time-scale and scales the time-step proportionally. This provides a better candidate for a truly general-purpose gravitational time-step criterion than the usual prescription derived from the gravitational acceleration, which does not respect the equivalence principle, breaks down when $\boldsymbol {a}=0$ , and does not obey the same dimensional scaling as the true time-scale of orbital motion. We implement the tidal time-step criterion in the simulation code GIZMO, and examine controlled tests of collisionless galaxy and star cluster models, as well as galaxy merger simulations. The tidal criterion estimates the dynamical time faithfully, and generally provides a more efficient time-stepping scheme compared to an acceleration criterion. Specifically, the tidal criterion achieves order-of-magnitude smaller energy errors for the same number of force evaluations in potentials with inner profiles shallower than ρ ∝ r^-1 (I.e. where $\boldsymbol {a}\rightarrow 0$ ), such as star clusters and cored galaxies. For a given problem these advantages must be weighed against the additional overhead of computing the tidal tensor on-the-fly, but in many cases this overhead is small.

The evolution of self-gravitating clouds of isothermal gas forms the basis of many star formation theories. Therefore it is important to know under what conditions such a cloud will undergo monolithic collapse into a single, massive object, or will fragment into a spectrum of smaller ones. And if it fragments, do initial conditions (e.g. Jeans mass, sonic mass) influence the mass function of the fragments, as predicted by many theories of star formation? In this paper we show that the relevant parameter separating monolithic collapse from fragmentation is not the Mach number of the initial turbulence (as suspected by many), but the infall Mach number M_infall∼ √{G M/(R c_s^2)}, equivalent to the number of Jeans masses in the initial cloud N_J. We also show that fragmenting clouds produce a power-law mass function with slopes close to the expected -2 (i.e. equal mass in all logarithmic mass intervals). However, the low-mass cut-off of this mass function is entirely numerical; the initial properties of the cloud have no effect on it. In other words, if M_infall≫ 1, fragmentation proceeds without limit to masses much smaller than the initial Jeans mass.