© The Authors 2026

1 Introduction

At least half of all stars, including our Sun and almost all massive stars, formed in clustered systems of hundreds to thousands of objects (Lada & Lada 2003; Zinnecker & Yorke 2007; Bressert et al. 2010; Adams 2010; Megeath et al. 2016, 2022; Adamo et al. 2020). The prevalence of clustered star formation is a direct consequence of the hierarchical nature of the interstellar medium (ISM) from which stars form (Elmegreen 2010). The ISM shows substructures on all spatial scales (Scalo 1985; Elmegreen & Falgarone 1996; Williams et al. 2000; Molinari et al. 2010): molecular clouds (≥ 1 pc), filaments (∼0.1−1 pc), clumps (∼0.3− 1 pc), and cores (≲ 0.1 pc). The fragmentation of these structures leads to the observed distribution of young stars, which are typically arranged in groups or clumpy subclusters (Kuhn et al. 2014; Zhou et al. 2024b). Only a fraction of these young systems ultimately remain as stellar bounded clusters after they merge into smoother distributions (Kruijssen 2012; Grudić et al. 2018; Li et al. 2019; Krumholz & McKee 2020; Krause et al. 2020; Zhou et al. 2025). Star clusters form in dense clumps that are typically located at the intersection of multiple filaments (Myers 2009; Schneider et al. 2012; Motte et al. 2018; Zhou et al. 2022). These dense clusters contain one or more high-mass stars (M ≥ 8 M_⊙), whose stellar winds, ionization, heating, radiation pressure, and ultimately supernova explosions play a fundamental role in shaping the evolution of the interstellar medium, from planetary systems to the entire galactic ecosystem (Portegies Zwart et al. 2010).

The formation of dense clusters is a complex theoretical problem. It is governed by the interplay of gravitational collapse, turbulence, and magnetic fields in an environment that is continuously altered by radiative and mechanical feedback from the embedded high-mass stars still in formation (Mac Low & Klessen 2004; Dale & Bonnell 2011; Krumholz et al. 2014; Krumholz & McKee 2020; Li et al. 2014; Girichidis et al. 2020; Padoan et al. 2020; Grudić et al. 2021; Guszejnov et al. 2021; Vázquez-Semadeni et al. 2019, 2025). In addition to these processes, the final structure and properties of clusters are influenced by other factors. These include the initial cloud conditions and morphology (Girichidis et al. 2011), which may also be related to the local galactic environment (Smith et al. 2020); the possible coalescence of protostellar cores in the early evolutionary phases (Dib 2023); and the hierarchical merging where distinct groups of pre-existing star-forming cores or young subclusters assemble (Maschberger et al. 2010).

Current models of high-mass stars and cluster formation are divided into categories that mostly differ by the dynamical state of the parent structure that fragments and where the material that ends up in the protostars comes from. One category is composed of static models that assume that high-mass star formation starts in substructures created by the supersonic turbulence present in molecular clouds (McKee & Tan 2003), with possible contributions from the magnetic field. The material is initially gathered in these substructures and builds up the entire mass reservoir for the star formation process. In this category, there are different scenarios distinguished by the structure type that is initially assembled: cores for the core-fed models, for example the turbulent model of McKee & Tan (2003), or the entire clump for the in situ formation scenario introduced by Longmore et al. (2014). In contrast, another category posits a more dynamical scenario where star formation is coeval to a further accumulation of gas in the high-density structures: cores and clumps. This additional contribution increases the reservoir of material available for accretion onto young protostars. Longmore et al. (2014) refer to this category as conveyor belt mode to recall that the mass of star-forming material is never assembled at a single time (see also the discussion in Krumholz & McKee 2020), but is continuously transferred. The competitive accretion model introduced by Bonnell et al. (2001) falls into this category since it prescribes that cores gather their masses from the clump local environment, with a mass growth that is more efficient for objects close to the bottom of the gravitational well (Bonnell & Bate 2006). While many models agree that flows funnel material from low to high densities to assemble massive stars, they differ on what drives the local accretion flow. For example, it is driven by gravity in the global hierarchical collapse (GHC) model (Ballesteros-Paredes et al. 2018; Vázquez-Semadeni et al. 2019, 2025)), while by the large-scale random (turbulent) velocity field in the inertial flow scenario of Padoan et al. (2020).

Over the years, several studies have been dedicated to characterizing the initial conditions of star formation by observing young stellar clusters (Gutermuth et al. 2009; Kuhn et al. 2014; Parker & Alves de Oliveira 2017; Parker 2018; Dib et al. 2018). Although they are considerably young, these systems are already influenced by their dynamical evolution, which does not allow one to trace back where and under which initial conditions the stars originally formed. This leaves severe uncertainties about how the time-dependent process of fragmentation started and how it proceeds during the early evolution of the system, before most of the gas is dispersed, questions that can be answered only by observing the distribution of star-forming cores.

Modern millimeter interferometers, such as the Atacama Large Millimeter/submillimeter Array (ALMA), the NOrthern Extended Millimeter Array (NOEMA), and the SubMillimeter Array (SMA) have achieved the angular resolution and sensitivity required to unveil the fragmentation properties in the dense clumps, which are the crowded and heavily extincted sites where high-mass stars and clusters form, typically located a few kiloparsec from the Sun (Zhang et al. 2014, 2015; Csengeri et al. 2017; Sanhueza et al. 2017; Beuther et al. 2018; Fontani et al. 2018; Svoboda et al. 2019; Sadaghiani et al. 2020; Moscadelli et al. 2021; Palau et al. 2021; Liu et al. 2020, 2022; Morii et al. 2021, 2023; Traficante et al. 2023; Avison et al. 2023; Ishihara et al. 2024; Liu et al. 2024; Zhou et al. 2024a). In recent years, an important observational effort has been carried out to study these massive clumps, identifying and characterizing the embedded population of star-forming fragments in order to provide important observational constraints to the models of high-mass stars and cluster formation. Several surveys focus on clumps in specific evolutionary phases, either very young and quiescent systems embedded in infrared dark clouds (Zhang et al. 2014; Sanhueza et al. 2017; Svoboda et al. 2019; Sanhueza et al. 2019; Anderson et al. 2021; Morii et al. 2023; Rigby et al. 2024) or systems in more advanced stages that are associated with infrared bright emission (Beuther et al. 2018; Liu et al. 2020, 2024; Xu et al. 2024). Other surveys selected their targets to cover a wide range of evolutionary phases (Avison et al. 2023; Traficante et al. 2023; Ishihara et al. 2024), but observed only a limited number of systems. The ALMA-IMF Large Program (Motte et al. 2022) surveyed 15 massive Galactic protoclusters near the Sun (d<5.5 kpc) for a detailed analysis of the core mass function in different environments and evolutionary stages. These observations have built up a rich, but rather heterogeneous dataset of observed properties of the fragments hosted in about two hundred massive clumps.

Although observed core masses are an important reference point for star formation models, the spatial distribution of cores offers an additional and robust diagnostic with critical information on how clusters assemble and the relevance of various physical mechanisms in their formation process. Fragmentation leaves a different signature imprinted on the separations of cores in a turbulent or gravity-dominated medium, which is assumed in different categories of models, characterized by different values for the Jeans length (Mac Low & Klessen 2004; Vázquez-Semadeni et al. 2025). Several observations have reported separations that are close to or smaller than the average thermal Jeans length of the hosting structure (clump or cloud; Palau et al. 2015; Kainulainen et al. 2017; Beuther et al. 2018; Sanhueza et al. 2019; Morii et al. 2024; Ishihara et al. 2024; Beuther et al. 2025). However, turbulent fragmentation is not completely ruled out as an additional contribution from turbulent support is required to explain observations in other systems (Rebolledo et al. 2020; Avison et al. 2023) and in nearby clouds (Ishihara et al. 2025). In general, distinct fragmentation patterns with strong core alignments, clustered arrangments with cores randomly distributed, or low fragmentation and rather symmetrical spatial distributions are observed, whether the dominant energetic contribution in the system is the magnetic field, turbulence, or gravity, respectively (Tang et al. 2019). Anisotropic accretion flows, whose presence may be traced by small-scale subfilaments, provide additional fragmentation, increasing the number of cores and their location in the system. In addition, dynamical effects, such as the merging of subclusters (Bonnell et al. 2003), clump global gravitational collapse (Vázquez-Semadeni et al. 2017), and tidal interactions between the cores (Zhou et al. 2024a; Li 2024), modify the overall structure and spatial distribution of young systems. Consequently, variations in these structural properties across different evolutionary stages can be used to determine the presence and impact of such dynamics. Observations suggest that the typical separation between cores is shorter in more evolved systems (Beuther et al. 2021; Traficante et al. 2023; Ishihara et al. 2024). Moreover, the core distribution is also found to be more compact in the later stages (Xu et al. 2024) and in regions with higher star formation rates (Dib & Henning 2019). These two pieces of evidence both support the presence of a strong gravitational pull, a signature that a global collapse, as predicted by the GHC theory (Vázquez-Semadeni et al. 2025), is active in massive clumps.

Additional information is provided by looking at the location of the massive members of the clustered system. Although young stellar clusters often show the presence of mass segregation (MS), with high-mass objects that are most likely to lie close to the cluster center (Allison et al. 2009; Parker & Alves de Oliveira 2017; Parker 2018; Dib & Henning 2019; Jia et al. 2025), it is unclear whether those patterns are connected to their initial formation or are due to their subsequent dynamical evolution (see also Portegies Zwart et al. 2010). A similar statistical analysis carried out on the dense core population shows that MS may be present since the early stages (Sanhueza et al. 2019; Sadaghiani et al. 2020; Nony et al. 2021; Morii et al. 2024; Xu et al. 2024), which is indicative that competitive accretion may be a dominant path of their formation (Bonnell & Bate 2006).

Although several indications were drawn from the available dataset, the diverse selection criteria of the targets, the different spatial resolutions and sensitivity levels achieved by the various surveys, and the low number of observed clumps limit the statistical analysis and direct comparison of the observed fragment population. The ALMA Cycle 7 large program ALMA evolutionary study of high-mass protocluster formation in the GALaxy (ALMAGAL) (2019.1.00195L, PI: S. Molinari, P. Schilke, C. Battersby, and P. Ho) was designed to mitigate these problems by observing an unprecedented large sample of 1013 clumps covering a wide range of physical conditions and Galactic locations (Molinari et al. 2025). The ALMAGAL program provides an extended and homogeneous dataset to study the fragmentation process in dense and massive clumps. In this paper we analyzed the spatial distribution of the cores identified in the ALMAGAL fields by Coletta et al. (2025) to provide a robust characterization of where the cores are located in massive and potentially cluster-forming clumps, what their typical separation is, how they are distributed, and how these properties vary across the early evolutionary phases. This work is a complement to the comparison between the physical properties of the ALMAGAL core population, specifically the level of fragmentation, core masses, and efficiency, and the large-scale clump parameters presented by Elia et al. (2026).

The structure of this paper is as follows. We discuss the details of the ALMAGAL program in Sect. 2, summarizing the observational strategy for the survey and the properties of the ALMAGAL clumps and the extracted core population used in the analysis. We present the results of the analysis of the spatial distribution of ALMAGAL cores in Sect. 3. Specifically, in Sect. 3.1 we describe the general morphological properties of the groups of cores, defining where they are distributed. In Sect. 3.2, we determine the separations between the cores and compare them with the clump averaged properties in Sect. 3.3, focusing on how the separations vary during different evolutionary phases. We discuss how they relate to the expectation of the Jeans fragmentation theory in Sect. 3.4. In Sect. 3.5, we present additional statistical descriptors to illustrate how cores are distributed and whether mass segregation is present in systems, checking if they show a variation across the different evolutionary phases. Finally, we discuss our results within the framework of the generally adopted models in Sect. 4.

2 Data

2.1 The ALMAGAL clump sample

The ALMAGAL large program (2019.1.00195.L, PI: Molinari, Schilke, Battersby, Ho) observed with ALMA in Band 6 (∼220 GHz, ∼1.4 mm) 1013 dense massive clumps located close to the Galactic Plane, and selected from the Hi-GAL (Molinari et al. 2010, 2016b; Elia et al. 2017, 2021) and RMS surveys (Lumsden et al. 2013), to investigate their fragmentation properties. The criteria adopted for the target selection are described in Molinari et al. (2025) and are based on thresholds on the Herschel continuum emission, the clump mass, M, and the average surface density, Σ_cl, determined from the Herschel Hi-GAL data (Molinari et al. 2016b; Elia et al. 2017, 2021). An initial analysis of the ALMA ACA data provided new estimates for the clump heliocentric distances, derived from the radial velocities measured from the spectral lines of multiple molecular tracers (Benedettini et al., in prep.). The ALMAGAL clumps, with physical properties updated thanks to these new measurements, have masses M between ∼100 and 12 000 M_⊙, bolometric luminosity, L, between ∼50 and ∼5 × 10⁵ L_⊙, and are located in the Galaxy at heliocentric distances, d_cl, from ∼1 up to 14 kpc, with only one clump lying closer than 1 kpc (Molinari et al. 2025; Elia et al. 2026). Notwithstanding these updates, the observed ALMAGAL clumps effectively sample the entire evolutionary path of massive Galactic clumps that are likely to host the formation of high-mass stars and clusters. In fact, they have clump-averaged surface densities, Σ_cl, ranging from ∼0.1 up to ∼10 g cm⁻², which covers the thresholds suggested by theoretical modeling and observational constraints for the formation of high-mass stars (Krumholz & McKee 2008; Kauffmann & Pillai 2010; Baldeschi et al. 2017, see also the discussion in Beuther et al. 2025). Moreover, we note that the ALMAGAL sample includes clumps with a ratio of their bolometric luminosity to their mass, L/M, which is uniformly distributed between ∼0.5 and ∼50 L_⊙/M_⊙, and includes cases with a L/M down to ∼0.05 and up to ∼400 L_⊙/M_⊙, although highly evolved phases are less sampled, as there are fewer systems with L/M>100 L_⊙/M_⊙ compared to the other evolutionary stages. Since the L/M ratio is a well-established indicator of the evolutionary stage of clumps (Molinari et al. 2016a), the ALMAGAL clumps provide a complete and rather homogeneous sampling of the evolution of massive clumps with active star formation. We present in Fig. 1 a thorough comparison between the L/M ratio, the surface density, Σ_cl, and the heliocentric distance, d_cl, of the ALMAGAL targets. There are 25 clumps located at distances greater than 9 kpc after updating their distances that are not shown in the left and central panels for clarity reasons, but the distribution of their properties is similar to the ones presented. Figure 1 shows that there is no relationship between the ratio L/M and either Σ_cl or d_cl, indicating that target selection did not introduce any strong bias due to an eventual correlation between these two quantities. In contrast, the ALMAGAL sample presents a relation between Σ_cl and d_cl, introduced as a consequence of the selection criteria on the clump mass.

Fig. 1 Distribution of the L/M ratio and average surface density, Σcl, as function of the clump heliocentric distance for the ALMAGAL sample, and comparison of these two properties. Each panel includes the 2D density distribution computed dividing the intervals of L/M, Σcl, and dcl into bins of width 0.3 dex, 0.2 dex, and 500 pc, respectively. The solid light gray lines indicate the contour levels of the 2D density distribution derived with kernel density estimation (KDE) corresponding to 10, 30, 50, 70, 95% of the sample.

2.2 ALMAGAL observation strategy and continuum images

The ALMAGAL program adopted an observational strategy of single pointings centered on the emission peaks identified in the Herschel Hi-GAL 250 μm images (Molinari et al. 2016b; Elia et al. 2017). Each ALMAGAL clump has been observed with a single pointing observation with ALMA, both with the ACA 7 m and the main 12 m antenna arrays. Observations have been carried out using two different configurations for the main array to image all ALMAGAL clumps with a similar physical resolution of ∼1000 au (Sánchez-Monge et al. 2025), although they are located over a wide range of heliocentric distances (Molinari et al. 2025). The C 2+C 5 configuration has been adopted to observe clumps with d_cl ≲ 4.7 kpc, achieving a resolution of ∼0.3^′′, while the C3+C6 configuration for clumps with d_cl ≳ 4.7 kpc for a resolution of ∼0.15^′′. These configurations cover the long and short baselines observed by the main array which, combined with the ACA compact array, allowed the recovery of the flux over all spatial scales up to ∼29^′′, corresponding to the maximum recoverable scale of the ACA array in band 6. In the following discussion, we refer to the group of clumps observed with the two configurations as near and far groups for the targets observed with C 2+C 5 and C 3+C 6 configurations, respectively. Although the two groups were well separated in d_cl when the initial target selection was made, this distinction has been lost with the update of the kinematic distances (see Benedettini et al., in prep.). However, this revision has a marginal impact on the statistical results, as most of the targets in each group are still located in a limited interval of d_cl equal to 2 ≤ d_cl ≤ 4.7 kpc and 4.7 ≤ d_cl ≤ 7.5 kpc, for near and far, respectively (see also Sect. 3.2.2).

The detailed description of data reduction, calibration, and imaging is presented in Sánchez-Monge et al. (2025). The continuum images obtained from the combination of the data from all the arrays have a circular field of view (FoV) with a diameter of ∼36^′′, limited by the primary beam size of the 12 m antenna and the cut to an attenuation of 0.3 applied during the image processing (Sánchez-Monge et al. 2025). In each image, the noise level is rather homogeneous in the central region of the FoV with an average value of ∼0.05 mJy/beam (Sánchez-Monge et al. 2025), but naturally increases towards the edges due to the antenna primary beam attenuation that causes a decrease in sensitivity by a factor of two at angular distances greater than ∼13.8^′′ from the phase center.

2.3 ALMAGAL core catalog

Coletta et al. (2025) extracted compact sources in ALMAGAL continuum images using a modified version of the CuTEx photometric package (Molinari et al. 2011). The entire catalog of identified compact sources is composed of 6348 objects detected at 5 σ in 844 clumps. The nature of these sources can be assessed by their photometric and physical properties. All of these sources have a centrally peaked intensity profile, reflecting the high detection efficiency of the CuTEx algorithm for this morphology and its low sensitivity to extended or flat-profile emission. These sources span sizes between ∼500 and ∼5000 au and have mean densities greater than 10⁶ part cm⁻³ (Coletta et al. 2025). However, most of them have sizes between 1000 and 2500 au (see Fig. 13 in Coletta et al. 2025) and quite high densities between 1.4 and 8 · 10⁷ part cm⁻³. These properties are characteristic of centrally condensed star-forming cores whose innermost regions have undergone collapse. Hence, it is probable that a large fraction of the ALMAGAL sources already host at least one protostellar object. Following this interpretation, we refer to all sources in the ALMAGAL catalog as cores throughout this paper, with the understanding that they likely represent the inner, high-density regions of larger, more extended structures.

We adopted the positions provided by Coletta et al. (2025) to analyze how the cores are distributed in the fields and to measure their relative separations, with the aim of characterizing the fragmentation process in high-mass clumps. Of the 844 systems, 130 clumps host only one core and are not suitable for our analysis. We recognize that the reliability of our analysis and the diagnostics adopted for the spatial distributions are highly dependent on the number of cores. We consider as our full sample all the systems with N_core ≥ 4, excluding the 200 clumps with two or three hosted cores, which are 110 and 90 clumps, respectively (Coletta et al. 2025). The inclusion of these systems has a negligible impact on the results discussed in the paper, and the properties of their separations are consistent with the systems with 4 ≤ N_cores ≤ 6. We adopted the full sample to discuss the distribution of core separation and the comparison with the Jeans fragmentation theory, as systems with a low number of cores may correspond to phases when fragmentation is at the beginning. In contrast, to study the global structural properties of the spatial distribution of the cores, i.e., the type of distribution, mass segregation, we restricted our analysis only to systems with a high degree of fragmentation (N_core ≥ 8 or higher).

In summary, the analyzed dataset comprises 514 clumps (N_core ≥ 4), containing a total of 5728 cores, greatly increasing the statistics available in the literature. However, we identify subsamples with different cuts on the minimum fragmentation level for the various aspects of the analysis of the core spatial distribution. Table 1 summarizes the sample sizes of these subsamples.

3 Results

3.1 The spatial distribution of ALMAGAL cores

The FoV of the ALMAGAL continuum images allows us to probe the fragmentation of massive clumps over circular regions with linear diameters of ∼0.3−1.4 pc in clumps located between 2 and 8 kpc, which are the bulk of the ALMAGAL sample. Visual inspection reveals that the identified cores are arranged in highly diverse configurations. The observed patterns depend on three main factors: the degree of fragmentation (number of cores), the presence (or absence) of extended continuum emission, and the morphology of that emission. The variety in morphologies disclosed by the ALMAGAL data is partially illustrated by the examples presented in Figs. 2 and 3.

We observe a difference between systems where multiple cores are detected and those that are poorly fragmented. In the first case, cores are often grouped together and located in a restricted area of the observed field, typically found on a major patch of extended continuum emission with an irregular morphological appearance. Although there are cases where the emission has a circular shape or modest elongation, it is not uncommon to find cases in which it appears quite elongated (see discussion in Molinari et al. 2025). The distribution of cores typically follows the most prominent continuum emission features, but they are also found outside those boundaries. Cores are also found lying in substructures with a highly filamentary appearance, typical of molecular clouds and present at all spatial scales, as highlighted by Herschel observations (André et al. 2014; Schisano et al. 2014, 2020; Könyves et al. 2015; Dib et al. 2020; Schneider et al. 2022; Hacar et al. 2023). Occasionally, these subparsec filaments are found converging toward the main emission patch as in hub-filament systems (Myers 2009; Hacar et al. 2018; Motte et al. 2018) (see the bottom left panel of Fig. 2). In other fields, the extended continuum emission splits into several well-separated regular patches, each of which hosts a group of associated cores. This is indicative that massive clumps can have an additional higher hierarchical level of substructuring (see the lower panels of Fig. 3), whose presence introduces a subclustering of the core population.

In contrast, systems characterized by a low level of fragmentation do not show a well-defined pattern in the core distribution. These cores are typically sparse in the observed field and do not present a strong signature of grouping or clustering. However, we identified cases where the cores are aligned along specific directions, which extend over several arcseconds. These aligned patterns sometimes match one or a few extended filaments crossing a fraction of the FoV. Even in the cases lacking the detection of extended emission, it may be possible that the hosting filament is undetected for sensitivity reasons and that a signature of its presence is only revealed by the alignment of the cores.

In general, ALMAGAL observations show that the cores are associated with almost every extended region of continuum emission. The high angular resolution reveals that these regions have a wide variety of appearances, as discussed in Molinari et al. (2025). As a consequence of the association, the spatial distribution of star-forming cores shows a similar diversity, likely due to the interplay of different physical factors, such as gravitational force, turbulence, magnetic fields, and feedback mechanisms (Mac Low & Klessen 2004; Krumholz 2006; Li et al. 2014; Shima et al. 2017; Tang et al. 2019; Sanhueza et al. 2025). All of these effects contribute to shape the fragmentation in clumps, affecting both the number and the distribution of cores during their evolution. In this paper, we refer to any group of cores identified in the ALMAGAL clumps with N_cores ≥ 4 as a “cluster”, even when they are limited in number, smaller than the members of typical embedded star clusters.

In the following sections, we introduce quantitative metrics and diagnostic tools to characterize the variety of observed patterns described above. Our aim was to answer the questions of where the cores are formed, how they are distributed, and what clues of the fragmentation process are left imprinted in their spatial distribution. To this end, we analyzed the size and morphology of the observed core distributions, and then focused on the mutual separation between cores. This is an important observable of the fragmentation process, which we related to the clump-averaged properties and compared to the Jeans fragmentation theory. Finally, we evaluated additional diagnostics capable of defining how cores are distributed and how such a distribution changes as a function of the evolutionary phase of the entire system. In addition, we note that a first rough classification of the observed variety of ALMAGAL clusters can be drawn using the convex hull polygon¹(CH) applied to the cores. A simple classification derives from the comparison of the number of vertices of the $\mathrm{C}\mathrm{H},N_v^{hull}$, with the number of cores located well inside the polygon boundaries, $N_{in}^{hull}$ (see Figs. 2 and 3 for reference). About ∼37% of the ALMAGAL clumps studied in this work (190 out of 514) have more cores within the CH boundaries than at the edges $(N_{in}^{hull}>N_v^{hull})$, a pattern that defines a well-localized group of objects with a surrounding sparse population, which is similar to a star cluster. On the contrary, there are 64 clumps (∼12%) where all detected cores form the CH polygon ($N_{in}^{hull}=0$). These systems have a small number of cores, from four to six sparse cores, and typically lack substantial associated extended continuum emission that would be expected if a large amount of high-density material is present. However, most ALMAGAL systems (∼51%) are between these two extreme patterns, since they show only a few cores inside the CH polygon. These systems can be young precursors of a clustered system where fragmentation is still at the beginning, or they can be a population of sparse objects derived from the local fragmentation of the hosting molecular cloud. Nonetheless, it is worth to evaluate the core distribution even in this latter case, as clusters are also recognized to form through mergers of group of cores initially distributed throughout the entire cloud (Chen et al. 2021; Dobbs et al. 2022).

Fig. 2 Examples of the distribution of cores found in four ALMAGAL clumps overlapping on the dust continuum map at 1.4 mm from which they were extracted. These fields show examples of clustered systems with circular and elliptical patterns and cases of aligned cores. The contour levels correspond to 4n × σ level and σ equal to the noise level measured on the image. The positions of the cores from the ALMAGAL catalog extracted by Coletta et al. (2025) are shown with gray star markers. In each field we indicate the cluster geometrical center and its radius (see the definition in the text) with a black cross and a dark blue circle, respectively. We also show the convex hull polygon (cyan segments) and the corresponding best-fitting ellipse (black dashed line) adopted for the morphological characterization of the system. The solid thick black segments are the MST edges determined with Prim’s algorithm (Prim 1957), connecting all the cores in the field.

3.1.1 Center and size of the clusters of ALMAGAL cores

In general, the ALMAGAL cores are located in a portion of the observed FoV, although there are several cases where they spread over the entire, with cores that are also found near the edges. The centers of the core clusters closely match the center of the observed ALMA fields. We adopted two methods to determine the cluster center: (i) the geometrical and (ii) the mass-weighted average of the core positions. These methods produce similar statistical results, although these positions may differ for individual systems. The ALMAGAL cluster centers are mostly located within 3^′′ from the center of the ALMA observed fields, which corresponds to the photometric center of the Hi-GAL clumps as traced in the Herschel 250 μm maps (Molinari et al. 2016b; Elia et al. 2021). The offset between these positions is significant only in a small number of systems: as a reference, it exceeds 6^′′, which is equal to ∼1/3 of the FoV, only in 35 fields (∼7% of the investigated sample).

We estimated the radius of the cluster, R_cluster, as the distance of the farthest core from the cluster center (Cartwright & Whitworth 2004). These radii are smaller than 13^′′ in 384 systems, i.e., 75% of the investigated sample. This region, corresponding to ∼0.72 × FoV, is roughly coincident with the fraction of the images where the noise level is more homogeneous (Sánchez-Monge et al. 2025). However, there are several systems where cores are also detected close to the FoV’s edge; most of those clumps are systems with d_cl<4 kpc. It is likely that for these systems, the single pointing mode of the ALMAGAL observations is sampling only the inner region of the clump. This is a critical selection effect that must be taken into account.

The top panel of Fig. 4 shows the distribution of the cluster radius measured by converting R_cluster into a linear size, $R_{cluster}^{ph}$, using the updated heliocentric distances (Molinari et al. 2025, Benedettini et al., in prep.). The distribution for the full sample shows two key features: a peak close to the median value of ∼0.2 pc, and a notable drop at $R_{cluster}^{ph}\approx 0.32\mathrm{p}\mathrm{c}$. To determine the impact of the FoV bias on the true physical distribution, we selected the 347 clumps with d_cl ≥ 3.7 kpc. As shown in the bottom panel of Fig. 4, for those systems, the observed area is sufficiently large to identify clusters with radius $R_{cluster}^{ph}\ge 0.32\mathrm{p}\mathrm{c}$. In this sample, the peak at ∼0.2 pc disappears, revealing a distribution that is rather uniform in the interval $0.1\le R_{cluster}^{ph}\le 0.32\mathrm{p}\mathrm{c}$. However, the cut-off at 0.32 pc remains even in this selected sample, suggesting that it is a true physical upper limit to the cluster size. We find only 36 systems with $R_{\text{cluster}}^{\text{ph}}\ge 0.32\mathrm{pc}$ (out of the 302 clumps with 3.7 ≤ d_cl ≤ 6.5 kpc), a number that is sufficiently small to not represent a significant statistic. Our conclusion is that the clumps with 3.7 ≤ d_cl ≤ 6.5 kpc provide the most robust representation of the true cluster size distribution in a generic massive dense clump. As a consequence, it is quite likely that there would be cores located outside the observed area in the 168 fields with N_core ≥ 4 and d_cl ≤ 3.7 kpc, meaning that the measurements of the core separations can be potentially incomplete in these systems, while they offer a better snapshot of the central region of the clump.

Our estimates of $R_{cluster}^{ph}$ can be affected by uncertainties. First, this definition of cluster radius depends on the effective membership of the farthest core, which is not always reliably confirmed. Secondly, there is a decrease in sensitivity towards the image edges, meaning that cores with low fluxes located at the margin of the FoV may remain undetected. It is interesting to compare these estimates with the size of their hosting clump, which defines where there is an overdensity of material that is potentially available for star formation.

Ideally, it would be useful to compare $R_{cluster}^{ph}$ with the extension of the continuum emission measured on the same ALMAGAL data. However, as discussed by Molinari et al. (2025), ALMAGAL interferometric data are affected by spatial filtering, and in several fields, they may miss a significant flux fraction. This makes ALMAGAL data unsuitable for reliably tracing the full, extended emission of the clump, even if they reveal the complexity of the clump inner structure. For this reason, we adopted the clump sizes, R_cl, determined from Herschel 250 μm data and defined as the circularized half-width at half-maximum of the elliptical Gaussian fit with CuTEx (Molinari et al. 2016b; Elia et al. 2021). This definition is certainly limited and a simplification of the complexity of the clump structure. However, our goal here is not to define an exact physical boundary, which is also not well defined in molecular clouds. Rather, our aim is to use a consistent and robust proxy for the region containing the bulk of the clump mass. We note that ALMAGAL single pointings are sufficiently wide to probe a wider area than the Herschel beam at 250 μm (with a full-width half maximum of ∼18^′′).

The size of the clumps ranges from ∼0.1 to ∼0.9 pc, and most of them have R_cl ≤ 0.35 pc, as shown in Fig. 5, which is comparable to our measurements of the cluster radius. The distribution of the ratio of $R_{cluster}^{ph}$ to R_cl reported in Fig. 5 indicates that the cores are typically located between ∼0.3 and ∼1.2 × R_cl. Hence, the formation of cores takes place in any fraction of the hosting clump, and may even extend to its close outskirts, but rarely beyond the overdensity that has formed in the molecular cloud. Although it is possible that additional cores may reside farther away from the center of the clump, they must have low mass to remain undetected towards the edge of the FoV. As the sensitivity in this area drops by a factor of ∼3, these undetected cores would have masses below ≤ 0.6 M_⊙, adopting the completeness limit indicated by Coletta et al. (2025).

Finally, we found that our measurements of $R_{cluster}^{ph}$ do not show any correlation with the clump surface density Σ_cl or the L/M ratio as shown by Fig. 6. A similar result is found for the ratio $R_{cluster}^{ph}/R_{cl}$.

Fig. 3 As in Fig. 2. These fields show examples of cores distributed over filamentary features and in well-separated local substructures, both patterns that are sometimes found in ALMAGAL continuum images.

Fig. 4 Top panel: distribution of the cluster radius measured from the average core positions for the sample of 514 ALMAGAL clumps with Ncore ≥ 4 (in gray) and the subsample composed 347 clumps with dcl ≥ 3.7 kpc (magenta). Bottom panel: clump heliocentric distance vs. measured cluster radius. The dashed line shows the physical size of the ALMAGAL FoV, while the dotted line indicates the area where the sensitivity is constant within a factor of two, i.e., angular sizes ∼13.5′′. Cluster radii exceeding the FoV sizes are possible in systems with a substantial offset between the cluster and the image center.

Fig. 5 Distribution of the ratio of the radius of the core cluster to the clump radius (gray histogram). The dot-dashed line, which refers to the top and right axes, shows the distribution of the clump radius derived from Herschel Hi-GAL data at 250 μm from Elia et al. (2021).

3.1.2 Morphology of the clusters of cores

The definition of cluster radius introduced in the previous section provides a robust measurement for spherical systems, for which the projected circular area $A=\pi R_{cluster}^2$ is a valid approximation for the size of the entire cluster. Nonetheless, cores are often distributed over elongated or irregular regions, for which the circular area gives a considerable overestimate of the cluster size (see, e.g., the bottom panels of Fig. 2). The CH polygon allows us to better estimate the cluster area and determine the occurrence of spherical, elongated, and irregular clusters. This method is widely adopted in studies of embedded clusters in star-forming regions (Schmeja & Klessen 2006; Gutermuth et al. 2009; Dib & Henning 2019) and it was recently adopted in Molinari et al. (2025) to also characterize the extended continuum emission in the ALMAGAL fields that shows significant morphological diversity.

We note that the observed variety in the spatial distribution of cores is not unexpected, since similar results are also found when observing stellar systems, of which the ALMAGAL cores are the precursors. In fact, Kuhn et al. (2014) observed different appearances in young star systems, with the most recurrent morphologies represented by chains of subclusters, clumpy ones, and isolated groups. The elongation of the cluster is a simple metric that can identify these possible patterns. To measure the cluster elongation, we fitted an ellipse to the area covered by the CH, measuring its orientation, PA_CH, and the ratio between the major and minor axes, e_CH. A different method was proposed by Schmeja & Klessen (2006), who introduced the quantity $$\xi =\frac{R_{cluster}}{R^{\mathrm{C}\mathrm{H}}},$$(1) where R_cluster is the circular cluster radius and R^CH is the convex hull radius, which is defined as the radius of the circle with an area equivalent to CH, A^CH. We note that Schmeja & Klessen (2006) adopted an additional correction factor $1/(1-N_v^{hull}/N_{cores})$) to statistically recover the real area on which the sources are distributed (Hoffman & Jain 1983) (see also the discussion in Appendix A of Parker (2018)). However, such a factor was found to be reliable only when ∼200 points are available (Ripley & Rasson 1977), while in our case it would introduce large discrepancies since $N_v^{hull}$ may be very close to N_cores. Therefore, we computed R^CH without applying this factor, although we also report the CH area normalized by it in Table 2. The distributions of the two estimators e_CH and ξ are similar, with similar average values. The distribution of ξ is narrower as ξ has a smaller dynamic despite the irregular shape of the CH. The most common shape of the cluster is a slightly elongated ellipse, with a median e_CH ∼2.2[ξ ∼2.3]. While we found systems in which the cluster has a rather circular shape (e_CH ≤ 1.3), they represent a minority with only 46 cases (∼9%) of 514. Other systems have strong elongations as high as ≈ 20−30, depending on the adopted estimator.

Approximately circular shapes may result from a distribution where cores are uniformly distributed in the system, but also from patterns where a few filamentary features are converging towards a central hub and are symmetrically distributed over different directions, as in the top left panel of Fig. 2. Figure 7 shows the distribution of the measured elongations for the clusters classified into the three groups described in Sect. 3.1. Similar results are obtained if the classes are defined in terms of the fragmentation level. High elongations are measured in systems with a limited number of cores and include the cases where they align along specific directions or are distributed over extended filamentary structures. Clustered systems typically have smaller elongation, with a median and interquartile range equal to $\xi =1.9_{-0.3}^{+0.4}(e_{\mathrm{CH}}=1.9_{-0.4}^{+0.6})$. As they are representative of systems with a high degree of fragmentation, this means that their cores are typically distributed over half of the circular area ascribed to the cluster, in slightly flattened structures like the ones shown in the bottom left panel of Fig. 2. However, it is worth noticing that a fraction of ∼7−16% (depending on the adopted estimator) of the 190 systems classified as clustered have elongations greater than three. While these high values suggest rather elongated morphologies, they also derive from peculiar distributions of the cores with highly asymmetric patterns. An example is shown in the bottom right panel of Fig. 2, where a group of cores is aligned along a direction converging towards an extended and rather circular emission causing a large elongation due to their separation from the main system.

Fig. 6 Cluster radius as a function of the clump-averaged properties: L/M (left panel) and surface density Σcl (right panel), color-coded by the clump distance. The black lines indicate the contours of the 2D probability density distribution corresponding to 10, 25, 50, 75, and 95% levels.

Fig. 7 Whisker plots showing the distribution of the two measured estimates of the cluster elongation for the different classes of ALMAGAL clumps identified by the number of cores inside the CH polygon: sparse, where no further core lies inside it; intermediate, where there are fewer cores than the ones defining the polygon itself; clustered, where there are more cores than the ones at polygon vertices. The left panel presents the elongation eCH, defined as the ratio of the major to the minor axis of the best fitting ellipse to the CH polygon, while the right panel refers to the ratio ξ of the circular area of the cluster to the area of the CH polygon. The first and third quartiles, and the 10th and 90th percentiles of the distribution are indicated by the boxes and the segment extremes, respectively, while the central segment corresponds to the median.

3.2 Core separations

3.2.1 Methods

Several mathematical tools are available to characterize the spatial distribution of a set of localized points, which are representative of the members of a cluster. In particular, an important metric to measure is the mutual separations between those points, whose distribution may be related both to the fragmentation scales and to the degree of clustering (Cartwright & Whitworth 2004; Schmeja et al. 2008; Sánchez & Alfaro 2009; Allison et al. 2009; Parker 2018; Dib & Henning 2019). Those tools are based on the computation of the distance matrix δ(i, j), a symmetric matrix whose elements are the angular distances between the points taken pairwise. Given N sources, the distance matrix δ(i, j) includes all possible $\frac{(N-1)\times N}{2}$ connections and fully characterizes the relative separations for all elements of the cluster. In this work, we do not require the complete information provided by δ(i, j), as, to describe the length of fragmentation, it is sufficient to relate each core to its direct neighbor. Those separation are usually computed with the minimum spanning tree (MST) method.

The MST is a graph that connects all cores of the system and corresponds, among all the possible connecting paths, to the set of segments with the shortest total length. It is composed of (N−1) segments, named MST edges, that are a subset of the δ(i, j). A simple way to estimate the typical separations between the cores is the statistical analysis of the distribution of the MST edges.

We computed the MST with an IDL implementation of the Prim’s algorithm (Prim 1957), identifying in the distance matrix δ(i, j) the edges l_i that compose the MST. Prim’s algorithm starts from a single position, adopted as a vertex, and identifies the shortest among all the possible connections with the adjacent sources, then proceeds iteratively until all the sources are finally connected, avoiding to form closed loops in the path. Examples of the computed MSTs are presented in Figs. 2 and 3. We identified a total of 5214 edges, l_i, measured in all 514 ALMAGAL fields analyzed in this study.

We followed two different approaches in the analysis of such a large ensemble² of separations. On the one hand, we independently analyzed the distribution of the separations in each clump. We measured statistical quantities, such as the median and mode of distribution, to provide an estimator for the characteristic core separation in each system. On the other hand, we consider the entire ensemble as a diagnostic of the distributions of core separations, relating all the measured separations to the clump physical properties. In such a case, the statistical estimators are determined over subsamples obtained by stacking all the separations measured in clumps with similar values of a specific property, such as the L/M ratio, ignoring that those clumps are distinct structures where possible different external conditions may be present. These two approaches allow for a well-rounded and complete description of the statistical information present in our measurements and to determine the possible link to the clump-averaged properties.

With the first approach, we provided estimates for the core separations that are clump-averaged, similar to other physical properties provided for clumps, such as density or temperature. These estimates also highlight the information provided by poorly fragmented systems, which would be hidden in the entire ensemble, as they are numerically outweighed by the contribution of clumps with several cores. Instead, with the entire ensemble, we provide an overall snapshot of the fragmentation process and relax the assumption of a unique characteristic core separation for clump. This analysis method provides a better statistical characterization thanks to the size of the sample and a more robust handling of outlier cores, whose presence would affect the clump-averaged estimators, especially for systems with a low degree of fragmentation.

3.2.2 The distribution of core separations

We divided the entire ensemble according to the set of array configurations adopted for the ALMAGAL observations, i.e., the near and far groups, respectively. The reason for this choice is twofold. First, the FoV is resolved by these two configurations in a different number of synthesized beams, meaning that these two datasets could potentially provide a substantially different density of point-like sources. In more detail, the FoV of the combined ALMAGAL continuum images is covered by ∼2000 and ∼8100 synthesized beams with the C 2+C 5 and C 3+C 6 configurations, respectively. Secondly, this division allows us to control the impact of the clump distance, as these two groups are composed of clumps located mainly in two different ranges of heliocentric distances. The two subsamples are composed of 4012 and 1716 cores, identified in 330 and 184 clumps, and provide 3682 and 1532 separations, respectively, in the near and far groups. The discrepancy in the number of cores is mainly due to the sensitivity limit of the survey as discussed by Coletta et al. (2025).

The right panels of Fig. 8 show the distributions of the projected linear separations between the cores given by the MST edges, l_i, obtained with the revised heliocentric distances for the ALMAGAL clumps. Most of these separations (∼70%) exceed the spatial resolution of the data by a factor of ∼3. In general, they range from ≲ 1000 to ∼100 000 au (corresponding to ≲ 0.005−0.5 pc), with the median of the two distributions equal to ∼6600 and ∼7400 au for near and far sample, respectively. Although the two distributions have a similar median and first quartile, they do not derive from the same underlying distribution, as indicated by the Kolmogorov-Smirnov (KS) test rejecting this null hypothesis with a p ≲ 10⁻⁵. The major difference between the two distributions is given by the large tail present in the far sample due to an increase in the number of separations with l_i > 15 000 au, caused by the wider physical area covered by the FoV in those clumps. Those separations produce a secondary bump in the distribution at ∼30 000 au.

We present the projected linear separations as a function of the clump distance in the left panels of Fig. 8 to compare them with the spatial resolution and the FoV of the observations. Clumps with heliocentric distances greater than 9 kpc are omitted in the plot for better visual clarity, but are included in the analysis. We found separations between cores ranging from the beam up to half the FoV’s size, indicating that the enlargement of the observed area includes additional cores located at the outskirts of the system that are more widely separated. The inclusion of these cores introduces a weak trend in the median and mode of the distribution, calculated in bins of width 500 pc, as also shown in Fig. 8. The mode increases steadily with the distance for both samples and varies from ∼3000−3500 au at a distance of 2 kpc up to ∼7000−7500 au for d_cl ≈ 4 kpc. A similar trend, but with slightly larger values, is found for the median of l_i, which varies from ∼4000 to ∼9000 au and ∼7000 to ∼8000 au for the near and far group, respectively. The median and mode of the distributions of l_i have the largest variations for distance d_cl<3.7 kpc, where we most likely observed only the inner region of the cluster, as discussed in Sect. 3.1.1. This suggests that typical separations in the central portion of these clusters may be shorter, with differences in the estimated characteristic lengths in the densest regions of the clump. The gradual shift of the median towards larger values is ascribed to the inclusion of the cores that are more widely separated. The variations of the mode are explained by the fact that the inclusion of cores with larger separations extends the interval of the data analyzed with the KDE method. As the interval increases, the KDE computation adopts a wider bandwidth, with the effect of shifting the position of the peak of the continuous distribution to larger values. Finally, the distributions are affected by the sensitivity limit of the observations, since in the clumps located closer to the Sun it is possible to identify fainter cores. The ALMAGAL catalog is complete for objects M ≳ 0.2 M_⊙ (Coletta et al. 2025), but less massive cores are regularly detected in clumps with d_cl ≤ 4 kpc. To deal with this effect, we discuss separately the near and far groups. Moreover, we also checked whether our results holds when selecting only the clumps located in narrow distance intervals.

Fig. 8 Left panels: linear separation between ALMAGAL cores as a function of the distance of the ALMAGAL clump for near (top panel) and far configuration sample (bottom panel). The gray shaded areas indicate the ranges of the limiting sizes set by the beams of the ALMAGAL images. The black and colored solid lines (blue and green) indicate the median and the mode of the distribution of the separations, respectively, while the colored area traces the distribution width, defined at the 68% peak intensity. The dashed and dotted lines indicates one-half of the observed field and the full field, respectively. The vertical lines show the adopted distance bins with the relative number of the separations. The 20 separations lying below the spatial resolution are found where the beam is strongly elongated, and are larger than the beam minor axis. In addition, a limited number of separations are not included in the figure as they exceed the selected intervals: 180 in the near sample (162 of which have dcl ≥ 9 kpc), and six in the far sample. Right panels: overall distribution of the projected linear separations in the near (top) and far (bottom) sample, respectively. The thick colored lines (blue and green) show the probability distribution function estimated with the KDE method. The gray shaded area indicates the interval spanned by the circularized beams in the corresponding ALMAGAL images.

3.2.3 Deprojection of the measured separations

The measured separations correspond to the projection in the plane of the sky of the intrinsic 3D separations between cores hosted in the clump. The distribution of 3D separations can be statistically derived by evaluating the effect of the projection on any generic vector randomly oriented in the 3D space. The details of the calculations are reported in Appendix A. On average, the projection along the line of sight of a generic vector oriented in the 3D space reduces its intrinsic length by a factor equal to π/4 (see Eq. (A.4)). As a first approximation, we corrected the measured separation between each core pair by this factor to scale to an average deprojected value, which we consider as the intrinsic separation $l_i^{depr}$. In such a case, the distributions of $l_i^{depr}$ retain the shape of the projected separations as we simply applied a scaling factor between the two quantities.

3.2.4 Characteristic separations between cores

The variety of observed patterns is reflected in different distributions of the separations between the hosted cores in each ALMAGAL field. To provide a representative characterization of these distributions, we determined the minimum, maximum and median of the separation l_i, which we report in Table 3. We also included the mode of the distribution to indicate the most frequently observed separation, which is particularly relevant for systems where one or more groups of cores are surrounded by a more widespread population. In these systems the median would provide a bias estimator, which may not be representative of the typical observed separation. In each clump we estimated the continuous probability density function $P^k\left(l_i^k\right)$ of the separations of the hosted cores using the Gaussian KDE method implemented in the scipy library. We identified the peak position of $P^k\left(l_i^k\right)$ to determine the mode, l_mode, and the spread of the distribution, the latter defined as the width of P^k (l_i) at 68% of the intensity peak, avoiding any assumption of the shape of the PDFs, which we observed to be often quite asymmetric. These quantities are measured for the 437 clumps hosting more than 5 cores, which provide at least four measurements for the core separations, and are included in Table 3. Statistically deprojected values in the k-th clump, $l_i^{depr,k}$, can be obtained by scaling by the constant factor 4/π (see the discussion in Sect 3.2.3).

Our estimates of the characteristic, clump-averaged, separation are shown in Fig. 9. They vary from ∼1500 to ∼60 000 au, and values larger than ∼15000 au are found mainly in clumps with a low level of fragmentation. A decrease in the average separations in highly fragmented systems is expected for simple geometrical reasons, but we verified that the observed trend is not caused by the reduction of the available area per source. In fact, we compared with the expected decrease in separations when sources are uniformly distributed over a field with a fixed angular size and found that the observed averages are sensibly smaller than these expectations.

However, we note that the clump-averaged estimates provide limited information about the distribution of separations in a clump. Individual distributions $P^k\left(l_i^k\right)$ have a diversity of shapes, often presenting multi-modal patterns, features that are not easily captured by a single averaged quantity. A thorough statistical analysis and classification of the distribution $P^k\left(l_i^k\right)$ will be conducted in a future study. Here, we compared l^{median, k} and l^{mode, k} in Fig. 9 to report the occurrence of asymmetric distributions. These two quantities agree within 25% for 316 systems (∼72%), and represent cases where P^k(l^k) is rather symmetric around a central value equal to the measurements provided. However, noticeable differences between the two estimators are not rare. In systems where they differ, the distributions of core separations P^k(l^{depr, k}) are sensitively skewed and asymmetric. Quantitatively, the median separation is larger than the mode in 96 systems (∼22% of the sample), with such a pattern found in both near and far systems, with 59 and 37 cases, respectively. The opposite case, with a median smaller than the mode, is rarer and is found in only 25 systems (≤ 5% of the sample).

Asymmetric distributions can be caused by the complex morphology of the 3D clump structure and the spatial distribution of cores. However, it is important to note that it is extremely unlikely to obtain multi-modal distributions for the projected separations starting from a single peaked distribution for the intrinsic 3D separations. One possibility is that these patterns are due to the limited degree of fragmentation observed, which could not be sufficient to properly sample the true distribution. Alternatively, they may also be the signature of different fragmentation conditions in groups of ALMAGAL clumps, with diverse characteristic separation lengths. We observed cases where a fraction of the core population is located farther away than the main group of more tightly packed cores. Indications that multiple characteristic lengths may be present in the fragmentation of massive clumps were already found in the literature (Zhang et al. 2021).

Fig. 9 Relation between the median of the core separations in each clump and the mode determined with the KDE method. The colors indicate fields with different number of cores, ranging from 4 to 50 in logarithmic scale. The dashed line indicates the 1: 1 relation, while the dotted lines indicate the range of variation of ± 25%.

3.3 Relation between the core separations and the clump properties

The core positions directly probe the fragmentation of the clump as they indicate where the dense material collapsed. The distribution of patterns and core separations reflects the diversity of initial conditions or evolutionary paths that are observed in the ALMAGAL systems. This wide observed diversity has also been highlighted by Elia et al. (2026), who related the properties of the core population to the ones of the hosting clumps finding a significant scatter in all their relations, interpreted as a result of different fragmentation conditions (Elia et al. 2026). Nonetheless, they also identified statistical trends between some fragmentation metrics. For example, they found that the number of cores, N_cores, increases with Σ_cl and the evolutionary indicators L/M ratio and T_bol. Similarly, here we discuss how the separations between cores relate to the clump-averaged properties, using the values derived from the Herschel Hi-GAL data (Molinari et al. 2025) (see Sect. 2.1).

3.3.1 Correlation with the clump evolution

The formation of cores is a consequence of the collapse of the dense material, which proceeds over timescales that are a few times the free-fall time of the gas. The degree of fragmentation changes with time, as more clump material may become gravitationally unstable, depending on the local conditions, and cores evolve into stars. The local conditions vary over the entire clump evolution because of the internal dynamics of the clump or the feedback of the newly formed stars. Changes in the separations between cores during the clump evolution, if present, are useful to trace the variations of these conditions, determining how fragmentation proceeds.

We show in the left panels of Fig. 10 the normalized density distributions of the deprojected separations measured in clumps divided into three main classes according to the value of the ratio L/M³, which is usually adopted as a diagnostic for the evolutionary state of the clump (Molinari et al. 2016b) (see Sect. 2.1). The classes were chosen to separate the clumps in the youngest phases of their evolution, which we considered as the systems with a ratio L/M<1 L_⊙/M_⊙, from those that are highly evolved and host multiple protostars, which we selected here as those with a ratio L/M>10 L_⊙/M_⊙. These values match the intervals identified by Molinari et al. (2016a) for different phases of clump evolution, and allow us to also compare our results with those of the surveys ASHES (Sanhueza et al. 2019), CORE (Beuther et al. 2018), and ASSEMBLE (Xu et al. 2024), each of which observed systems in a specific evolutionary phase. The two classes of young and evolved systems are composed of 129 and 209 clumps, respectively, leaving 176 clumps in the intermediate one (see Table 4). Figure 10 shows that the distribution of the young class, with L/M<1 L_⊙/M_⊙, and the evolved one, with L/M>10 L_⊙/M_⊙ are substantially different. The peak and median decrease to smaller values in the population of evolved clumps, which is also characterized by a narrower distribution. The Kolmogorov–Smirnov and Mann-Whitney U tests highlight the differences of these distributions, since they both return a p-value ≪ 10⁻⁷, rejecting the null hypothesis that these two datasets derive from the same underlying population. Specifically, in young systems with L/M<1 L_⊙/M_⊙, the near and far distributions have median values equal to ∼12 000 and ∼18 000 au (∼0.06−0.085 pc), respectively. However, these values describe only partially the distributions presented in Fig. 10, which have complex shapes. The distribution peaks estimated with the KDE method correspond to ∼14 000 and ∼34 000 au, respectively. These values overlap with the separations equal to ∼14 000−35 000 au measured in the 39 protoclusters observed by the ASHES survey (Sanhueza et al. 2019; Morii et al. 2024). We note that while the ASHES survey maps a larger angular area than ALMAGAL’s single pointings for near clumps, the two surveys cover a similar physical area in the more distant far sample, for which the comparison is less biased by the observed area. The distribution of the latter sample is broad and negatively skewed, with a peak that matches the upper values reported by Morii et al. (2024). An additional comparison is provided by the 12 young starless clumps observed by Svoboda et al. (2019) where the median of the measured projected separations is 0.083 pc, corresponding to a deprojected value of ∼22 000 au, which is marginally larger than the median values reported above, but still compatible with the ones measured in the youngest systems of our ALMAGAL sample.

In the evolved clumps with L/M>10 L_⊙/M_⊙, the separations decrease to median values of ∼7300−7700 au (∼0.035− 0.038 pc). These values match the results reported by the CORE (Beuther et al. 2018) and ASSEMBLE (Xu et al. 2024) surveys, which targeted the high-mass and high-luminosity clumps characterized by 10 ≤ L/M ≤ 1000 L_⊙/M_⊙ and that are in the later phases of their evolution. The median of the separations measured in the 20 clumps of the CORE survey is equal to ∼6000 au, with a interquartile interval varing from ∼4000 au to ∼15 000 au, once their measurements are corrected for the projection effect (see Sect. 3.2.3). Even if this value is marginally smaller than our measurements, it is considered compatible given the widths of the distributions presented in Fig. 10. The typical separation measured in the 11 high-mass clumps of the ASSEMBLE survey is equal ∼0.035 pc (∼7200 au) which is in agreement with our results.

A considerable fraction of the ALMAGAL clumps have a L/M ratio between the two values discussed above. In these systems, we find a different distribution of the separations for the clumps of the near and far sample. In the near sample, the distribution is similar to that of the evolved systems, the opposite holds for the far sample, where it resembles that of the young clumps. We measured median values between ∼9000 and ∼15 000 au, with the spread influenced by the wide area probed in the far sample. However, the KS test rejects the null hypothesis that each of the two sets of distributions derives from a unique underlying distribution with p ≲ 10⁻⁷ and p ≲ 0.05 for the near and far sample, respectively. We mainly attribute the different shapes of the two distributions to the different degrees of fragmentation of their clump population. The near sample is dominated by highly fragmented clumps with N_cores>8, which strongly impact the shape of the distribution since their separations contribute for more than ∼80% of the total. Conversely, the far sample consists mainly of poorly fragmented clumps (N_cores ≤ 8), whose separations are a relevant (∼65%) fraction of the distribution presented in Fig. 10. In summary, the different shapes are produced by the different relative contributions of poorly and highly fragmented systems combined with the evidence that separations are typically larger in clumps with N_cores ≤ 8 than in those with a higher degree of fragmentation. A minor, secondary effect may be caused by the larger physical area observed in the far group, which favors the measurements of wider separations.

Only a few clumps with 1<L/M<10 L_⊙/M_⊙ were observed before ALMAGAL, a large fraction of them were obtained by the SQUALO (Traficante et al. 2023) and the DIHCA survey (Ishihara et al. 2024). A typical separation of ∼7000 au is reported from the 30 clumps of the DIHCA survey, although this includes systems with higher L/M, since their sample covers ∼1 to 50 L_⊙/M_⊙. The average separation and, in general, the entire distribution closely match the results from the more evolved clump populations, similar to what is found for the near sample of the ALMAGAL survey. In contrast, the comparison with the SQUALO survey is severely limited since only 3 clumps (out of 13 surveyed in total) fall within this interval of L/M (Traficante et al. 2023). Traficante et al. (2023) reported a low level of fragmentation, measuring minimum separations between 13 000 and 29 000 au. Although these values should be considered as a lower limit because of the projection effect, these wide separations are nonetheless comparable with our results from the far sample systems.

In summary, our statistical analysis reveals an evolutionary trend in core separation. When ALMAGAL clumps are grouped by their L/M ratio, we find that the cores in more evolved systems are generally closer to one another than those in younger systems. Evidence of a relationship between core separations and the evolutionary diagnostic L/M was already proposed by Traficante et al. (2023), who compared the separations in the 13 young and moderately evolved targets of their SQUALO survey to those of the 20 high-mass and high-luminosity clumps of the CORE survey (Beuther et al. 2018). Similar indications were also found by including the 30 moderately evolved systems (∼1<L/M<50 L_⊙/M_⊙) of the DIHCA survey (Ishihara et al. 2024) and comparing the 11 massive, luminous, and evolved clumps of the ASSEMBLE survey (Xu et al. 2024) with the 39 high-mass young clumps observed by the ASHES survey (Sanhueza et al. 2019; Morii et al. 2024).

The sample size and the uniform coverage across the L/M parameter space of the ALMAGAL survey provide robust statistical significance to this trend. It allows us to move beyond a coarse three-group comparison and analyze in detail the variation as a continuous function of L/M. The right panels of Fig. 10 show the entire ensemble of separations as a function of the clump L/M ratio. Although the separations $l_i^{depr}$ have a wide scatter for any given ratio L/M, we find a steady decreasing trend traced by the median of the distributions computed in equally spaced bins in logarithm of L/M. This behavior is not caused by the weak correlation found with distance. In fact, it is still present when selecting systems in the two distance intervals 3.7 ≤ d_cl ≤ 4.7 kpc and 4.7 ≤ d_cl ≤ 6.5 kpc, for which the statistical values of the separations do not show any dependence on distance (see Sect. 3.2.2), despite the minor quantitative differences emerging due to the smaller sample sizes. Apparently, the two samples near and far show different trends. In the near sample, we found that the median has a major decrease from ∼15 000 to ∼10 000 au occurring until L/M ≈ 1 L_⊙/M_⊙, followed by a minor decrease to 7000−7800 au for increasing values up L/M ≈ 100 L_⊙/M_⊙. In the far sample, the median is typically larger in the early phases, with a value of ∼20 000 au, and shows a steeper and more steady decrease starting at L/M ≈ 1 L_⊙/M_⊙. These differences are due to the smaller sample size, lower level of fragmentation, and wider observed area in clumps of the far sample, as discussed above. The shift of the median to larger values is caused by the increase in the statistical weight of clumps with a limited degree of fragmentation whose cores are more distantly separated, as discussed above. It is interesting to compare the trend traced by the mode of the distributions. The right panel of Fig. 10 shows that in the two samples the modes are coincident for L/M ≳ 1 L_⊙/M_⊙. Instead, we observe a match between the values of the near sample and the shapes of the 2D density contours for L/M<1 L_⊙/M_⊙. These shapes suggest a potential bimodality in the distribution, with one of the peaks that follows the trend found in the near sample. More generally, the distributions of the separations in the far sample are all broad for clumps with L/M ≤ 10 L_⊙/M_⊙. Moreover, the profiles of the KDE distributions shown in the left panel of Fig. 10 present two local maxima, one approximately at ∼40 000 au and another that shifts from ∼10 000 down to ∼7000 au. The maximum at large separations is built mainly, but not exclusively, by systems with a low degree of fragmentation, whose cores are more widely distributed. In fact, about ∼72% of the separations greater than ∼20 000 au of the far sample are measured in clumps with N_cores ≤ 8. To verify that these maxima correspond to a statistically significant bimodality, we ran the Silverman (1981) test, following the implementation described in Beraldo e Silva et al. (2021, see their Appendix A), applied to the distributions in different intervals of L/M. In a few intervals of L/M, we were able to assess that the corresponding distribution of separations has more than one significant peak. However, we were unable to extend this to the entire L/M range due to the limited sample size in a few intervals of L/M. In general, we found that the null hypothesis that at least two peaks are present can be assessed only with a probability Prob > 0.8, although it becomes significant (Prob > 0.95) over limited intervals of L/M.

Finally, we present in Fig. 11 the same analysis using clump-averaged separations, showing only the median separation as a function of the L/M ratio since the mode has a similar trend. Similarly to the case of the entire ensemble, the averaged quantities also show a gradual decrease in the separation between the cores with L/M, as we determined a Spearman correlation coefficient equal to ρ_s=−0.37 [−0.38] for the median [mode] with a p-value ≪ 10⁻⁹. These values indicate that a correlation is indeed present; however, the Pearson correlation coefficient of ρ_p=−0.26[−0.28] suggests that it is not a single-linear relation over the entire interval of L/M. Figure 11 also indicates that the observed relationship is connected to the degree of fragmentation, which also statistically increases with the L/M ratio (see the discussion in Elia et al. 2026). As discussed for the case of the entire ensemble, the relationship changes for L/M ∼1 L_⊙/M_⊙, above which it becomes rather flat. This value sets a boundary above which ∼90% of the fields are intermediate or clustered, according to the classification introduced in Sect. 3.1, or, alternatively, have a moderate to high degree of fragmentation, as they mostly host more than eight cores.

Fig. 10 Left panels: density distribution of the deprojected core separations measured in the all clumps with L/M ≤ 1 (orange), 1 ≤ L/M ≤ 10 (green), and L/M ≥ 10 L⊙/M⊙ (purple) with the continuous distribution computed with the KDE method overlapping with the colored continuous line. Clumps of the near and the far are shown in the top and bottom panel, respectively. The colored vertical thick lines show the median of these distributions. Right panels: Deprojected core separations as function of the clump L/M for the systems in the near (top panel) and the far (bottom panel) subsample, respectively. The thick black lines indicate the median of the distributions of the core separations computed in bins of L/M of width 0.3 dex. The thick colored lines, in blue and in green for the near and far subsample, respectively, indicate the peak of these distributions, while the colored area traces their spread defined as the width at 68% of peak intensity. The dashed lines show the variation of the mode (in light blue and green) and of the median computed over subsamples having 3.7 ≤ dcl ≤ 4.7 kpc and 4.7 ≤ dcl ≤ 6.5 kpc. The numbers at the top give the size of the subsamples from which the distribution in each bin is computed. The colored contours trace, from the outside inward, fractions of the sample lying outside that rise from 15% to 95%, with a step of 10% for the near (in blue) and far (in green) subsample, respectively. The median and the peak derived for the near sample are also indicated in the bottom panel for comparison purposes.

Fig. 11 Relation between the clump-averaged separations between the cores provided by the median and the L/M ratio, color-coded to indicate the number of cores detected in the clump. To support the visualization, different sizes are adopted to identify systems 4<Ncores<8 (low), 8 ≤ Ncores<16 (intermediate), and Ncores ≥ 16 (high). The different symbols refer to the clumps of the near (circle) and far (diamond) sample. The solid lines show the median of the distribution of the separations computed in bins of L/M of width 0.4 dex for the near (blue) and far (green) sample, respectively, while the black dashed line corresponds to the median of all the clumps with Ncores ≥ 4.

3.3.2 Relation with clump surface and volume density

The ALMAGAL clumps cover systems with surface densities, Σ_cl, varying from ∼0.1 to ∼10 g cm⁻², which allow us to determine whether there is a connection between the places where the material collapsed and the averaged density of the clump. Density is a key quantity in the fragmentation theory and, in the case of Jeans fragmentation, sets the scales that ultimately shape the separation between cores. Figure 12 shows the distribution of all the separations as a function of the clump surface, Σ_cl, and volume number densities, n_cl for the near and far samples. The clump-averaged Σ_cl and n_cl are computed with the assumption of spherical geometry for the clump as $${\Sigma }_{cl}=\frac{M_{cl}}{\pi R_{cl}^2}\text{and}{n}_{cl}=\frac{3M_{cl}}{4\pi \mu m_H{R}_{cl}^3},$$(2) where M_cl and R_cl are the clump mass and radius, respectively, m_H is the atomic mass and μ is the mean molecular weight per free particle adopted to be equal to 2.37, which corresponds to a molecular gas with solar composition⁴ (Kauffmann et al. 2008).

Quantitatively, the typical separation between the cores steadily decreases from ∼25 000 au to ∼8000 au in clumps with an increasing Σ_cl from ∼0.3 g cm⁻² to ∼5 g cm⁻², and a drop to ∼6000 au for Σ_cl ≥ 5 g cm⁻², where the relation appears to be flattening. A decreasing trend is also found when the separations are compared to the clump volume density n_cl, with a substantial decrease between ∼10⁴ and ∼10⁶ cm⁻³ (see the lower panel of Fig. 12). Both relations show a wide scatter of the measured separations, and weak linear relations appear only from the statistical estimators, i.e., the median and the mode. We checked for possible power-law relations between l^depr versus Σ_cl and n_cl by fitting to the data with the functional forms $l^{depr}\propto {\Sigma }_{cl}^{\alpha }$ and $l^{depr}\propto n_{cl}^{\beta }$. The derived exponents of these fits are reported in Table 5 and indicate that a weak correlation exists between the separations and these quantities. The empirical trends traced by the medians closely match the result of the fits derived using all the separations in each of the two samples and show a tight power-law behavior for n_cl ≥ 10⁵ cm⁻³.

Fig. 12 Deprojected core separations against the clump surface (top panel) and volume (bottom panel) density. The distributions of the core separations are computed in bins of Σcl and ncl of width 0.15 and 0.2 dex, respectively. The peaks and spreads of these distributions, defined as the width at 68% peak intensity, are indicated by lines and area for the near (blue) and the far (green) samples, respectively. The medians of the distributions for these two sets are indicated as dashed and dotted black lines, respectively. The best fitting power-law relation computed over the entire set of separations and whose results are reported in Table 5 are indicated as light gray and dark gray lines for the near and far sample, respectively.

3.4 Comparison with the Jeans fragmentation

In this section, we compare the observed distribution of core separations with the expectations of the Jeans fragmentation theory for a non-magnetized, homogeneous, isothermal and self-gravitating clump. Under these conditions, there is a critical mass (and a corresponding length) above which the medium is not able to remain stable and therefore collapses (Jeans 1902). These critical values depend only on the density, n, and temperature, T (see Appendix B). Although these conditions refer to an idealized case which may not be representative of real star-forming clumps, the analytical prescription given by Jeans theory offers a simplified model against which the observed fragmentation in the ALMAGAL clumps can be compared. Specifically, we are looking for whether thermal fragmentation is sufficient to statistically describe the separations between cores or whether additional contributions, such as turbulent support, are needed. Since core masses are affected by higher uncertainties than core separations, we focused our analysis on the comparison between the latter and the thermal Jeans lengths, while we report the results for the masses in Appendix C.

3.4.1 The Jeans length and mass for ALMAGAL clumps

We computed the thermal Jeans length, ${\lambda }_J^{th}$, and mass, $M_J^{th}$, for the ALMAGAL clumps with Eq. (B.2) and Eq. (B.3), adopting as temperature the one of the dust, T_d, derived from Herschel/Hi-GAL data and reported in Table A.1 in Molinari et al. (2025). At the typical densities of ALMAGAL clumps (n>10⁴ cm⁻³), gas and dust are expected to be thermally coupled (Goldsmith 2001). This has been observationally confirmed, with relatively small uncertainties discussed below, in a sample of 1086 Hi-GAL clumps (Merello et al. 2019). In these clumps, on average dust temperatures were found to be consistent with temperatures measured from NH₃ emission (Merello et al. 2019), considered a reliable temperature tracer for the gas phase (Juvela et al. 2012; Guzmán et al. 2015). The values of ${\lambda }_J^{th}$ and $M_J^{th}$ for each ALMAGAL clump are reported in Table 3. The thermal Jeans length ranges from ∼3500 to ∼60 000 au(∼0.017−0.3 pc), with a median value of ∼14 500 au(∼0.07 pc), and Jean masses from ∼0.2 to 14 M_⊙, with a median equal to 0.89 M_⊙.

The uncertainties of these quantities are determined by those on the temperature and density, which depend on both their direct measurements and on the assumptions adopted for their estimates. The fit to the Herschel data typically provides small uncertainties for T_d, quantified as ∼2 K (Elia et al. 2017; Svoboda et al. 2019). This value translates into a different contribution for the young and evolved systems of the ALMAGAL sample, characterized by a typical temperature of ≈ 10 and ≈ 40 K, respectively. The uncertainties in ${\lambda }_J^{th}$ and $M_J^{th}$ vary from 2% to 10% and 6% to 30%, with the highest values for younger systems. Although the dust and gas temperature are found to be in agreement on average, their relationship presents a measurable scatter whose amplitude quantifies how well we can assess the thermal coupling between these two phases. Typical discrepancies between T_g and T_d are of the order ∼25% (see Fig. 2 in Merello et al. (2019)), where the distribution of log (T_g/T_d) has a width of ∼0.2 dex. These discrepancies introduce an additional random contribution of ∼13% and ∼38% to the uncertainties in ${\lambda }_J^{th}$ and $M_J^{th}$, respectively. In addition to the random scatter, T_d and T_g show a systematic offset in different stages of the clump evolution, which may lead to an underestimation or overestimation of ${\lambda }_J^{th}$ and $M_J^{th}$. Ammonia-derived gas temperatures are found to be systematically ∼20% higher than T_d in young systems such as prestellar clumps (Merello et al. 2019), where star formation activity is absent, while no significant offset is found in evolved systems. However, those evolved phases also have other temperature diagnostics as additional molecular tracers become available. The gas temperatures T_g estimated from these diagnostics are generally higher than T_d derived from the fitting of the spectral energy distribution of the Herschel data (Molinari et al. 2016b; Lin et al. 2022). Specifically, Molinari et al. (2016b) estimated gas temperatures from CH₃C₂H(2−1) emission finding that they are about ∼40% higher than dust ones. Similar discrepancies were also found by Lin et al. (2022), who instead analyzed the CH₃CN and CH₃OH emission lines in a few clumps with L/M ≥ 10 L_⊙/M_⊙, although they also report extreme cases associated with HII regions, where estimated T_g are as high as ∼130 K which is a few times higher than T_d.

The discrepancy between gas and dust temperature may be explained by different reasons. Although fitting the dust continuum provides an accurate determination of T_d, this emission may include contributions from the foreground or background. In addition, it is not always clear whether the molecular line emission probes the entire clump or only a fraction of it (Mininni et al. 2025), in particular in more evolved systems where temperature gradients are certainly set by the hosted protostars. However, the systematic difference between T_d and T_g indicates that the clump-averaged ${\lambda }_J^{th}$ and $M_J^{th}$ may be different from our estimates in different evolutionary stages. Nonetheless, the expected average increase of these quantities does not exceed ∼20% and ∼40% for ${\lambda }_J^{th}$ and $M_J^{th}$, respectively. At the same time, we consider an additional random contribution of ∼10% and ∼30% to the uncertainties in ${\lambda }_J^{th}$ and $M_J^{th}$ induced by the typical scatter of ∼20% between the two estimates T_g and T_d.

Unlike the case of temperature, the average volume density n_c is not directly measured, but is typically inferred from measurements of the clump mass and radius. Assuming an uncertainty of 50% in M_cl, 20% in the measured source size and 15% in the clump distance, the uncertainty in n_cl is about ∼60%, which results in ∼30% for ${\lambda }_J^{th}$ and $M_J^{th}$. These estimates are obtained by assuming a uniform distribution of the material, which may not be realistic due to the density stratification present in clumps (Beuther et al. 2002; Williams et al. 2005; Palau et al. 2014; Beuther et al. 2024). In these systems, the density profile has a typical power-law shape, i.e., $\rho ={\rho }_0{\left(\frac{R}{R_0}\right)}^{-p}$ with ρ₀ and R₀ the density at the reference radius R₀. The local density in the inner region of size R^′ is therefore enhanced by a factor ${\left(\frac{R^{\text{'}}}{R_0}\right)}^{-p}$, which depends on the value of p. Observations of highmass clumps indicate that p varies in the interval 1<p<2, with an average value of p ≈ 1.6 (Beuther et al. 2002; Williams et al. 2005; Palau et al. 2014; Beuther et al. 2024, see also the discussion in Beuther et al. 2025 for a summary of the recent measurements). The observed variation appears to be connected to the clump evolutionary state. In fact, Lin et al. (2022) reported an increase from p ≈ 1 to ∼1.5 in their sample, and Gieser et al. (2022) found slightly shallower profiles with p ≈ 1.6 in the younger regions, whereas they also measured p=2.0 ± 0.2 in the evolved massive clumps of the CORE surveys (Gieser et al. 2021; Beuther et al. 2024).

In principle, the density profile and the exponent p could be constrained for each clump by fitting the Hi-GAL observations. However, this would require a careful and detailed radiative transfer modeling to properly account for the temperature structure, which is beyond the scope of this work. Our aim here is to quantify the amplitude of the uncertainty that the assumed density structure introduces into our estimates for Jeans properties. In the case of ALMAGAL clumps, fragmentation is effectively taking place in a region defined by the cluster radius, $R_{cluster}^{ph}$ (see Sect. 3.1 and Table 2). Due to density stratification, the average density is higher in the central regions, which results in a decrease for ${\lambda }_J^{th}$ and $M_J^{th}$ by a factor $C^{\frac{p}{2}}$, where $C=\left(R_{cluster}^{ph}/R_c\right)$. The median value of C is ∼0.46, with quartiles of ∼0.36 and ∼0.59, which translate into an expected ${\lambda }_J^{t{h}^{\text{'}}}=0.{68}_{-0.1}^{+0.09}\times {\lambda }_J^{th}$, for the young systems with p ≈ 1, and ${\lambda }_J^{t{h}^{\text{'}}}=0.{46}_{-0.1}^{+0.13}\times {\lambda }_J^{th}$ for evolved systems where p ≈ 2.

In general, our assumptions in T_g and n_c introduce systematic shifts in the average estimates of ${\lambda }_J^{th}$ and $M_J^{th}$, that go in opposite directions. The amplitudes of these shifts are similar, and, on average, compensate; therefore, we did not apply any correction. However, since their contribution leads to an increase in the spread of the uncertainty, in our analysis we considered a typical random error on ${\lambda }_J^{th}$ of at least ∼50%.

3.4.2 ALMAGAL core separations compared to Jeans length

Figure 13 shows the comparison of the two clump-averaged estimators for the core separation, the median l^median and the mode l^mode (see Sect. 3.2.2), with the estimated thermal Jeans length, ${\lambda }_J^{th}$. In these plots, we mark the level of fragmentation and L/M of the clump with the symbol size and the color, respectively. A large fraction of ALMAGAL clumps have core separations that are consistent with the thermal Jeans length ${\lambda }_J^{th}$. The mode [median] separations range from ∼0.07[0.1] × to $\sim 6.5[11.0]\times {\lambda }_J^{th}$. The distributions of the ratio $l^{median}/{\lambda }_J^{th}$ and $l^{mode}/{\lambda }_J^{th}$ are centered on values close to 1 and cover the narrow intervals included in Fig. 13 with dark and light gray shaded areas representing the interquartile and the 10th to 90th percentile interval, respectively. Quantitatively, the median and the mode are centered around 0.91 × and $0.76\times {\lambda }_J^{th}$, respectively. The first and third quartiles are equal to 0.56 ×−1.60 × and $0.5\times -1.26\times {\lambda }_J^{th}$ for the median and the mode, respectively. The 10th-90th percentile ranges extend from 0.29 × to 2.10 × for the mode and 0.36 × to $2.53\times {\lambda }_J^{th}$ for the median. Both the clump-averaged median and mode reveal a similar behaviour with the thermal Jeans length.

Although on average the core separations are in agreement with the thermal Jeans length, we observe that the systems with different L/M are distributed in different regions of Fig. 13, indicating that the ratio $R=l^{depr}/{\lambda }_J^{th}$ varies with L/M and, therefore, changes during the clump evolution. The trend is shown in Fig. 14, where we present a panoramic view of how the ratio R varies with L/M for all separations of the ensemble. We found a decreasing trend which is the result of two effects: the decrease in the core separations l^depr discussed in Sect. 3.3.1, and a steady increase in ${\lambda }_J^{th}$ with L/M. We note that the increase in ${\lambda }_J^{th}$ is mainly driven by the general warming of the clump material by the radiation of newly born stars. Although it also depends on density n_cl, this has a marginal impact when averaging over bins of L/M ratio. In fact, our sample includes clumps with density spanning over the entire range of values in each of those bins (see Fig. 1). More generally, the selected ALMAGAL clumps do not show an increase in their density with evolution (see Sect. 2.1).

Given the uncertainties in ${\lambda }_J^{th}$, our result is that the observed separations are consistent with the thermal Jeans length computed with the average properties derived from Herschel observations for all clumps with the L/M ratio between 1 and 100 L_⊙/M_⊙. Figure 14 shows an average systematic decrease in R with L/M and that the core separations become significantly smaller than the thermal Jeans lengths in the evolved systems with L/M ≫ 50 L_⊙/M_⊙(L/M ≫ 100 L_⊙/M_⊙. in the case of near sample). This trend for L/M ≳ 1 L_⊙/M_⊙ is mainly driven by an increase in ${\lambda }_J^{th}$ and not by a decrease in the core separations. The observation of sub-Jeans separations in these evolved stages does not rule out gravity-dominated thermal fragmentation. In these phases, protostellar feedback could create local density enhancements that exceed the mean clump density used for our ${\lambda }_J^{th}$ calculation. Therefore, Jeans-like fragmentation of these locally denser regions would naturally produce cores with separations smaller than our clump averaged estimates.

It is interesting to look in detail at the distribution of R for young and poorly fragmented clumps with L/M<0.5− 1 L_⊙/M_⊙. The distributions have a median R ∼1.5, but the shape of the density contours suggests a bimodality, which we confirmed with the statistical tests described in Sect. 3.3.1. Although the two peaks for the far sample systems are connected to the bimodal distribution of l_i, they are also present in the near systems. One of the two peaks of the distributions is centered at R ∼1, indicating a population of cores which are separated by a thermal Jeans length, while the second peak is shifted at higher values, which are centered at about R ∼3.

Fig. 13 Comparison between the clump-averaged estimators for the core separations, specifically the median (top panel) and mode (bottom panel) of the distribution in each clump, and the computed Jeans length. The fragmentation level in the clump is indicated by the symbol sizes, with increasing size going from low (4 ≤ Ncores<8) and intermediate (8 ≤ Ncores<16), up to high (Ncores ≥ 16) core numbers. The dashed lines are drawn as a reference and correspond to the separation equal to ${\lambda }_J^{th}$ (thick line) and 2 × and $0.5\times {\lambda }_J^{th}$ (thin lines), respectively. The gray shaded boxes indicate the amplitude of the distributions of the ratio $l_i^{depr}/{\lambda }_J^{th}$ (see text), and correspond to the interquartile range (dark gray) and the 10th to 90th percentile range (light gray), respectively.

3.5 How are cores distributed in the fields?

In previous sections, we analyzed where cores are located in star-forming clumps and their relative separations to discuss how fragmentation takes place in potentially cluster-forming systems. How cores are distributed in the systems and the identification of potential different mass-dependent clustering of the cores offer additional information on the cluster formation process. To such aim, we analyzed two additional diagnostics, the Q parameter and the mass segregation ratio, Λ_MSR, which are typically used to characterize young proto-clusters (Allison et al. 2009; Cartwright & Whitworth 2004; Parker 2018; Dib & Basu 2018), and were also recently estimated for groups of cores (Dib & Henning 2019; Sanhueza et al. 2019; Sadaghiani et al. 2020; Avison et al. 2023; Morii et al. 2024; Xu et al. 2024).

3.5.1 The Q parameter

Statistical analysis provides several diagnostics capable to characterize the spatial distribution of young stellar groups. Cartwright & Whitworth (2004) introduced the Q parameter as a tool to quantify the degree of substructures in the distribution of the members of young stellar clusters and, if present, to identify cases where their radial density profiles are centrally concentrated (Cartwright & Whitworth 2004; Parker & Goodwin 2015; Jaffa et al. 2017). The Q parameter is defined as the ratio of the mean edge length of the cluster MST, m̄, and the correlation length of the cluster, s̄, both with normalization factors defined to remove any dependence on the number of members and the size of the cluster. In formula these are $$\overline{m}=\sum _{i=1}^{Nc-1}\frac{l_i}{\sqrt{N_{cores}\pi R_{cluster}^2}},$$(3) where l_i are the MST edges and R_cluster is the size of the cluster defined as in Sect. 3.1.1 and $$\overline{s}=\frac{{\overline{l}}_i}{R_{cluster}},$$(4) with l_i defined as the average separation between the cluster members. The ratio Q=m̄/s̄ allows us to distinguish between two specific spatial configurations, one where the members are radially distributed with a central concentration and the other where they have fractal type of subclustering. Cartwright & Whitworth (2004) showed that systems with members distributed in smooth and centrally concentrated patterns have Q > 0.8, with values that increase for higher concentrations represented by higher exponents p in the number density profile N(r) ∝ r^−p. In contrast, systems whose members are organized in a fractal type of subclustering, described with a fractal dimension D, have Q<0.8, and smaller values of Q correspond to systems with a higher level of subclustering (lower values of D).

Although this diagnostic is typically used in the statistical analysis of young star clusters (Parker 2014; Gregorio-Hetem et al. 2015; Dib et al. 2018), it has only recently been applied to the distribution of star-forming cores, finding, for example, that Q correlates with star formation activity (Dib & Henning 2019). Specifically, Q increases with the density of star formation rate in low-mass star-forming clouds (Dib & Henning 2019). Moreover, the comparison between the core population observed by the ASHES (Sanhueza et al. 2019; Morii et al. 2024) and ASSEMBLE surveys (Xu et al. 2024) shows that different spatial arrangements are found depending on the evolutionary stage of the clump. Most of the ASHES clumps (∼70%), that are young and with L/M ≲ 2 L_⊙/M_⊙, have Q ≲ 0.8, which characterizes subclustered systems. In contrast, higher values of Q are typically measured in the more evolved clumps of the ASSEMBLE survey that have L/M ≳ 10 L_⊙/M_⊙ (Xu et al. 2024). The presence of a weak correlation between the Q parameter and the ratio L/M has been proposed to reflect a time evolution effect, where the signature of a gradual transition from a subclustered pattern to a centrally condensed configuration is a consequence of gravitational collapse (Xu et al. 2024). We tested the statistical robustness of the correlation Q versus L/M with the ALMAGAL sample, which largely extends the sample size of the ASHES and ASSEMBLE surveys from which it was initially derived.

The Q parameter measured in the ALMAGAL clumps with N_cores ≥ 4 shows values ranging from 0.4 to 0.95, with a median equal to 0.65, and half of the systems having 0.6<Q<0.75. These values would suggest that the ALMAGAL sample is composed mainly of clumps where the cores are distributed in subclusters, if the classification criteria introduced by Cartwright & Whitworth (2004) were applied. However, recently Avison et al. (2023) questioned the diagnostic usefulness of the Q parameter when applied to a small number of sources (N_cores< 10). They simulated clusters composed of small numbers of members (N_cores=5 and 10) to characterize the observations of the TEMPO survey, where they found between 2 and 15 fragments (Avison et al. 2023). Most of the ALMAGAL clumps with these conditions (N_cores<10) would be classified as highly subclustered as ∼80% of them have Q ≤ 0.6. To quantitatively assess whether the Q parameter can be used as a useful diagnostic in the ALMAGAL fields that have a wider interval of the number of cores than TEMPO survey (N_cores from 2 to 49), we extended the analysis of Avison et al. (2023). This analysis is presented in Appendix D and is also applicable to the ASHES and ASSEMBLE fields, which host 8 to 37 and 13 to 37 fragments, respectively (Morii et al. 2024; Xu et al. 2024).

Our result is that the Q parameter retains a diagnostic power in probabilistic terms only for systems with N_cores ≥ 10. We defined two intervals of Q where it is most likely that the cores follow a centrally concentrated (Q ≳ 0.8 for high power-law exponent p) or subclustered distribution (Q ≲ 0.65, lower fractal dimension), with a value of Q=0.75 that approximately separates these two regimes. These results revise the initial impression on the ALMAGAL systems and indicate that both spatial distribution patterns are present in the 226 ALMAGAL systems with N_cores ≥ 10. In fact, adopting the tentative limit of Q=0.75, there are 155 fields (∼68%) with Q<0.75 and 71 fields (∼32%) with Q ≥ 0.75, and similar proportions are found if higher cuts on N_cores are applied. However, our simulations show that the statistical distributions of Q for the two classes are not strongly separated, but they instead overlap. This led us to assess the classification with sufficient reliability only for the 60 fields where Q<0.65 and for the 38 where Q>0.8. Although we reduce the sample size with this choice, we were still able to confirm that both patterns are observed in massive dense clumps.

We show in Fig. 15 the values of the Q parameter as a function of the clump L/M ratio for the 226 ALMAGAL fields with N_cores ≥ 10, and we included measurements from the ASHES pilot (Sanhueza et al. 2019) and the ASSEMBLE survey (Xu et al. 2024). The 39 fields with L/M ≤ 2 L_⊙/M_⊙ expand the sample provided by the ASHES survey for the early evolutionary phases (Sanhueza et al. 2019; Morii et al. 2024), and indicate that these early phases have cores that are mainly distributed in fractal-type subclusters, with 13 fields showing Q<0.65. Fields with a clear signature of a concentrated distribution are rare in these young clumps, and only four of them can be robustly classified as centrally concentrated (Q>0.8). Although the occurrence of systems with a more pronounced centrally concentrated distributions (i.e., Q>0.8) is certainly higher in the later stages of their evolution, mild or high subclustered distributions of cores are still present at later stages, as Q<0.65 is found in 27 systems out of 124 with L/M>10 L_⊙/M_⊙.

We do not find a significant relationship between Q and L/M, on the contrary, we observe a large scatter of Q for any interval of L/M, which is not connected to the level of fragmentation or to other clump-averaged properties, for example the clump surface density. It is unclear what the reasons are behind this large scatter, which may not be only caused by statistics but also by additional factors not considered in our analysis, as, for example, the surrounding local environment. A weak dependence of Q with the L/M ratio may be identified, but the Pearson and Spearman correlation coefficients equal to ρ_p=0.18 and ρ_s = 0.185 with p-values of 0.0067 and 0.0053 are significantly lower than the values reported by Xu et al. (2024): ρ_p=0.60 and ρ_s= 0.57. This indicates that if a transition between spatial distribution patterns exists in massive clumps, it may require additional factors than simply the internal clump evolution.

Fig. 14 Left panel: comparison of the ratio R of the deprojected separations to ${\lambda }_J^{th}$ and the clump L/M ratio for the near (top panel) and far sample (bottom panel), respectively. The black dashed line shows the median of the distributions of R computed in bins of L/M equally spaced in logarithm of width 0.4 dex. The thick colored lines and the colored area trace the mode and the spread (defined as the width at 68% of the peak) of these distributions, adopting the blue and green colors for the near and far sample, respectively. The numbers at the top give the size of the subsamples adopted to compute the distribution in each bin. The gray shaded area shows the interval 0.5 ≤ R ≤ 1.5, corresponding to an uncertainty on the ${\lambda }_J^{th}$ of a factor of 2. The colored contours trace, from the outside inward, fractions of the sample lying outside that rise from 15% to 95% with a step of 10% for the near (in blue) and the far (in green) subsample, respectively. Right panel: normalized density histogram and estimated probability density distribution estimated with the KDE method (colored continuous lines) of the ratio $R=l_i/{\lambda }_J^{th}$ measured in the three subsamples selected by different intervals of L/M tracing clumps in different evolutionary phases: young (L/M ≤ 1 L⊙/M⊙ in orange), intermediate (1 ≤ L/M ≤ 10 L⊙/M⊙ in green), and evolved (L/M ≥ 10 L⊙/M⊙ in purple).

Fig. 15 Relation between the Q parameter and the L/M ratio of ALMAGAL clumps which host at least ten cores. The two dashed lines values Q=0.65 and Q=0.8 defining the intervals for which it is most likely that the hosted cores follow a fractal type of subclustering (Q ≤ 0.65) and radial distribution with a central concentration (Q ≥ 0.8), as indicated by the Monte Carlo simulations in Appendix D. The dotted line points to the Q=0.75, which separates the two regimes with a lower statistical significance. Estimates from the ASHES (dark gray; Sanhueza et al. 2019) and ASSEMBLE surveys (light gray; Xu et al. 2024) are included, and are shown as star markers.

3.5.2 Mass segregation in the ALMAGAL fields

A further clue on the formation and early evolution process can be extracted by also taking into account the mass of the cores. Observations of evolved star clusters reveal that they are often mass segregated, with massive objects that are more concentrated towards the cluster centers than the low-mass ones (Hillenbrand & Hartmann 1998). This observational signature is often interpreted as a consequence of the dynamical evolution caused by close encounters between pairs of cluster members that lead the system to a global state of energy equipartition. As a result, the less massive cluster components rise to higher orbits around the gravitational center, whereas the more massive cluster components move more slowly, falling into lower orbits. However, young clusters may have been born with a primordial MS set by the formation process itself. In fact, primordial MS with respect to the center of the gravitational potential is a natural outcome of the competitive accretion model for high-mass stars (Bonnell & Bate 2006).

Evidence of MS has already been found in observations of young cluster systems in both high-mass and intermediate- to low-mass environments (Chen et al. 2007; Gennaro et al. 2011; Parker 2018; Plunkett et al. 2018; Dib et al. 2018; Dib & Henning 2019; Sadaghiani et al. 2020; Nony et al. 2021). This indicates that MS is not produced merely by long-term dynamical evolution, but it can be already established in the early phases. However, it remains unclear whether it is primordial, i.e., set in the initial formation stages, or the result of a rapid, early dynamical evolution of the system of cores. For instance, weak or no MS signature is found in the 39 young clumps of the ASHES survey (Sanhueza et al. 2019; Morii et al. 2023). In contrast, supporting the idea that it develops over time, clear segregation is observed in the more evolved systems, with 8 of 11 clumps in the ASSEMBLE survey presenting clear signatures of MS (Xu et al. 2024).

We searched for MS signature in ALMAGAL systems using the mass segregation ratio Λ_MSR (Allison et al. 2009; Parker 2018), which compares the total length of the MST of the N_MST most massive cores with the typical one of the entire protocluster. The value of Λ_MSR for N_MST cores is given by $${\mathrm{\Lambda }}_{MSR}\left(N_{MST}\right)=\frac{⟨l_{random}⟩}{l_M}\pm \frac{{\sigma }_{l,random}}{l_M},$$(5) where 〈 l_random〉 is the average of the total length of the MST derived by randomly selecting a total of N_MST cores from the cluster, while l_M is the total length of the MST of the N_MST most massive cores. We computed the average 〈 l_random〉 and the dispersion σ_{l, random} by randomly selecting the cores 1000 times as in Morii et al. (2023) and Xu et al. (2024). If MS is present, then the N_MST most massive cores are more concentrated, and the total length of their MST would be shorter than the one built from randomly selected cluster members, which will give Λ_MSR ≥ 1. In contrast, Λ_MSR would be ≈ 1 if the cores had a similar spatial distribution regardless of their mass.

We adopted the core masses estimated from Coletta et al. (2025) and computed Λ_MSR(N_MST) in the 296 ALMAGAL fields hosting at least eight cores to limit the uncertainties due to their small number. This threshold corresponds to the minimum number of detections in the ASHES field for which the ratio Λ_MSR was also estimated (Morii et al. 2024). If we do not take into account any uncertainty in the core mass, we find that a signature of mass segregation is revealed in 171 fields, which have a MSR profile that differs from unity for more than the random uncertainty σ_{l, random}/l_M (it exceeds unity in 131 fields). We adopted a MonteCarlo approach to take into account the uncertainty on the core masses. For each system, we generated 500 unique realizations where every core mass was drawn from a Gaussian distribution centered on its value, M_core, with a standard deviation of 20%. We computed the ${\Lambda }_{MSR}^j\left(N_{MST}\right)$ with the relative error for each of these 500 synthetic realizations, fully mapping the range of possible outcomes given the assumed mass uncertainties. The fraction of realizations where ${\Lambda }_{MSR}^j\left(N_{MST}\right)$ exceeds unity over the dispersion σ_{l, random} defines the probability that mass segregation is present in the system. We classified a system as robustly mass segregated if this fraction is higher than 70%, a condition found in 110 fields (∼37% of the total). This considerable fraction suggests that cores are already segregated in mass prior to dispersal of the gas present in the clump, and raises the question of whether segregation is set at the initial formation stage or is a consequence of the early evolution, which requires a comparative analysis between systems at different evolutionary stages.

Given their large numbers, we do not show all MSR profiles. Instead, we characterized each profile with metrics indicating average value, peak, and extent, reporting those quantities in Table 6. We included the mass segregation integral (MSI), $I_{\mathrm{\Lambda }}^{MSR}$ introduced by Xu et al. (2024) to quantify the fractional area of the MRS profile that exceeds unity. This quantity allows the comparison of the MSR profiles of different systems and is defined as $I_{\mathrm{\Lambda }}^{MSR}=\frac{\sum _{i=2}^{N_{MST}^{max}}{\mathrm{\Lambda }}_{MSR,i}}{N_{MST}^{max}-1}$, where Λ_{MSR, i} is the value of the profile Λ_MSR(N_i) for the most massive cores N_i. Moreover, we determined the maximum of the profile, ${\mathrm{\Lambda }}_{MSR}^{max}(N_{MST}^{max})$, and the largest number of cores, $N_{MST}^{crit}$, where it still remains above unity. These two quantities extend the information provided by the MSI, as they indicate how much more concentrated the most massive cores are than any other randomly selected group and the largest number of cores where the spatial distribution is found different.

Figure 16 shows the distribution of clumps with N_cores, ≥, 8 as a function of the ratio L/M for systems with and without MS signature, indicating how frequently MS is found during the entire evolutionary sequence of massive clumps. We divided the sample of 296 clumps into three different groups based on the L/M ratio: L/M<2 L_⊙/M_⊙, 2<L/M ≤ 10 L_⊙/M_⊙, and L/M>10 L_⊙/M_⊙. The first and latter intervals are also sampled by the targets of the ASHES and ASSEMBLE surveys, allowing us to extend their results with larger sample sizes. These groups consist of 59, 83, and 154 clumps and have an average number of cores equal to 12, 15 and 17, respectively. The 110 systems where mass segregated is robustly confirmed divide in 10, 25, and 75 systems in each of these three groups. All of those systems have a different MSI distribution, and we rejected the null hypothesis that they derive from the same underlying one with the Kolmogorov–Smirnov and Mann-Whitney U tests with confidence levels set by p_KS<0.029 in the two cases with L/M<2 L_⊙/M_⊙ and 2<L/M<10 L_⊙/M_⊙, and with a more robust assessment of p_KS, p_U ≲ 10⁻⁵ in the more numerous group with L/M>10 L_⊙/M_⊙.

The occurrence of MS in the young clumps of the ALMAGAL sample (∼17% with 10 out of 59) is similar (and compatible within the uncertainty of the counting statistic) to the one measured in the ASHES sample (13% with 5 out of 39). These systems have a ${\mathrm{\Lambda }}_{MST}^{max}$ between ∼2.5 and ∼7 with an average value of ∼5, indicating that when MS is present, even in these early phases, it usually has a moderate to high impact on the distribution of the cores. One peculiar system shows a very high value of ${\mathrm{\Lambda }}_{MST}^{max}\sim 20$, but visual inspection of the system indicates that it is closely associated with multiple mid-infrared bright sources and, thus potentially affected by the strong feedback of nearby massive young stellar objects. Since there is only one case, we consider it an outlier. The fraction of clumps with MS steadily increases in more advanced stages of evolution. It first rises to 30% (25 out of 83) and then reaches ∼49% (75 out of 154) for 2<L/M<10 L_⊙/M_⊙ and L/M>10 L_⊙/M_⊙, respectively. An indication of this effect was already suggested by Xu et al. (2024) who found that 8 out of 11 clumps with L/M>10 L_⊙/M_⊙ are mass segregated, but the large sample size provided by the ALMAGAL large program allows us to robustly and clearly assess this trend, also measuring the increase for intermediate evolutionary stages.

Furthermore, we also found that its signature affects a larger number of cores and becomes stronger in more evolved systems. Figure 17 shows the distributions of the maximum value of the mass segregation ratio, ${\mathrm{\Lambda }}_{MST}^{max}$, and of $N_{MST}^{crit}$ in the three aforementioned samples. The intensity of segregation is always greater than ∼2.5, with most systems having Λ_MSR>3, which is considered a high level of MS (Pouteau et al. 2023). The median value of ${\mathrm{\Lambda }}_{MST}^{max}$ is ∼5 in less evolved clumps, and increases to ∼6.4 and ∼6.2 for the intermediate and more evolved subsamples, respectively. This suggests that the cores with higher mass become more concentrated as the clump evolves.

Segregation generally affects only two or three most massive cores of the clusters, as ∼80% of the system where segregation is present have $N_{MST}^{crit}\le 3$, with 54 and 34 systems where $N_{MST}^{crit}=2$ and 3, respectively. While the evolution leads to an average increase in the number of cores (Elia et al. 2026), at the same time we found that the number of cores showing a tighter spatial distribution also increases and may involve up to the five most massive cores in the system (only two systems show $N_{MST}^{crit}\ge 6$).

In conclusion, although MS may be present in the earliest stages, it is about three times more common in evolved clumps (L/M>10 L_⊙/M_⊙). Segregation is more pronounced in systems in more evolved phases, when cores with masses of a few M_⊙ in ∼1000−2000 au are commonly found in the system (see Fig. 18 in Coletta et al. 2025). The gravitational concentration of the entire clump or the gas accretion toward its center has been suggested to be responsible for the observed increase in MS (Xu et al. 2024). We do not find any connection between the signature of MS and the clump surface density Σ_cl or the cluster radius R_c, as it would be expected if gravitational collapse was present. However, it is important to note that Herschel derived surface densities might not be sufficient to highlight the increase in density of the clump central region because of their poor resolution. Further analysis of the dynamics and the gas distribution in clumps with high signatures of MS is necessary using molecular line tracers at high angular resolution to verify whether clump material is efficiently transferred toward central regions and the most massive cores.

Fig. 16 Distribution of the ALMAGAL clumps hosting at least 10 cores in different bins of the ratio L/M. The colored histograms show the systems were a signature of MS is found, with different colors referring to the different intervals of L/M ratio, i.e., orange L/M ≤ 2 L⊙/M⊙, green 2<L/M ≤ 10 L⊙/M⊙ and purple L/M>10 L⊙/M⊙. The solid black line, which refers to the right y-axis scale, indicate the relative fraction of systems with MS with respect to the observed ones in bins of L/M, while the dot-dashed line shows the average occurrence fraction of MS systems in the three intervals of L/M.

Fig. 17 Boxplots showing the distribution of maximum value of the mass segregation ratio, ${\mathrm{\Lambda }}_{MSR}^{max}$, divided for the different values $N_{MST}^{crit}$ that identify the interval where the MSR profile is statistically higher than unity, in the 110 ALMAGAL clumps where signatures of MS are found. Those clumps are divided into three groups according to the L/M ratio composed of 59, 83, and 154 systems. The rectangular area of the boxplots shows the interquartile intervals, while the whiskers represent the 10th to 90th percentile intervals; the diamond symbols indicate the cases outside this interval.

4 Discussion

ALMAGAL observations reveal that fragmentation in high-mass clumps is widely diversified, with star-forming cores spatially arranged across a broad spectrum of configurations, from sparse, mildly concentrated, up to highly clustered. This observed heterogeneity, which mirrors the wide scatter found for core masses and fragmentation degrees discussed by Coletta et al. (2025) and Elia et al. (2026), does not preclude a deeper analysis. Rather, systematic analysis of spatial distributions and core-to-core separations provides a powerful diagnostic tool for the physical mechanisms that govern where and how fragmentation takes place. Our data show that most of the small-scale (∼1400 au) cores are confined within the inner clump region (≤ 0.35 pc). Additional cores forming at larger distances are rare or, if present, are likely low-mass objects undetectable in ALMAGAL fields, which would be compatible with their formation in low-density environment. This result might favor an in situ formation scenario for clusters (Dib et al. 2010). However, it does not contradict the more dynamic conveyor belt-type models (Krumholz & McKee 2020; Vázquez-Semadeni et al. 2019) that involve anisotropic rather than monolithic contraction, with material being fed into the clump along large-scale flows. The Global Hierarchical Collapse (GHC) model, an example of these dynamic models, predicts gravity-driven accretion along extended filaments to supply material to cores and hubs (Vázquez-Semadeni et al. 2019, 2025). GHC simulations show that low-mass cores may form along these peripheral filaments (Gómez & Vázquez-Semadeni 2014; Vázquez-Semadeni et al. 2017), a scenario consistent with observed cores found along filaments connected to high-density clumps (Myers 2009; Kirk et al. 2013; Peretto et al. 2014; Hacar et al. 2017; Chen et al. 2019; Álvarez-Gutiérrez et al. 2024). The ALMAGAL cores found on the outskirts of the clumps may be formed in similar flows whose underlying filaments remain undetected because of the limited sensitivity and coverage. In this context, the presence of aligned cores in several ALMAGAL systems can be an additional footprint of the presence of these flows. In fact, numerical simulations show that the cluster structure retains indications of the native distribution and organization in structures of the star-forming material (Vázquez-Semadeni et al. 2017).

Two significant results are found from the analysis of the separations between cores in the ALMAGAL high-mass clumps. The separation is on average of the order, or smaller, of the pure thermal Jeans length over almost the entire early evolution of the clump sampled by the ALMAGAL survey. This evidence confirms and extends the results of previous studies based on a limited number of clumps (Palau et al. 2015, 2018; Beuther et al. 2018; Liu et al. 2020; Zhang et al. 2021; Morii et al. 2024; Ishihara et al. 2024). Separations close to the thermal Jeans length are predicted in the case that the fragmentation process is driven mainly by gravitational contraction (Palau et al. 2014; Vázquez-Semadeni et al. 2025). Nonetheless, we found indications of systems in which the separations between the cores are significantly larger than the thermal Jeans length, ${\lambda }_J^{th}$, by an average factor of three. One possible explanation is that there is a contribution due to turbulence, which provides additional support and increases the Jeans length (Mac Low & Klessen 2004). The turbulent Jeans length can be computed by replacing in Eq. (B.2) the sound speed c_s with the gas velocity dispersion σ_v given by ${\sigma }_v=\sqrt{c_s^2+{\sigma }_{nth}^2}$, where ${\sigma }_{nth}^2$ are the non-thermal motions present in the gas. The typical properties of the ALMAGAL clumps in the earliest phases have a median value of c_s ≈ 0.2 km s⁻¹; therefore, to explain the observed ratio R ∼3 it would be sufficient to have a moderate non-thermal velocity dispersion σ_nt ≈ 0.6 km s⁻¹, which corresponds to Mach number ℳ ≈ 5, evaluated as $\mathcal{M}=\sqrt{3}{\sigma }_{nth}/c_s$. To compare it with the average conditions present in the entire hosting molecular cloud, we matched the position of the ALMAGAL clumps with the clouds identified in the SEDIGISM survey using the ¹³CO(2−1) emission (Duarte-Cabral et al. 2021). We found that the typical σ_v of the matched clouds is ∼1.8 km s⁻¹, with a distribution with interquartile interval from ∼1.3 to ∼2.6 km s⁻¹. These values indicate a remarkable contribution of non-thermal velocity dispersion at the cloud scale, which would produce significantly larger separations between cores than the observed ones.

Another possibility is that different fragmentation regimes, each with its own characteristic length, are present in a fraction of ALMAGAL clumps, which would also explain the asymmetries found in the distribution of separations P(l_i). These regimes may be due to different contributions of the magnetic field, local turbulence, and density distribution (Fontani et al. 2018), be the result of the assembly of several subclusters (Morii et al. 2023), or the fragmentation in the filaments connected to the accretion flow expected in the GHC model (Vázquez-Semadeni et al. 2019). A more robust interpretation of our observations requires a one-to-one comparison with dedicated simulations such as these developed within the Rosetta Stone project which model the formation and temporal evolution of massive clumps (Lebreuilly et al. 2025), produce post-processed maps at multiple wavelengths (Tung et al. 2025), and create mock observations that closely resemble ALMA data with high fidelity (Nucara et al. 2025).

In addition to the comparison with the Jeans fragmentation theory, an important result emerging from the ALMAGAL clusters of cores is the progressive statistical decrease exhibited by the separations with the L/M ratio, which can be interpreted as evidence that cores are more tightly packed in more evolved systems compared to the younger ones. This evidence has been advocated in support of the GHC model for high-mass and cluster formation as it naturally derives from the global contraction of the systems (Traficante et al. 2023; Xu et al. 2024; Vázquez-Semadeni et al. 2025). Our analysis shows that typical separations decrease by about a factor of three, with a major decrease from ∼20 000 down to ∼10 000 au taking place at L/M ≈ 1 L_⊙/M_⊙. Above this value, our data show only marginal variations, and the observed trend does not clearly reflect a decrease in core separations due to the ongoing global collapse of the clump. If a dynamical shrinking of the cluster is occurring, it is not the dominant effect traced by the median of the population of separation measured in ALMAGAL fields. The value L/M ≈ 1 L_⊙/M_⊙ is reached when the luminosity of already formed protostars becomes relevant for the entire clump (Molinari et al. 2008), and corresponds to an increase in the degree of fragmentation observed in these systems. The relation of the degree of fragmentation with evolution is discussed by Elia et al. (2026) who reported a steady increase on average in the number of cores over the entire interval 0.1 ≤ L/M ≤ 100 L_⊙/M_⊙. The increase in fragmentation is also presented in Fig. 11, which also shows the average variation of the core separations. We note that the observed decrease in these separations follows the transition from clumps characterized mainly by a limited number of cores (N_cores ≤ 8) to more crowded systems.

We attempted to draw an overall scenario from ALMAGAL observations. Our results indicate that fragmentation occurs in multiple substructures of the clump, which are generally widely separated by ∼15 000−30 000 au(∼0.07−0.14 pc), measured in systems where fragmentation is at the onset. These substructures host a limited number of already collapsed cores, but also fragments that did not yet have enough time to grow sufficiently to reach a high contrast over the parent structure to be detectable, as proposed by Louvet et al. (2024) to justify the observed small fragmentation present in young systems. As the clump evolves, additional cores collapse and become detectable in the ALMAGAL fields with a sufficiently high significance. These newly collapsed cores form in the midst of the original population, but a minor fraction of these objects also form in the outskirts of the clump, even though less material is available. Those effects produce a decrease in the observed average separation while the additional cores form by collapsing material.

Nonetheless, it is notable that systems with cores distributed following a smooth and radially concentrated profile start to be observed at the same L/M value where the typical core separations decrease. In addition to the increase in fragmentation level, the emergence of those systems may be the signature of a gradual merging and dispersion of the initial, fractal distribution in subclusters, that characterizes the young ALMAGAL clumps. This outcome is consistent with predictions from dynamical-type models where conveyor belt-like processes are active (Longmore et al. 2014; Krumholz & McKee 2020; Vázquez-Semadeni et al. 2025). According to these predictions, time evolution causes an increase in the degree of concentration of cores as the system shrinks and dynamically relaxes, as proposed by Xu et al. (2024). However, while our results indicate that smooth, centrally concentrated distributions occur more frequently at later stages, we do not observe a tight correlation with evolution that would robustly support this scenario over the entire sample. Instead, we still observe that subclustering signatures are present in a considerable fraction of systems with high values of the L/M ratio. As this pattern is often observed also in young embedded stellar clusters (Kuhn et al. 2014), this suggests the presence of mechanisms that work against the merging of these substructures or that are able to slow the process down for long timescales. The scatter of Q values and, more generally, the observed heterogeneity of the fragmentation in the ALMAGAL sample questions whether a unique path of evolution, like the one described above, can be traced for all the massive dense clumps observed. Additional factors that are not taken into account in this analysis may influence the evolution of the spatial distribution of cores and lead to the scatter of Q values for each given interval of the L/M ratio. Low values of the Q parameter are also measured when the main system increases its radial concentration, but additional substructures form in the outskirts of the field and become statistically significant. Moreover, it has been found that the intensity of the internal mechanical feedback affects the final structure and distribution of stars in clusters. Numerical simulations including protostellar jets have typically higher values of Q over time than those without any feedback or with only the radiative one (Verliat et al. 2022). In addition to the evolutionary effect, the variations in measured Q may reflect differences in the initial conditions with which the region had initially assembled and fragmented, as pointed out by Dib & Henning (2019) from the analysis of nearby low-mass star-forming regions. As the ALMAGAL fields cover different regions in the Galaxy, the clumps may be formed with different physical conditions imprinted by the large scale environment. In summary, our results are compatible with a dynamical evolutionary scenario, which defines the evolutionary path of a subset of massive clumps, and we cannot exclude the possibility that dynamical effects such as gravitational collapse, system shrinking, and dynamical relaxation are active in massive clumps.

5 Conclusion

We analyzed the spatial distribution of the 5728 dense cores extracted by Coletta et al. (2025) in the 514 ALMAGAL clumps with at least four detections. We provided a quantitative and statistically robust characterization of the fragmentation process in massive and potentially cluster-forming clumps, determining where the cores are located, the typical length of their reciprocal separations, and how they are distributed. Our measurements compose a comprehensive set of observational constraints for the fragmentation process against which to compare the predictions of different theories of cluster formation. We focused our analysis on how the properties of the statistical distribution change during the evolution of the systems before gas dispersal.

Our results are as follows:

• The ALMAGAL clumps show a wide variety of fragmentation patterns, with cores distributed over regions that extend from ≲ 0.1 up to ∼0.8 pc, which closely match the inner size of the hosting clump as determined from Herschel observations. Although cores are also observed in the outskirts of the clump, they are either rare or are low mass (and thus often undetected), as expected from formation in a low-density environment. We introduced a tentative classification based on the components and elongation of the convex hull polygon defined from the core positions, finding that several configuration patterns are present in the sample. We classified them into aligned, sparse, circular, elongated, and subclustered arrangements, with the most common pattern being moderately elongated, with aspect ratio ∼2−3, although we note that highly elongated systems are also present;

• The projected separations between the cores measured with the MST method show a broad distribution with a median value of ∼6600 au and ∼7400 au for clumps that are located at heliocentric distances 2<d_cl<4.5 kpc and 4.5<d_cl< 7.5 kpc, respectively. Smaller median separations are measured in systems with d_cl ≤ 3.7 kpc, suggesting that there are different characteristic separations in the inner regions of the clumps, which are expected to be denser. Evidence that a single characteristic length for a clump may not completely describe the observed fragmentation is also found via a comparison of the medians and modes of the distribution of the core separations, which differ from each other in 121 systems (∼24% of the analyzed sample). These differences are due to a local subclustering or to the presence of a more widely spaced population of cores in the outskirts of the system that followed a different fragmentation pattern;

• We identified a statistically significant shift in the distribution of core separations with the L/M ratio, which can be interpreted as a variation during clump evolution. The observed shift varies significantly in the early phases (L/M ≲ 2 L_⊙/M_⊙) when it drops from ∼20 000 au to ∼10 000 au for 2 ≤ L/M ≤ 10 L_⊙/M_⊙, and displays smaller changes for more evolved phases dropping to ∼7000 au for systems with L/M ≥ 10 L_⊙/M_⊙. Although the observed trend can be statistically described by a steady decrease in these two intervals, the larger variation appears to be produced by an increase in the degree of fragmentation in more evolved systems, as more cores are detected;

• We found a good agreement between the core separations and the thermal Jeans lengths ${\lambda }_J^{th}$ computed from the clump-averaged physical parameters, with a median ranging between ∼0.5 and ∼2 times ${\lambda }_J^{th}$. This agreement remains when the separations are divided according to the clump L/M. Nevertheless, we identified a small population of cores whose separations are centered at $\sim 3\times {\lambda }_J^{th}$, which are still smaller than the prediction of turbulent fragmentation of the hosting molecular cloud. These cores are mostly hosted in clumps with low levels of fragmentation and in early phases of evolution, with L/M ≲ 0.5−1 L_⊙/M_⊙. However, our results are compatible with thermal Jeans fragmentation being active in massive clumps. Although the masses are affected by large uncertainties and systematics, we note that the typical mass of cores is compatible with the thermal Jeans mass only in young and pristine clumps, where we observe that the core mass distribution is centered at $M_{\text{core}}\sim 0.6{\mathrm{M}}_J^{th}$;

• We did not find robust evidence of a progressive transition of the structure of the spatial distribution of the cores with evolution. The analysis using the Q parameter indicates that the cores are often distributed in a fractal type of subclustering in the early stages characterized by L/M ≤ 1 L_⊙/M_⊙. Although we found several clusters that are most likely distributed with radially concentrated profiles in the evolved clumps with L/M ≥ 10 L_⊙/M_⊙, subclustering is still often detected in those phases. The existence of these systems suggests that, in addition to the clump average properties, there are other factors that play a role in shaping the structure of the cluster over the evolution;

• Signatures of mass segregation are found in 110 fields out of the 296 ALMAGAL clumps hosting more than eight cores. The occurrence of mass segregation is lower (∼17%) in young systems, but steadily increases with the evolutionary stage to ∼30% in intermediate-evolved clumps with 2 ≲ L/M<10 L_⊙/M_⊙ and up to ∼48% in the more evolved systems in the ALMAGAL sample with L/M ≥ 10 L_⊙/M_⊙.

Data availability

Tables 2, 3, and 6 are available at the CDS via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/707/A221.

Acknowledgements

We are very grateful to the anonymous referee for constructive comments that helped clarify a number of points in the paper. The Teams at INAF-IAPS, Heidelberg, Paris-Saclay and University of Bologna acknowledge financial support from the European Research Council via the ERC Synergy Grant “ECOGAL” (project ID 855130). R.S.K. also acknowledges financial support from the Heidelberg Cluster of Excellence (EXC 2181 – 390900948) “STRUCTURES”, funded by the German Excellence Strategy, and from the German Ministry for Economic Affairs and Climate Action in project “MAINN” (funding ID 50 OO2206), and is also grateful to the Harvard Radcliffe Institute for Advanced Studies and Harvard-Smithsonian Center for Astrophysics for their hospitality and support during his sabbatical. A.S.-M. acknowledges support from the RyC2021-032892-I grant funded by MCIN/AEI/10.13039/501100011033 and by the European Union ‘Next Generation EU’/PRTR, as well as the program Unidad de Excelencia María de Maezto CEZ2020-001058-M, and support from PID2023-144675NB-I00 (MCI-AEI-FEDER, UE). Part of this research was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration (80NM0018D0004). R.K. acknowledges financial support via the Heisenberg Research Grant funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under grant no. KU 2849/9, project no. 445783058. P.S. was partially supported by a Grant-in-Aid for Scientific Research (KAKENHI Number JP22H01271 and JP23H01221) of JSPS. P.S. was supported by Yoshinori Ohsumi Fund (Yoshinori Ohsumi Award for Fundamental Research). S.W. has been supported by the Collaborative Research Centre 1601 (SFB 1601) funded by the Deutsche Forschungsgemeinschaft (DFG). L.B. gratefully acknowledges support from the ANID BASAL project FB210003. Y.-W. T. is supported by NSTC grants 112-2112-M-001-066- and NSTC 114-2112-M-001-042. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2019.1.00195.L. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

References

Appendix A Probability distribution of the projected separation

The measurable separations between cores are the projection on the plane of the sky of the real 3D distances between pair members. To derive these 3D distances, we can only rely on statistical analysis that, once applied to the distribution of projected separations, allows us to infer the 3D distance between cores. This approach is based on the estimation of the probability distribution density of the projected distances L_p, p(L_p) of a pair of objects that are randomly placed in the 3D space. As a first step, we consider only the effect of the projection on these two objects with fixed separation L, without considering any distribution of their separation. Assuming two points in the 3D space with a fixed separation L, the projected separation observed along a generic direction is equal to L_p=L sin i, where i is the inclination angle of the vector connecting the two points and the line of sight (see Fig. A.1). For a randomly oriented connecting vector between the two points, the probability of observing a separation L_p is equal to the area of the annulus dS=2 π L_p L di normalized to the surface of a sphere of radius L, in other words: $$dP\left(L_p\right)=\frac{dS}{4\pi L^2}=\frac{L_p}{L}\frac{di}{2}=\frac{\sin idi}{2}.$$(A.1) Differentiating L_p=L sin i and doing some algebra yields the following: $$di=\frac{d{L}_p}{L\cos i}=\frac{1}{L}\frac{d{L}_p}{\sqrt{1-{\sin }^{2}i}}.$$(A.2) Replacing Eq. A.2 in Eq. A.1, we obtain the probability density distribution function p(L_p), $$p\left(L_p\right)=\frac{dP}{d{L}_p}=\frac{L_p}{L\sqrt{L^2-L_p^2}},$$(A.3) which is normalized to unity as ${\int }_0^{L}p\left(L_p\right)d{L}_p=1$.

The expected value 〈 L_p〉 of the projected separation L_p is derived from the probability density distribution function p(L_p) as $$\langle L_p\rangle ={\int }_0^L{L}_{p}p\left(L_p\right)d{L}_p=\frac{\pi }{4}L$$(A.4) providing the average factor that link the intrinsic 3D separation to the projected one in the case of random orientation for a pair of objects. A similar estimate has been provided by Ishihara et al. (2024), against which this computation provides a statistical distribution. More generally, in the case that the 3D distance between the two points is not fixed, but has a probability distribution function P^′(L), then the projected separations will be distributed as $$p\left(L_p\right)=\int \frac{1}{L}\frac{P^{\text{'}}\left(L\right)L_{p}dL}{\sqrt{L^2-L_p^2}},$$(A.5) whose formulation has the advantage of describing how to statistically deproject the projected probability distribution function P^′(L) that can be derived from the observed separations. Equation A.5 can be inverted with numerical methods, or by assuming a shape for P^′(L). An important aspect of Eq. A.5 is that it is extremely unlikely for a unimodal 3D distribution P^′(L) to produce a multimodal 2D projected distribution p(L_p).

Fig. A.1 Scheme of projection of a pair of randomly oriented objects in 3D space

Appendix B Thermal Jeans analysis

The linear analysis of the stability of a non-magnetized, homogeneous, isothermal, and self-gravitating clump shows that there is a critical length above which density perturbations are always unstable against gravitational contraction (Jeans 1902; Mac Low & Klessen 2004). This critical length, and the mass contained in a sphere of this diameter, are usually referred to as the thermal Jeans length, ${\lambda }_J^{th}$, and mass, $M_J^{th}$, which, in the absence of turbulent support, depend on the density, ρ_cl, and the internal sound speed, c_s as $${\lambda }_J^{th}={\left(\frac{\pi c_s^2}{G{\rho }_{cl}}\right)}^{\frac{1}{2}}\text{and}{M}_J^{th}=\frac{{\pi }^{\frac{5}{2}}}{6}\frac{c_s^3}{\sqrt{G^3{\rho }_{cl}}},$$(B.1) where G is the gravitational constant.

Expressing the sound speed in terms of the temperature of the clump and assuming a mean molecular weight μ=2.37, which corresponds to a molecular gas with solar composition (Kauffmann et al. 2008; see also our Sect. 3.3.2), these relations become $${\lambda }_J^{th}\approx 4300\mathrm{a}\mathrm{u}(\frac{\mathrm{T}}{10\mathrm{K}}{)}^{\frac{1}{2}}(\frac{\mathrm{n}}{{10}^6{\mathrm{c}\mathrm{m}}^{-3}}{)}^{-\frac{1}{2}}\text{and}$$(B.2) $$M_J^{th}\approx 0.28{\mathrm{M}}_{\odot }(\frac{\mathrm{T}}{10\mathrm{K}}{)}^{\frac{3}{2}}(\frac{\mathrm{n}}{{10}^6{\mathrm{c}\mathrm{m}}^{-3}}{)}^{-\frac{1}{2}},$$(B.3) where n is the number (volume) density in the system.

Fig. C.1 Distribution of the ratio $M_{core}/M_J^{th}$ observed in the ALMAGAL clumps. Light blue: Entire population of 6218 cores; dark blue: Population of cores in clumps with Ncores ≥ 4.

Appendix C ALMAGAL core masses versus thermal Jeans mass

The estimated thermal Jeans masses in the ALMAGAL clumps range from ∼0.3 to ∼14 M_⊙, a range of masses that is fully observable by the sensitivity of the survey (Molinari et al. 2025), and that is above the completeness limit of ∼0.2 M_⊙ of the ALMAGAL core catalog (Coletta et al. 2025).

Together, the distribution of the ratio $M_{core}/M_J^{th}$ presents a wide spread over about two orders of magnitude, but is heavily skewed towards low values of the ratio. The median value is $M_{core}/M_J^{th}\sim 0.3$ and the first and third quartiles are equal to 0.1 and 0.75, respectively, indicating that most of the ALMAGAL cores have masses considerably lower than $M_J^{\text{th}}$. Specifically, 4987 of 6218, equal to ∼80% of the entire sample, have $M_{core}/M_J^{th}$. Even taking into account the large uncertainties in the core mass and $M_J^{th}$, the result is that a large fraction of cores with sub-Jeans masses are present. Similar results were also found in other surveys, like, for example, the ASHES survey, where a large population of sub-Jeans mass cores has also been detected (Sanhueza et al. 2019; Morii et al. 2023).

Figure C.2 shows the density distribution of the ratio $M_{core}/M_J^{th}$ as a function of Σ_cl and L/M. We do not observe any strong correlation between the clump surface density Σ_cl and the ratio ${\mathrm{M}}_{\text{core}}/{\mathrm{M}}_J^{\text{th}}$. The large majority (∼70%) of the detected cores are identified in clumps with 0.6 ≲ Σ_cl ≲ 2.5 g cm⁻² whose distribution of $M_{core}/M_J^{th}$ does not depend on the clump surface density Σ_cl. They all have a wide distribution that extends from 0.02 to ∼10, with a typical peak around ∼0.1−0.3.

It is more interesting to analyze the distribution of the ratio ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}$ as a function of the L/M ratio. While the ratio ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}$ still shows a general wide spread, we find a broad general trend with the L/M ratio. It is important to note that for clumps in the early stages of evolution, we identify a better agreement of the core masses with the computed thermal Jeans mass. In fact, the distribution of ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}$ for cores in clumps with L/M ≲ 1 L_⊙/M_⊙ peaks at ∼0.6, and shows a relatively narrow distribution with the first and third quartiles equal to 0.3 and 1.33. This agreement is not due to the sample size, as the distribution is derived from 1186 cores, equivalent to ∼20% of the entire sample. On the other hand, for the core separations, we observe that the distribution ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}$ shifts its peak towards lower values and broadens with increasing values of L/M. The peak gradually shifts to ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}\sim 0.1$ for clumps with L/M ≳ 50 L_⊙/M_⊙. Moreover, we note that while the peak shifts, the distribution becomes progressively more right skewed, presenting a considerable tail at high values of the ratio, and increase its width. This observed behavior is due to the combination of two effects. First, we measured an increase in the estimated ${\mathrm{M}}_J^{th}$ with L/M in the sample of ALMAGAL clumps, which shifts the distribution of the ratio ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}$. The increase in ${\mathrm{M}}_J^{th}$ is mainly driven by a general warm-up of the clump material in systems with higher L/M, due to the heating from the already formed protostars. The clump-averaged volume density, n, has a minor impact on the estimated ${\mathrm{M}}_J^{th}$, since the observed ALMAGAL clumps have a distribution of n that spans the same interval independently of L/M. Secondly, Coletta et al. (2025) found that the core mass distribution extends to higher masses in the clumps with L/M, which explains the broadening of the distribution of the ratio in these systems. Significant differences between M_core and ${\mathrm{M}}_J^{th}$ are measured in clumps with L/M>2 L_⊙/M_⊙, where the average ratio ${\mathrm{M}}_{core}/{\mathrm{M}}_J^{th}$ is found to be ≲ 0.3, indicating that the newly formed cores have masses smaller than the thermal Jeans mass. One possible explanation for the abundance of subJeans masses may reside in the increase in the number of cores with evolution. The formation of additional fragments increases the strength of the tidal forces present in the cluster. Those forces are suggested to also promote gravitational instabilities, with the effect of reducing the critical mass required for local collapse (Zavala-Molina et al. 2023)

Appendix D Simulations of the Q parameter in regimes with low number of fragments

The Q parameter introduced by Cartwright & Whitworth (2004) is usually adopted to characterize the spatial distributions of the components of clusters, since its value distinguishes between systems where the elements are centrally concentrated and those where they are arranged in subclusters, which are well represented by a fractal distribution (Cartwright & Whitworth 2004; Parker 2018; Dib & Henning 2019). Cartwright & Whitworth (2004) used synthetic clusters to identify the intervals of the values of the Q parameter that characterize these different spatial distributions. They found that Q ≃ 0.8 separates between the two aforementioned distributions: centrally concentrated clusters have Q>0.8, while subclustered systems have Q<0.8. The threshold value of Q ≃ 0.8 identifies systems whose elements are uniformly distributed.

However, as correctly pointed out by Avison et al. (2023), this classification was defined for star clusters composed of a large number of components, from several tens to hundreds of elements. On the contrary, the typical cluster of cores observed with ALMA is generally composed of a smaller number of sources per field, up to a few tens in the case of the ALMAGAL survey. Poorly crowded clusters with 2 to 15 detected sources per field were found in the TEMPO survey by Avison et al. (2023), who were the first to investigate the reliability of the Q parameter applied to the case of small sample sizes. Their result was that Q parameter can lead to inconclusive results when N_cores<10, as it typically returns values below the threshold value of 0.8 independently of the intrinsic spatial distribution of the cluster elements. Although their results suggest that the Q parameter has a limited applicability to typical ALMA observations, they analyzed only two sets of synthetic clusters with N_core = 5 and 10, as these values represent the degree of fragmentation observed in the TEMPO survey. Here, we extend their work to the case of ALMAGAL clumps, where Coletta et al. (2025) detected between 2 and 49 sources, investigating the reliability of the Q parameter in these cases.

Fig. C.2 Density distribution plot of the core mass normalized to the thermal Jeans mass, MJ, presented as function of the clump L/M (right panel) and surface density, Σc (left panel). The contours identify the different areas for which, going from the external to the internal, a fraction of 10%, 25%, 50%, 75%, 95% of the entire distribution lies outside. The horizontal lines are reference values of the ratios: unity (dashed), 0.1 (dot-dashed), and 5 (dotted).

We followed the prescription of Cartwright & Whitworth (2004) to simulate synthetic clusters with centrally condensed and subclustered spatial distributions for their components. We simulated 3D clusters following the former type of distribution as spherical systems composed of N objects radially distributed with a number density N(r) ∝ r^−p. These are simply generated by creating N triples of polar coordinates (r, θ, φ) with the following formulae: $$r={\left(\frac{(3-p)R}{3}\right)}^{\frac{1}{3-p}},$$(D.1) $$\theta =arccos\left(\right(2\Theta )-1),$$(D.2) $$\varphi =2\pi \Phi .$$(D.3) Here R, Θ, and Φ are numbers extracted randomly with a uniform distribution between 0 and 1.

Simulations of subclustered systems are generated by creating fractal distributions with the box method described in Goodwin & Whitworth (2004) and also adopted by Cartwright & Whitworth (2004) and Cartwright et al. (2006). In short, these patterns are created by initially dividing the region of the cluster into eight cubes with each side equal to R_c/2, placing a first generation of parent seeds in the center of each of these cubes. Each cube is further split into 8 subcubes, whose centers are populated by a first generation of children, which may become fertile. Only fertile elements are kept in each cycle, and, in the next one, those children become parents, spawning a new generation. The iterative process creates the fractal distribution for the cluster that contains a larger or smaller number of substructures depending on the probability function adopted for the fertility. This probability is set by the fractal dimension D as 2^(D−3), allowing fewer children to mature for lower values of the fractal dimension D, and thus generating a system with more substructures. Moreover, to avoid extremely regular structures, we introduce a random displacement of amplitude 0.01 × R_c to the position of the spawning children. Cycles are iterated until a number of cluster components larger than the target N are generated, then we randomly select N of these components to define a single synthetic cluster realization.

We ran seven sets of simulations: three for each of two types of distribution sampling with p = 1, 2, and 2.9 (low, medium, and high concentrations) and with D = 1.5, 2.0, and 2.5, (high, mild, and low level of subclustering), respectively, and one where the components are uniformly distributed in the cluster, which is simply generated adopting p=0 (an approximately similar result is obtained with D=3). Each set is composed of a group of synthetic clusters computed varying N_cores from 4 to 30 with a step of 2. Our simulations are designed over a range of N_cores smaller than that of the ALMAGAL survey, since higher number of fragments rarely occurs in our observations. We computed 10⁵ synthetic realizations for each value of N_cores, projecting the 3D spatial distribution in two dimensions and then deriving the distribution of the Q parameter. The large number of realizations is chosen to carefully determine the probability distribution function of the Q parameter. To speed up the calculation of the 10⁵ realizations, we initially simulated 10⁴ systems with three times the intended number of sources N_cores, and then, for each system, we selected N_cores elements at random ten times.

The simulation campaign is composed of 1.1 M synthetic clusters, and the results are shown in Fig. D.1 where we report how the probability distribution functions (PDFs) of the Q parameter, represented as boxplots at the 10, 25, 75, and 90th percentiles, vary as a function of the number of cluster elements N for different types of spatial distribution. For each measured Q parameter, these PDFs could be theoretically adopted to assess the likelihood of a system being a concentrated cluster or a subclustered structure. However, we prefer to use these results to define the ranges of Q that are most likely associated with one or the other type of spatial distribution.

The plots in Fig. D.1 show that in the case of very small populations (N_core<8) the measured PDFs of different spatial distributions strongly overlap, confirming the result of Avison et al. (2023) that the Q parameter is highly unreliable in these systems. The situation improves for N_core ≥ 10, above which the PDFs of the two types of spatial distributions start to progressively shift one with respect to the other, although they still retain a certain degree of overlap. For this reason, we focus the discussion on the results of cases where N_core ≥ 10.

Fig. D.1 Boxplots representing the distributions of the Q parameter derived from the simulated clusters as a function of the number of cluster components. Each box indicates the median, the first, and the third quartile, while the whiskers extend from the 10th to the 90th percentile of the distribution of measured Q. The left panel compares the results of the simulations with uniform distribution, obtained with p=0, and with the lowest level of subclustering in the simulation campaign, corresponding to a fractal dimension D=2.5. The central panel presents the results for clusters with a radial distribution of their components, for the two cases with exponents p=1 and p=2. The right panel presents the results for clusters with a fractal distribution of their components with two different levels of subclustering, determined by a fractal dimension D=2.0 (mild subclustering) and D=1.5 (high subclustering). The reference values for Q discussed in the text are reported as horizontal lines, corresponding to Q=0.65 and Q=0.83 (dotted lines), 0.75 (dashed line), and the threshold Q=0.8 introduced by Cartwright & Whitworth (2004, thin line).

Fig. D.2 Same as Fig. D.1, but showing the direct comparison of the distributions of the Q parameter computed from the sets of synthetic clusters with the components distributed with a radial power law and with a fractal distribution. The comparison is shown for sets with low, mild, and high concentration degrees and level of subclustering.

These plots show that the criteria introduced by Cartwright & Whitworth (2004) cannot be directly applied to the case of small sample sizes. In particular, the threshold limit Q ≈ 0.8, which should also identify the case with uniform distribution, does not reflect the results obtained from our simulation campaign. Our simulations of the uniform distribution indicate that the median of the probability distribution functions slightly increases from ∼0.73 to 0.76 when the number of cluster components increases. More generally, the PDFs in this limit case are confined in an interval equal to 0.65 ≤ Q ≤ 0.82. Values close to ∼0.8 are only found in the upper tails of the distribution, for values greater than the third quartile, while the average of the medians is equal to Q=0.75. The simulation with uniform distribution is a limit case for the two types of spatial distribution analyzed; therefore, its statistical properties can be adopted as a reference to discuss the results of the other simulations. If we look at the median of the PDFs, we notice that the clusters with centrally concentrated distribution have Q ≳ 0.75, while those with mild and high levels of subclustering have Q<0.75.

However, adopting a defined value as a threshold limit does not take into account that PDFs have widths that extend over wide intervals of Q, and often overlap for different simulated spatial distributions. Figure D.2 shows the PDFs derived from the simulation campaign, comparing the sets of simulations for increasing degrees of concentration or subclustering. The two different spatial distributions may both produce systems with Q values between 0.65 ≤ Q ≲ 0.8, but outside this interval it can be determined with a reasonably high confidence whether the cluster components are distributed as a radial power-law or as subclusters. Simulations with subclustering have PDFs whose 90th percentile is always lower than 0.8, and systems with Q>0.8 are unlikely to have a fractal distribution. Similarly, since the 10th percentile of the probability distribution functions of the clusters with radial distribution is always higher than 0.65, we expect that when Q<0.65 is measured, these systems are most likely to present a level of subclustering.

All Tables

All Figures

Fig. 1 Distribution of the L/M ratio and average surface density, Σcl, as function of the clump heliocentric distance for the ALMAGAL sample, and comparison of these two properties. Each panel includes the 2D density distribution computed dividing the intervals of L/M, Σcl, and dcl into bins of width 0.3 dex, 0.2 dex, and 500 pc, respectively. The solid light gray lines indicate the contour levels of the 2D density distribution derived with kernel density estimation (KDE) corresponding to 10, 30, 50, 70, 95% of the sample.

In the text

Fig. 2 Examples of the distribution of cores found in four ALMAGAL clumps overlapping on the dust continuum map at 1.4 mm from which they were extracted. These fields show examples of clustered systems with circular and elliptical patterns and cases of aligned cores. The contour levels correspond to 4n × σ level and σ equal to the noise level measured on the image. The positions of the cores from the ALMAGAL catalog extracted by Coletta et al. (2025) are shown with gray star markers. In each field we indicate the cluster geometrical center and its radius (see the definition in the text) with a black cross and a dark blue circle, respectively. We also show the convex hull polygon (cyan segments) and the corresponding best-fitting ellipse (black dashed line) adopted for the morphological characterization of the system. The solid thick black segments are the MST edges determined with Prim’s algorithm (Prim 1957), connecting all the cores in the field.

In the text

	Fig. 3 As in Fig. 2. These fields show examples of cores distributed over filamentary features and in well-separated local substructures, both patterns that are sometimes found in ALMAGAL continuum images.
In the text

Fig. 4 Top panel: distribution of the cluster radius measured from the average core positions for the sample of 514 ALMAGAL clumps with Ncore ≥ 4 (in gray) and the subsample composed 347 clumps with dcl ≥ 3.7 kpc (magenta). Bottom panel: clump heliocentric distance vs. measured cluster radius. The dashed line shows the physical size of the ALMAGAL FoV, while the dotted line indicates the area where the sensitivity is constant within a factor of two, i.e., angular sizes ∼13.5′′. Cluster radii exceeding the FoV sizes are possible in systems with a substantial offset between the cluster and the image center.

In the text

	Fig. 5 Distribution of the ratio of the radius of the core cluster to the clump radius (gray histogram). The dot-dashed line, which refers to the top and right axes, shows the distribution of the clump radius derived from Herschel Hi-GAL data at 250 μm from Elia et al. (2021).
In the text

	Fig. 6 Cluster radius as a function of the clump-averaged properties: L/M (left panel) and surface density Σcl (right panel), color-coded by the clump distance. The black lines indicate the contours of the 2D probability density distribution corresponding to 10, 25, 50, 75, and 95% levels.
In the text

Fig. 7 Whisker plots showing the distribution of the two measured estimates of the cluster elongation for the different classes of ALMAGAL clumps identified by the number of cores inside the CH polygon: sparse, where no further core lies inside it; intermediate, where there are fewer cores than the ones defining the polygon itself; clustered, where there are more cores than the ones at polygon vertices. The left panel presents the elongation eCH, defined as the ratio of the major to the minor axis of the best fitting ellipse to the CH polygon, while the right panel refers to the ratio ξ of the circular area of the cluster to the area of the CH polygon. The first and third quartiles, and the 10th and 90th percentiles of the distribution are indicated by the boxes and the segment extremes, respectively, while the central segment corresponds to the median.

In the text

Fig. 8 Left panels: linear separation between ALMAGAL cores as a function of the distance of the ALMAGAL clump for near (top panel) and far configuration sample (bottom panel). The gray shaded areas indicate the ranges of the limiting sizes set by the beams of the ALMAGAL images. The black and colored solid lines (blue and green) indicate the median and the mode of the distribution of the separations, respectively, while the colored area traces the distribution width, defined at the 68% peak intensity. The dashed and dotted lines indicates one-half of the observed field and the full field, respectively. The vertical lines show the adopted distance bins with the relative number of the separations. The 20 separations lying below the spatial resolution are found where the beam is strongly elongated, and are larger than the beam minor axis. In addition, a limited number of separations are not included in the figure as they exceed the selected intervals: 180 in the near sample (162 of which have dcl ≥ 9 kpc), and six in the far sample. Right panels: overall distribution of the projected linear separations in the near (top) and far (bottom) sample, respectively. The thick colored lines (blue and green) show the probability distribution function estimated with the KDE method. The gray shaded area indicates the interval spanned by the circularized beams in the corresponding ALMAGAL images.

In the text

	Fig. 9 Relation between the median of the core separations in each clump and the mode determined with the KDE method. The colors indicate fields with different number of cores, ranging from 4 to 50 in logarithmic scale. The dashed line indicates the 1: 1 relation, while the dotted lines indicate the range of variation of ± 25%.
In the text

Fig. 10 Left panels: density distribution of the deprojected core separations measured in the all clumps with L/M ≤ 1 (orange), 1 ≤ L/M ≤ 10 (green), and L/M ≥ 10 L⊙/M⊙ (purple) with the continuous distribution computed with the KDE method overlapping with the colored continuous line. Clumps of the near and the far are shown in the top and bottom panel, respectively. The colored vertical thick lines show the median of these distributions. Right panels: Deprojected core separations as function of the clump L/M for the systems in the near (top panel) and the far (bottom panel) subsample, respectively. The thick black lines indicate the median of the distributions of the core separations computed in bins of L/M of width 0.3 dex. The thick colored lines, in blue and in green for the near and far subsample, respectively, indicate the peak of these distributions, while the colored area traces their spread defined as the width at 68% of peak intensity. The dashed lines show the variation of the mode (in light blue and green) and of the median computed over subsamples having 3.7 ≤ dcl ≤ 4.7 kpc and 4.7 ≤ dcl ≤ 6.5 kpc. The numbers at the top give the size of the subsamples from which the distribution in each bin is computed. The colored contours trace, from the outside inward, fractions of the sample lying outside that rise from 15% to 95%, with a step of 10% for the near (in blue) and far (in green) subsample, respectively. The median and the peak derived for the near sample are also indicated in the bottom panel for comparison purposes.

In the text

Fig. 11 Relation between the clump-averaged separations between the cores provided by the median and the L/M ratio, color-coded to indicate the number of cores detected in the clump. To support the visualization, different sizes are adopted to identify systems 4<Ncores<8 (low), 8 ≤ Ncores<16 (intermediate), and Ncores ≥ 16 (high). The different symbols refer to the clumps of the near (circle) and far (diamond) sample. The solid lines show the median of the distribution of the separations computed in bins of L/M of width 0.4 dex for the near (blue) and far (green) sample, respectively, while the black dashed line corresponds to the median of all the clumps with Ncores ≥ 4.

In the text

Fig. 12 Deprojected core separations against the clump surface (top panel) and volume (bottom panel) density. The distributions of the core separations are computed in bins of Σcl and ncl of width 0.15 and 0.2 dex, respectively. The peaks and spreads of these distributions, defined as the width at 68% peak intensity, are indicated by lines and area for the near (blue) and the far (green) samples, respectively. The medians of the distributions for these two sets are indicated as dashed and dotted black lines, respectively. The best fitting power-law relation computed over the entire set of separations and whose results are reported in Table 5 are indicated as light gray and dark gray lines for the near and far sample, respectively.

In the text

Fig. 13 Comparison between the clump-averaged estimators for the core separations, specifically the median (top panel) and mode (bottom panel) of the distribution in each clump, and the computed Jeans length. The fragmentation level in the clump is indicated by the symbol sizes, with increasing size going from low (4 ≤ Ncores<8) and intermediate (8 ≤ Ncores<16), up to high (Ncores ≥ 16) core numbers. The dashed lines are drawn as a reference and correspond to the separation equal to ${\lambda }_J^{th}$ (thick line) and 2 × and $0.5\times {\lambda }_J^{th}$ (thin lines), respectively. The gray shaded boxes indicate the amplitude of the distributions of the ratio $l_i^{depr}/{\lambda }_J^{th}$ (see text), and correspond to the interquartile range (dark gray) and the 10th to 90th percentile range (light gray), respectively.

In the text

Fig. 14 Left panel: comparison of the ratio R of the deprojected separations to ${\lambda }_J^{th}$ and the clump L/M ratio for the near (top panel) and far sample (bottom panel), respectively. The black dashed line shows the median of the distributions of R computed in bins of L/M equally spaced in logarithm of width 0.4 dex. The thick colored lines and the colored area trace the mode and the spread (defined as the width at 68% of the peak) of these distributions, adopting the blue and green colors for the near and far sample, respectively. The numbers at the top give the size of the subsamples adopted to compute the distribution in each bin. The gray shaded area shows the interval 0.5 ≤ R ≤ 1.5, corresponding to an uncertainty on the ${\lambda }_J^{th}$ of a factor of 2. The colored contours trace, from the outside inward, fractions of the sample lying outside that rise from 15% to 95% with a step of 10% for the near (in blue) and the far (in green) subsample, respectively. Right panel: normalized density histogram and estimated probability density distribution estimated with the KDE method (colored continuous lines) of the ratio $R=l_i/{\lambda }_J^{th}$ measured in the three subsamples selected by different intervals of L/M tracing clumps in different evolutionary phases: young (L/M ≤ 1 L⊙/M⊙ in orange), intermediate (1 ≤ L/M ≤ 10 L⊙/M⊙ in green), and evolved (L/M ≥ 10 L⊙/M⊙ in purple).

In the text

Fig. 15 Relation between the Q parameter and the L/M ratio of ALMAGAL clumps which host at least ten cores. The two dashed lines values Q=0.65 and Q=0.8 defining the intervals for which it is most likely that the hosted cores follow a fractal type of subclustering (Q ≤ 0.65) and radial distribution with a central concentration (Q ≥ 0.8), as indicated by the Monte Carlo simulations in Appendix D. The dotted line points to the Q=0.75, which separates the two regimes with a lower statistical significance. Estimates from the ASHES (dark gray; Sanhueza et al. 2019) and ASSEMBLE surveys (light gray; Xu et al. 2024) are included, and are shown as star markers.

In the text

Fig. 16 Distribution of the ALMAGAL clumps hosting at least 10 cores in different bins of the ratio L/M. The colored histograms show the systems were a signature of MS is found, with different colors referring to the different intervals of L/M ratio, i.e., orange L/M ≤ 2 L⊙/M⊙, green 2<L/M ≤ 10 L⊙/M⊙ and purple L/M>10 L⊙/M⊙. The solid black line, which refers to the right y-axis scale, indicate the relative fraction of systems with MS with respect to the observed ones in bins of L/M, while the dot-dashed line shows the average occurrence fraction of MS systems in the three intervals of L/M.

In the text

Fig. 17 Boxplots showing the distribution of maximum value of the mass segregation ratio, ${\mathrm{\Lambda }}_{MSR}^{max}$, divided for the different values $N_{MST}^{crit}$ that identify the interval where the MSR profile is statistically higher than unity, in the 110 ALMAGAL clumps where signatures of MS are found. Those clumps are divided into three groups according to the L/M ratio composed of 59, 83, and 154 systems. The rectangular area of the boxplots shows the interquartile intervals, while the whiskers represent the 10th to 90th percentile intervals; the diamond symbols indicate the cases outside this interval.

In the text

	Fig. A.1 Scheme of projection of a pair of randomly oriented objects in 3D space
In the text

	Fig. C.1 Distribution of the ratio $M_{core}/M_J^{th}$ observed in the ALMAGAL clumps. Light blue: Entire population of 6218 cores; dark blue: Population of cores in clumps with Ncores ≥ 4.
In the text

Fig. C.2 Density distribution plot of the core mass normalized to the thermal Jeans mass, MJ, presented as function of the clump L/M (right panel) and surface density, Σc (left panel). The contours identify the different areas for which, going from the external to the internal, a fraction of 10%, 25%, 50%, 75%, 95% of the entire distribution lies outside. The horizontal lines are reference values of the ratios: unity (dashed), 0.1 (dot-dashed), and 5 (dotted).

In the text

Fig. D.1 Boxplots representing the distributions of the Q parameter derived from the simulated clusters as a function of the number of cluster components. Each box indicates the median, the first, and the third quartile, while the whiskers extend from the 10th to the 90th percentile of the distribution of measured Q. The left panel compares the results of the simulations with uniform distribution, obtained with p=0, and with the lowest level of subclustering in the simulation campaign, corresponding to a fractal dimension D=2.5. The central panel presents the results for clusters with a radial distribution of their components, for the two cases with exponents p=1 and p=2. The right panel presents the results for clusters with a fractal distribution of their components with two different levels of subclustering, determined by a fractal dimension D=2.0 (mild subclustering) and D=1.5 (high subclustering). The reference values for Q discussed in the text are reported as horizontal lines, corresponding to Q=0.65 and Q=0.83 (dotted lines), 0.75 (dashed line), and the threshold Q=0.8 introduced by Cartwright & Whitworth (2004, thin line).

In the text

	Fig. D.2 Same as Fig. D.1, but showing the direct comparison of the distributions of the Q parameter computed from the sets of synthetic clusters with the components distributed with a radial power law and with a fractal distribution. The comparison is shown for sets with low, mild, and high concentration degrees and level of subclustering.
In the text