ABSTRACT

The Orion Nebula Cluster (ONC) is the most massive region of active star formation within a kpc of the Sun. Using Gaia DR3 parallaxes and proper motions, we examine the bulk motions of stars radially and tangentially relative to the cluster centre. We find an age gradient with distance to the stars in the ONC, from 385 pc for the oldest stars, to 395 pc for the younger stars, indicating that the star-forming front is propagating into the cloud. We find an organized signature of rotation of the central cluster, but it is present only in stars younger than 2 Myr. We also observe a net infall of young stars into the centre of the ONC’s deep gravitational potential well. The infalling sources lie preferentially along the filament; on the other hand, outflowing sources are distributed spherically around the cluster, and they have larger velocity dispersion. We further propose a solution to a long-standing question of why the ONC shows a weak signature of expansion even though the cluster is likely bound: much of this expansion may be driven by unstable N-body interactions among stars, resulting in low-velocity ejections. Though analysing signatures imprinted on stellar dynamics across different spatial scales, these observation shed new light on the signatures of formation and evolution of young clusters.

1 INTRODUCTION

The gravitational fields of star-forming molecular clouds can influence the dynamics of both stars and gas over significant spatial scales. In the solar neighbourhood, the most massive of these star-forming regions is the Orion Complex. It contains a number of clusters, the most massive of which is the Orion Nebula Cluster (ONC), alone containing well over 4000 stars within the radius ∼10 pc of each other, of which ∼2000 are concentrated within a radius of ∼2 pc within the Trapezium (e.g. McBride & Kounkel 2019).

Despite the masses and densities involved in assembling some of the massive star-forming regions like Orion, so far there have been only a few significant signatures of gravitational pull of the young populations on their surroundings. Kuhn et al. (2019) have analysed a sample of 28 young clusters; they found that almost none of them appear to have rotational support, and most of them show a preference of expansion with a typical velocity of ∼0.5 km s−1. Only two clusters, M17 and NGC 6231, have shown some signature of contraction, albeit at a low level of significance. Most importantly, when examining multiple subclusters that are all part of the same population strung together along a filament, they found only a little evidence of convergent motion within them, concluding that they are not involved in a hierarchical assembly that would collect smaller star clumps into a larger cluster.

Over the years, most of the efforts have been concentrated on the dynamical state of the ONC based on its proper motions. The focus has been primarily aimed almost exclusively on the central cluster within it, rather than the full ‘head’ of the Orion A molecular cloud. Early works with motions of stars measured from the photographic plates by Jones & Walker (1988) have concluded that the cluster is most likely not in virial equilibrium, but, as the escape velocity they estimated is comparable or smaller than the velocity dispersion (σv) they measured, they considered the cluster to be only weakly bound. Hillenbrand & Hartmann (1998) were able to fit the spatial distribution of the stars in the ONC well using the King model profile, and they have argued that while the entire cluster may not necessarily be relaxed, its core is expected to be. They found that the stellar mass within the cluster is only 40 per cent of what is expected for virial equilibrium, but that the remaining gas may account for the remainder, and that it may yet form a sufficient number of stars for the cluster to remain bound at later stages in its life. Later, Dzib et al. (2017) have used radio interferometry to improve astrometric measurements of stars, and they found no indication of either preferred radial or tangential motion of stars, which would be expected of a virialized population. Finally, since the release of Gaia DR2, several works, such as Kounkel et al. (2018), Kuhn et al. (2019), or Getman et al. (2019), have re-evaluated the dynamics of the ONC, with overall conclusions being that the cluster does appear to be bound, although it may show a weak preference towards expansion, mainly along its outer edges.

In this paper, we examine the dynamical state of Orion, the gravitational influence that it has on its vicinity, and the expansion of the stars within it. To clarify the terminology, we refer to the entire ‘head’ as the ONC, and we would frequently use Trapezium to refer to the central cluster core. In Section 2, we present the data used and the general methods of our analysis. Section 3 presents the primary results of our study. In Section 4, we discuss the implications of our findings of an expansion signature in the ONC despite the cluster likely being bound. We conclude with a summary in Section 5.

2 DATA AND METHODS

2.1 Data selection

One of the challenges that the previous works have had in analysing the dynamical state of the ONC is the sheer number of stars contained within it, as it makes it difficult to fully visualize the velocity vectors of stars. As such, although the cluster appears to be virialized on average, some of the substructure may be less apparent when everything is viewed as a whole. More than that, the ONC contains multiple generations of stars, some younger than 1 Myr, some with ages >4 Myr. While the distributions of these stars of different ages overlap, they do have their own unique features, and this may influence the kinematics as well (Beccari et al. 2017). To mitigate both of these issues, we split up the sample of the stars in the ONC based on their age.

The ages of young stars can commonly be determined from the placement of their photometry on the HR diagram. However, some confusion can remain regarding the ages of unresolved binary stars which have a higher combined luminosity in comparison to single stars; without taking this into account, the ages that would be inferred photometrically will be underestimated. Uncertainty in extinction also has a significant effect on the precise age determination. Although it is still possible to observe trends and gradients in age across given populations (Beccari et al. 2017), such age bias can often be difficult to quantify.

Instead, we use log g as a proxy for age. As a young star evolves and its radius shrinks, its log g increases until it reaches main sequence. Such measurement is a more direct method of determining the age of a star, as it is significantly less biased by multiplicity. To date, only a few studies have produced reliable log g values for young stars that have sufficient accuracy and precision to use them as an age proxy However, recently, Olney et al. (2020) and Sprague et al. (2022) have performed analysis of APOGEE spectra to use data-driven techniques to measure Teff and log g that can be compared to the pre-main-sequence isochrones. Ages estimated from log g are not affected by extinction (as log g is derived from the shape of the spectral lines rather than from the luminosity; moreover, APOGEE spectra are in infrared, where the effect of the extinction is significantly smaller than in the optical regime). Similarly, log g measurements appear to be largely unaffected by unresolved binaries. Olney et al. (2020) have found that among the pre-main-sequence stars, there is no systematic difference in log g distribution of spectroscopic binaries (both SB1 and SB2) relative to the total sample of young stars. Kounkel et al. (2021) have found that there is a slight bias in log g distribution of the SB2s found among the more evolved main-sequence stars; this bias is ∼0.1 dex in magnitude – comparable to the magnitude of the typical uncertainties. However, even if this bias is present in younger stars, such as the sample in Orion, given that SB2s are relatively rare in comparison to all unresolved binaries, spectroscopic ages estimates should on average be more robust.

We do note that using log g rather than luminosity does not compensate over some of the other fundamental factors of the systematic uncertainty in deriving stellar ages, including the uncertainty in the absolute calibration of the stellar parameters, as well as the ability of the isochrone models in accurately representing said parameters. This can often be a challenge, particularly among the cooler stars (e.g. David et al. 2019).

We use the spectroscopic sample of stars observed by APOGEE towards the ONC (Kounkel et al. 2018), and following a crossmatch with Gaia EDR3 (Gaia Collaboration 2021), we further restrict it via

  • Spatial cuts in RA and Dec. 82° < α < 85°, −6.5° < δ < −3.5°.

  • Cut in parallax 2 mas < π < 3.5 mas, which translates to a distance cut of 285–500 pc.

  • Proper motions within 6 mas yr−1 of the mean of the cluster, with the mean motion in the local standard of rest1 μα, lsr = 0.8 mas yr−1 and μδ, lsr = 2.8 mas yr−1. This translates to a window ∼11 km s−1 or ∼5–6 times the velocity dispersion of the clusters, excluding any fast dynamically ejected runaways such as those in McBride & Kounkel (2019), or the bulk of the potential contamination. The APOGEE sample towards the ONC primarily contains bona fide members, with only a small fraction of likely field stars, as such, this cut eliminates only 5 per cent of stars.

  • Spectroscopically determined Teff < 6500 K and log g> 3 dex, restricting the sample to the stars from which spectroscopic ages can be reliably determined.

This produces a sample of 1612 stars, of which 893 are found towards the Trapezium.

The log g values derived via APOGEE Net may still have minor systematic differences relative to various isochrones, such as e.g. PARSEC (Marigo et al. 2017), particularly at the cooler end of Teff < 3800 K. The isochrones suggest that log g is expected to decrease with decreasing Teff for low-mass young stars, whereas, in a given cluster of a given age, the measured log g values from Olney et al. (2020) may appear to slightly increase with decreasing Teff mirroring a relation that is seen among the main-sequence stars. As such, applying theoretical isochrones directly results in a significant overestimation of the ages of cool stars in comparison to their hotter counterparts. Without precise empirical isochrones that could be used to characterize this dependence, we approximate a simple monotonic relationship between age and log g without a dependence on Teff, as this offers a better consistency within a given cluster in comparison to the theoretical isochrones.

We roughly separate the sample into five bins, and translate these log g cuts into approximate age ranges: log g< 3.6 dex (≲1 Myr, 125 stars), 3.6 < log g < 3.8 dex (1–2 Myr, 265 stars), 3.8 <log g < 4.0 dex (2–3 Myr, 480 stars), 4.0 <log g < 4.2 dex (3–5 Myr, 445 stars), and 4.2 < log g dex (≳5 Myr, 297 stars). The typical precision in log g in the sample is ∼0.1 dex (10th and 90th percentile of the formal error distribution are 0.037 dex and 0.16 dex). We note that these ages are only approximate labels to provide the readers with some intuition regarding the star formation history of the ONC, but fundamentally, these remain to be cuts solely in log g rather than age directly.

There are independent age estimates from photometry for some of the stars against which the spectroscopic ages can be evaluated, e.g. ages of the individual stars from Kounkel et al. (2018) derived using Teff and bolometric luminosity relative to the PARSEC isochrones (Marigo et al. 2017), or the ages based on Gaia and 2MASS photometry interpolated using Sagitta neural net (McBride et al. 2021). There are some systematic differences, both between photometric ages estimated form these two approaches, as well as relative to the spectroscopic age estimates used here. However, in both cases, there is a strong trend of older stars having systematically higher log g. As such, while there may be some cross-contamination of the neighbouring bins, there are clear differences between the selection of the oldest and the youngest age bins, regardless of the method used to select them. The overall results presented in this paper are robust against the age metric that is used.

2.2 Analysis

To examine the motion of stars in the ONC in the reference frame at the centre of the cluster, we derive a metric representing the relative radial orientation of motion, cos θ, where θ is defined as
(1)
where |$\operatorname{atan2}$| is the 2-argument arctangent that is capable of returning values across the full range of a circle, α and δ are the right ascension and declination (in decimal degrees) of the star in question, μα, lsr and μδ, lsr are the star’s proper motions in right ascension and declination with respect to the local standard of rest (Schönrich et al. 2010), expressed in milli-arcseconds per year, αCOM and δCOM are the right ascension and declination of the ONC’s centre of mass (αCOM = 83.8° and δCOM = −5.4°), and μα, COM, lsr and μδ, COM, lsr are the mean proper motions of ONC members in right ascension and declination with respect to the local standard of rest (μα, COM, lsr = 0.8 mas yr−1 and μδ, COM, lsr = 2.8 mas yr−1). The centre of mass position and the mean velocities are approximations derived based on the sample, and they are used to evaluate intracluster dynamics. The typical uncertainty in θ from the reported uncertainty in μα and μδ is ∼3°, increasing to ∼10° if we consider a reasonable uncertainty in the reference position and velocity of the cluster. However, the results outlined in this work do not strongly depend on these precise averages; they are robust against the uncertainties, and they remain consistent when other reasonable estimates are used instead.

We then segregate the selected members of the ONC into three categories based on their computed cos θ values: cos θ > 0.5 are the sources that are preferentially moving away from the cluster centre (566 stars), cos θ < −0.5 are those moving towards the cluster centre (606 stars), and −0.5 < cos θ < 0.5 are those moving tangentially around the cluster/rotating around it (440 stars). The latter sources can be further subdivided through tangential orientation, with sin θ < 0 as the sources preferentially moving clockwise (234 stars), and sin θ > 0 as those moving counterclockwise (206 stars).

3 RESULTS

To evaluate the distribution of velocities in the ONC, as a comparison, we generated a random sample of proper motions for all of the stars, drawing both μα and μδ from the 0.9 × 1.1 mas yr−1 Gaussian velocity dispersion of the cluster (consistent with the underlying distribution of sources in the real data), preserving the positions of each star. We processed these synthetic velocities in a similar manner, and applied similar cuts. We find that in a purely randomly drawn velocities, there is expected to be an approximately equal number of infalling, outflowing, and rotating stars. In comparison, real data for the entire ONC show that the fraction of sources that are rotating is suppressed, and, instead, there is an excess in the number of stars that are infalling (Fig. 1). Applying two-sided Kolmogorov–Smirnov (KS) test between the full distribution of all alignments between the real and the synthetic data shows that the alignment of stars relative to the cluster centre relative to the random is different at 2.5σ level. Furthermore, zooming in only on the central cluster, within 0.4° of the cluster centre, we find there to be an excess in the number of the stars expanding outwards from the centre. KS test shows that the distribution of velocities is different from random at 3σ level. The trends shown in Fig. 1 are not affected by the precise cut on the cos θ, e.g. changing the threshold from ±0.5 to other thresholds does not substantively affect the significance of these trends. Similarly, uncertainty in θ does not skew the underlying distribution.

Number of sources that are infalling into the ONC (−1 < cos θ < −0.5), those that are preferentially moving away (0.5 < cos θ < 1), and those that are moving tangentially around the cluster (−0.5 < cos θ < 0.5). Note that although the bin size extends beyond ±1, no objects are found beyond these limits, and the geometry of the cos function is such that these splits at ±0.5 are expected to divide a random sample into three bins of comparable sizes. The left-hand panel shows the data for the entire ONC, and the right-hand panel just zooms in on the sources found in the inner 0.4° of the central cluster. The blue histograms show the orientation of motion of the real stars, red – stars with randomly generated velocities drawn from a Gaussian distribution consistent with the velocity dispersion of the cluster. The uncertainties are estimated as through Poisson approximation of sqrt(n).
Figure 1.

Number of sources that are infalling into the ONC (−1 < cos θ < −0.5), those that are preferentially moving away (0.5 < cos θ < 1), and those that are moving tangentially around the cluster (−0.5 < cos θ < 0.5). Note that although the bin size extends beyond ±1, no objects are found beyond these limits, and the geometry of the cos function is such that these splits at ±0.5 are expected to divide a random sample into three bins of comparable sizes. The left-hand panel shows the data for the entire ONC, and the right-hand panel just zooms in on the sources found in the inner 0.4° of the central cluster. The blue histograms show the orientation of motion of the real stars, red – stars with randomly generated velocities drawn from a Gaussian distribution consistent with the velocity dispersion of the cluster. The uncertainties are estimated as through Poisson approximation of sqrt(n).

3.1 Infall

Figs 2 and 3 show the motions of infalling, expanding, and rotating sources in each of the five age bins. The stars in the oldest age bin (≳5 Myr) predate the formation of the central cluster as a whole. Although there are stars that fall into various cuts of cos θ, there is no organized expansion or rotation to speak of; they appear to primarily trace random motions of stars.

Proper motions of stars towards ONC, separated into five age bins, shown in separate rows, based on their log g, with the velocity vectors colour-coded based on the orientation relative to the centre: red for infall, blue for expansion, magenta for clockwise rotation, green for counterclockwise rotation. The length of each vector corresponds to the motion of a star over the next 0.2 Myr. Typical velocities are 1.2 mas yr−1 or 2.2 km s−1. The black dot shows the assumed centre of the cluster. The first column shows the full sample, second column shows only the sources that are falling towards the centre, third column shows only the sources outflowing from the centre, and the fourth column shows the sources that are preferentially moving tangentially around the centre.
Figure 2.

Proper motions of stars towards ONC, separated into five age bins, shown in separate rows, based on their log g, with the velocity vectors colour-coded based on the orientation relative to the centre: red for infall, blue for expansion, magenta for clockwise rotation, green for counterclockwise rotation. The length of each vector corresponds to the motion of a star over the next 0.2 Myr. Typical velocities are 1.2 mas yr−1 or 2.2 km s−1. The black dot shows the assumed centre of the cluster. The first column shows the full sample, second column shows only the sources that are falling towards the centre, third column shows only the sources outflowing from the centre, and the fourth column shows the sources that are preferentially moving tangentially around the centre.

Same as Fig. 2, but zoomed in on Trapezium. The length of the vectors is decreased to the distance covered in 0.1 Myr.
Figure 3.

Same as Fig. 2, but zoomed in on Trapezium. The length of the vectors is decreased to the distance covered in 0.1 Myr.

However, at the oldest (≳5 Myr) age bin, there does appear to be a significant infall of stars moving towards the Trapezium from south of the ONC, in the vicinity of NGC 1980 (we note that each panel shows the current proper motions of the stars of that age, not necessarily the motions at the time of their formation). This infall becomes even more prominent in the next age bin (3–5 Myr); furthermore, there does appear to be significant infall from the northern part of the ONC, NGC 1977. This infall continues to persist in the younger age bins as well, until the density of sources in these regions deplete, as NGC 1980 and 1977 had stopped actively forming stars a few Myr ago. The infalling sources constitute the bulk of sources found north and south of the Trapezium; however, there are next to no infalling sources from either east or west of it (Fig. 4, top row). While we do not have proper motions for the molecular gas, given the similarity between radial velocity of the gas and stars that is often observed (e.g. Da Rio et al. 2017), this suggests that the filament that has produced the ONC is contracting due to the gravitational potential of both the Trapezium and the entire cluster as a whole.

Top: map of the ONC showing fraction of sources in a given healpix that are preferentially falling into the centre of the ONC, expanding outwards, and those that are preferentially moving tangentially around it. Bottom: same as above, but the velocities for each star have been generated randomly from a Gaussian distribution representing the cluster, to highlight the differences between the real data and the null hypothesis.
Figure 4.

Top: map of the ONC showing fraction of sources in a given healpix that are preferentially falling into the centre of the ONC, expanding outwards, and those that are preferentially moving tangentially around it. Bottom: same as above, but the velocities for each star have been generated randomly from a Gaussian distribution representing the cluster, to highlight the differences between the real data and the null hypothesis.

Given that the central cluster in the ONC is several times more massive and more than an order-of-magnitude denser than either NGC 1980 or 1977 (e.g. Megeath et al. 2016), its potential provides the singularly most dominant force in the region. However, it is difficult to state definitively when the infall has begun: whether it has started when the cluster has accreted a large fraction of its mass, whether it was when it started forming its first stars, or whether it predates their formation and originated in the gas. However, given the dynamics of other stars in the Orion Complex (Kounkel et al. 2018), and the known sources of stellar feedback in the region (e.g. winds from θ1 Ori C), it is unlikely that anything other than gravity would be responsible for the observed infall in these regions.

In contrast, if we examine the sample with randomly generated proper motions (Fig. 4, bottom), the distribution of sources that would have orientation of their proper motions suggestive of either infall, outflow, or rotation is expected to be more homogeneous and uniform across the cluster. Comparing the distribution of fractions of stars that are infalling per healpix (i.e. distribution shown in Fig. 4, considering a set of values across all healpix within a given range of δ) in the real data versus the random sample via a two-sided KS test shows that the two populations are distinct at >4σ level. The difference becomes more pronounced when segregating the ONC into three portions: the top (δ ≳ −5°), middle (−5.8° ≲ δ ≲ −5°), and bottom (δ ≲ −5.8°). The resulting KS-test statistics are shown in Table 1.

Table 1.

KS test statistics comparing the fraction of sources that are expanding, infalling, or rotating, in the ONC compared to the random sample. The asterisks highlight the instances that appear to be most significantly discrepant at >3σ.

InfallExpansionRotation
Full1.31e-05 (4.3σ*)0.030 (2.1σ)0.110 (1.6σ)
Top0.00600 (2.7σ)0.180 (1.3σ)0.0870 (1.7σ)
Middle0.0775 (1.8σ)0.00255 (3.0σ*)0.361 (0.9σ)
Bottom7.01e-07 (5.0σ*)0.0362 (2.0σ)0.514 (0.6σ)
InfallExpansionRotation
Full1.31e-05 (4.3σ*)0.030 (2.1σ)0.110 (1.6σ)
Top0.00600 (2.7σ)0.180 (1.3σ)0.0870 (1.7σ)
Middle0.0775 (1.8σ)0.00255 (3.0σ*)0.361 (0.9σ)
Bottom7.01e-07 (5.0σ*)0.0362 (2.0σ)0.514 (0.6σ)
Table 1.

KS test statistics comparing the fraction of sources that are expanding, infalling, or rotating, in the ONC compared to the random sample. The asterisks highlight the instances that appear to be most significantly discrepant at >3σ.

InfallExpansionRotation
Full1.31e-05 (4.3σ*)0.030 (2.1σ)0.110 (1.6σ)
Top0.00600 (2.7σ)0.180 (1.3σ)0.0870 (1.7σ)
Middle0.0775 (1.8σ)0.00255 (3.0σ*)0.361 (0.9σ)
Bottom7.01e-07 (5.0σ*)0.0362 (2.0σ)0.514 (0.6σ)
InfallExpansionRotation
Full1.31e-05 (4.3σ*)0.030 (2.1σ)0.110 (1.6σ)
Top0.00600 (2.7σ)0.180 (1.3σ)0.0870 (1.7σ)
Middle0.0775 (1.8σ)0.00255 (3.0σ*)0.361 (0.9σ)
Bottom7.01e-07 (5.0σ*)0.0362 (2.0σ)0.514 (0.6σ)

We do note that ONC is not an isolated system, and that it is all that can also be affected by the gravity of the larger Orion A molecular cloud, as well as gravity of the Orion Complex as a whole, and these forces could also influence the dynamics of the ONC. However, on the scales of a few pc (which is the scale over which we observe the infall), self-gravity of the cluster would dominate over the gravity of the extended structure.

3.2 Expansion

Stars with age <5 Myr begin to trace the central overdensity associated with the Trapezium, and the younger they are, the more centrally concentrated they become. Around the cluster core (including along the OMC 2/3 filament), there is little to no ‘organized’ infall of stars towards the centre. However, as soon as the central cluster begins to develop, there is a significant signature of both rotation and expansion, overshadowing the infalling sources that are there. This is consistent with the simulations, where infall is easier to definitively observe at wider separations from a massive cluster than in its immediate proximity (Kuznetsova, Hartmann & Ballesteros-Paredes 2015).

The bulk of the expanding sources can be traced back directly to the cluster centre and they tend to be distributed spherically around it. However, the expanding sources dominate the sample around the outer edges on the east and the west side of the cluster due to this spherical distribution, as there is no overlap with the infalling sources that are only found along the filament (Fig. 4). That preference for the outer edges of the cluster to be dominated by expansion has previously been noted by Kounkel et al. (2018), Kuhn et al. (2019), and Getman et al. (2019).

We do note that at the ages of 3–5 Myr, the sources that appear to move away from the cluster centre have their motion preferentially due north, with very few sources moving due south, in contrast to the younger age bins where the expansion does appear to be more spherically symmetric, or in contrast to the sources moving towards the cluster centre or those moving tangentially around it in the similar 3–5 age bin which also appear to be more uniformly distributed. It is unclear why such an asymmetry may have arisen, as there is no obvious bias in targeting of such sources in the underlying survey. Alternatively, the mean cluster position and velocity were carefully selected using the cluster as a whole, but there may have been a slight evolution in their precise placement at earlier epochs.

3.3 Rotation

Examining the rotation around the Trapezium, within 0.4° of the cluster centre, in sources that are preferentially moving tangentially, there is little to no preference in the direction among the stars that are older than 2 Myr. However, there is a strong preference for clockwise motion in younger stars, by a factor of ∼2 among 1–2 Myr old stars, and by a factor of ∼3 among <1 Myr old stars (Fig. 5). This preference for clockwise motion in youngest stars does not depend on the precise point adopted as the cluster centre around which rotation is evaluated. However, the excess is strongest when evaluated in the centre of the Trapezium, and it weakens the further it is displaced.

Distribution of direction of motion of sources within 0.4° around the cluster centre that preferentially exhibit tangential motion. Orientation <0 shows the sources moving clockwise; >0 are those moving counterclockwise. Sources are separated into four age bins.
Figure 5.

Distribution of direction of motion of sources within 0.4° around the cluster centre that preferentially exhibit tangential motion. Orientation <0 shows the sources moving clockwise; >0 are those moving counterclockwise. Sources are separated into four age bins.

It is difficult to say, however, if the older stars used to have a preferred orientation that has since been washed out over their lifetime through the dynamical evolution of the cluster, with the younger stars being the only ones with some memory of it, or if the angular momentum of the cluster has evolved over time as the cluster grew in such a manner as to develop a semicoherent rotation only for the currently forming stars.

Rotation of the molecular gas in ONC has been previously predicted in the toy model of the formation of the cluster from Hartmann & Burkert (2007), and it is expected that the stars would inherit the angular momentum from the gas. A weak signature of rotation of the central cluster has also been reported by Theissen et al. (2022).

Among the youngest (<1 Myr) sources that are rotating around the cluster centre, the specific angular momentum can be most easily seen to be conserved, with tangential velocity inversely proportional to the radius at which these stars are found (Fig. 6). From this constant r × vrot, it is possible to estimate the specific angular momentum of ∼7.4 × 1020 cm2 s−1, with both r and v being 2D projections in the plane of the sky of these youngest sources rotating around the cluster. This specific angular momentum is a factor of ∼1.8 higher than the specific angular momentum of the Orion B molecular cloud (Hsieh et al. 2021). However, in the older stars (including 1–2 Myr), the coherence between vrot and r is no longer apparent, that is to say, v has a significant scatter relative to r which washes away any trends.

Tangential velocity as a function of radial distance of the stars with age <1 Myr that are preferentially rotating around the cluster centre. The fitted blue line follows the relation of v = 2/r. The uncertainties in vrot include propagated errors in proper motions, parallax, and angle; most are smaller than the symbol size.
Figure 6.

Tangential velocity as a function of radial distance of the stars with age <1 Myr that are preferentially rotating around the cluster centre. The fitted blue line follows the relation of v = 2/r. The uncertainties in vrot include propagated errors in proper motions, parallax, and angle; most are smaller than the symbol size.

3.4 Radial velocity component

Examining the distribution of the stars in the ONC in the plane of the sky and their proper motions allows for an incredible precision in inferring the dynamical state of the cluster, however, some leverage can be also gained via examining distance of the stars and the radial velocity. However, as even with Gaia EDR3, parallaxes can be very uncertain, resulting in a large spread in the inferred distance; we limit the sample only to the sources with σπ < 0.04 mas, which at the distance of the ONC translates to ∼6 pc, which is still considerable, as it is comparable to the size of the cluster in the plane of the sky, but, none the less, allows to resolve some structure along the line of sight. However, imposing this constraint on σπ significantly limits the sample to only 533 stars, particularly towards the central cluster due to large degree of nebulosity in the region degrading the quality of the parallaxes.

We examine the distance versus the vertical extent of the cluster in δ converted to physical units in Fig. 7, highlighting the velocities of the sources in these two respective dimensions. The assignment of infalling and outflowing sources is retained from the previous figures based on their plane of the sky velocities. On a first glance, the outflowing sources in this plane appear to be somewhat different than in Fig. 2; we note that this is primarily due to the strict σπ cut preferentially excluding sources in the region where expansion is most strongly apparent.

Same as Fig. 2, but showing distance from the Sun versus the position along the filament in δ converted to physical units, with vectors converted from the μδ and radial velocities. The length of the vectors corresponds to the distance covered in 0.4 Myr. Typical uncertainties in distance are ∼6 pc.
Figure 7.

Same as Fig. 2, but showing distance from the Sun versus the position along the filament in δ converted to physical units, with vectors converted from the μδ and radial velocities. The length of the vectors corresponds to the distance covered in 0.4 Myr. Typical uncertainties in distance are ∼6 pc.

When examining the full sample of infalling and outflowing sources only in the radial velocity space, where the cut on parallax quality is not necessary, we find that the stars that are infalling have a good agreement with RVs of the gas along all δ in the cluster. On the other hand, the outflowing stars not only often have a larger σv (Table 2, Fig. 8), they also may be offset from the gas. This is most strongly apparent at δ ∼ −5°, where the peak of the RV distribution of outflowing stars is blueshifted relative to the gas. The origin of this blueshifted stellar component has long since been questioned (Fűrész et al. 2008; Tobin et al. 2009; Kounkel et al. 2016). We can finally offer a partial explanation: radial velocities of these stars could have been dynamically processed by the central cluster. The mean velocity of the expanding stars is more comparable to the mean velocity of the central cluster; on the other hand, the gas north of the cluster is intrinsically redshifted. As such, the total distribution of the stars formed from the molecular gas in that region coupled with the stars that most likely originate in the central cluster but since have been scattered to the same line of sight becomes asymmetric relative to the gas.

Radial velocity of the 13CO molecular gas (Bally et al. 1987), and that of the stars in the sample, in the local standard of rest reference frame. The left-hand panel shows the full view of the cluster, and the right-hand panel is binned in five discrete slices along δ. The velocity distribution of the gas is shown in greyscale in the background of the left-hand panel, and as a black line in each tier of the panel on the right. Infalling and outflowing stars are shown in red and blue, respectively. On the right, each distribution is scaled relative to its peak. Typical uncertainty in RV is 0.1 km s−1.
Figure 8.

Radial velocity of the 13CO molecular gas (Bally et al. 1987), and that of the stars in the sample, in the local standard of rest reference frame. The left-hand panel shows the full view of the cluster, and the right-hand panel is binned in five discrete slices along δ. The velocity distribution of the gas is shown in greyscale in the background of the left-hand panel, and as a black line in each tier of the panel on the right. Infalling and outflowing stars are shown in red and blue, respectively. On the right, each distribution is scaled relative to its peak. Typical uncertainty in RV is 0.1 km s−1.

Table 2.

Radial velocity statistics.

InfallExpansionRotation
|$v_{med}^{a}$||$\sigma _{v}^{a}$|N*log gbvmedσvN*log gvrσvN*log g
−4.9° < δ < −4.5°13.41.4893.9711.73.0634.0913.12.0434.01
−5.3° < δ < −4.9°12.31.81113.9010.52.71553.9411.42.6773.90
−5.7° < δ < −5.3°9.61.91903.9310.03.42583.9010.02.22453.94
−6.1° < δ < −5.7°9.21.01394.128.61.6653.958.33.4404.02
−6.5° < δ < −6.1°9.60.9963.908.01.3313.937.62.1453.93
InfallExpansionRotation
|$v_{med}^{a}$||$\sigma _{v}^{a}$|N*log gbvmedσvN*log gvrσvN*log g
−4.9° < δ < −4.5°13.41.4893.9711.73.0634.0913.12.0434.01
−5.3° < δ < −4.9°12.31.81113.9010.52.71553.9411.42.6773.90
−5.7° < δ < −5.3°9.61.91903.9310.03.42583.9010.02.22453.94
−6.1° < δ < −5.7°9.21.01394.128.61.6653.958.33.4404.02
−6.5° < δ < −6.1°9.60.9963.908.01.3313.937.62.1453.93

Notes. aMean (lsr) radial velocity and radial velocity dispersion in km s−1 fitted as a Gaussian to the RV distribution in the slice, ignoring outlying wide wings to exclude likely spectroscopic binaries.

bMean log g of the stars in the slice.

Table 2.

Radial velocity statistics.

InfallExpansionRotation
|$v_{med}^{a}$||$\sigma _{v}^{a}$|N*log gbvmedσvN*log gvrσvN*log g
−4.9° < δ < −4.5°13.41.4893.9711.73.0634.0913.12.0434.01
−5.3° < δ < −4.9°12.31.81113.9010.52.71553.9411.42.6773.90
−5.7° < δ < −5.3°9.61.91903.9310.03.42583.9010.02.22453.94
−6.1° < δ < −5.7°9.21.01394.128.61.6653.958.33.4404.02
−6.5° < δ < −6.1°9.60.9963.908.01.3313.937.62.1453.93
InfallExpansionRotation
|$v_{med}^{a}$||$\sigma _{v}^{a}$|N*log gbvmedσvN*log gvrσvN*log g
−4.9° < δ < −4.5°13.41.4893.9711.73.0634.0913.12.0434.01
−5.3° < δ < −4.9°12.31.81113.9010.52.71553.9411.42.6773.90
−5.7° < δ < −5.3°9.61.91903.9310.03.42583.9010.02.22453.94
−6.1° < δ < −5.7°9.21.01394.128.61.6653.958.33.4404.02
−6.5° < δ < −6.1°9.60.9963.908.01.3313.937.62.1453.93

Notes. aMean (lsr) radial velocity and radial velocity dispersion in km s−1 fitted as a Gaussian to the RV distribution in the slice, ignoring outlying wide wings to exclude likely spectroscopic binaries.

bMean log g of the stars in the slice.

4 DISCUSSION

4.1 Distance to the ONC and receding star formation

Previously, there has been some contention regarding the orientation of the ONC’s parent filament in the plane of the sky: Kounkel et al. (2018) and Großschedl et al. (2018) found it to be largely flat, with a near-constant distance across its length, whereas Stutz, Gonzalez-Lobos & Gould (2018) and (Getman et al. 2019) found a significant variation in distance, with the central cluster being found at a larger distance than everything else in the ONC.

Separating the sources into different age bins based on their log g does show a peculiar trend in Fig. 7: the distance to the ONC does not appear to be constant as a function of age. As previously, we restrict the sample only to the sources with parallax uncertainty <0.04 mas, to maximize the resolution in distances. The median parallax with the corresponding reduced uncertainty of all the stars in each age bin is given in Table 3. It shows that the younger stars are found at increasingly larger distances. The two-sample KS test rejects the null hypothesis that the distribution of parallaxes in the oldest and in the youngest age bins originate from the same distribution with the P = 2.4 × 10−9, or ∼6σ.

Table 3.

Typical distance to the ONC as a function of age.

AgeParallaxDistance
(Myr)(mas)(pc)
1–22.5386 ± 0.0025393.92 ± 0.39
2–32.5567 ± 0.0019391.12 ± 0.29
3–52.5937 ± 0.0017385.55 ± 0.26
>52.5995 ± 0.0022384.68 ± 0.32
AgeParallaxDistance
(Myr)(mas)(pc)
1–22.5386 ± 0.0025393.92 ± 0.39
2–32.5567 ± 0.0019391.12 ± 0.29
3–52.5937 ± 0.0017385.55 ± 0.26
>52.5995 ± 0.0022384.68 ± 0.32
Table 3.

Typical distance to the ONC as a function of age.

AgeParallaxDistance
(Myr)(mas)(pc)
1–22.5386 ± 0.0025393.92 ± 0.39
2–32.5567 ± 0.0019391.12 ± 0.29
3–52.5937 ± 0.0017385.55 ± 0.26
>52.5995 ± 0.0022384.68 ± 0.32
AgeParallaxDistance
(Myr)(mas)(pc)
1–22.5386 ± 0.0025393.92 ± 0.39
2–32.5567 ± 0.0019391.12 ± 0.29
3–52.5937 ± 0.0017385.55 ± 0.26
>52.5995 ± 0.0022384.68 ± 0.32

Fig. 9 shows the interpolated distance of the stars at each age bin, and it accentuates the variance in distance further. The interpolation was performed by a neural net trained to predict distances of the stars as a function of their position in the cluster and their log g. The separation across different age bins appears to hold true across the entire filament, although given that younger stars are far more numerous in the Trapezium than elsewhere, this could result in the apparent asymmetry that Getman et al. (2019) have observed. On the other hand, if the sample selection is less sensitive to the youngest stars, this could produce a more uniform distance distribution as in Kounkel et al. (2018) and Großschedl et al. (2018).

Distance to the stars in the ONC separated in different age bins, with different lines showing the interpolated distance along the filament in each bin. The sources shown in the plot have parallax uncertainty <0.04 mas, typically ∼6 pc.
Figure 9.

Distance to the stars in the ONC separated in different age bins, with different lines showing the interpolated distance along the filament in each bin. The sources shown in the plot have parallax uncertainty <0.04 mas, typically ∼6 pc.

We note that stars with age <1 Myr may appear to break this trend of receding star formation, their average distance appears to be centred on 390 pc, closer than stars in 1–2 Myr age bin (Fig. 7). However, as they are most likely to be heavily extinguished (both due to being still embedded within the envelope of gas and due to having a protoplanetary disc, many without any optical emission), imposing strict parallax quality cuts, with σπ depending on G flux, may result in their census being particularly biased to sources sitting towards the front of the cloud.

The overall progression of distances from 5 Myr stars to 1 Myr stars itself is unlikely to be attributable to such a bias. Older stars tend to have somewhat lower luminosity (by ∼0.5 mag) due to having smaller radii, however the difference in the distance modulus between 385 and 395 pc is negligible (0.05 mag). Even considering the difference in extinction increasing with distance travelling through the cloud, it is difficult to explain older and younger stars not necessarily being co-located.

The stars (particularly those that are infalling, regardless of their age) tend to share a common radial velocity with the gas. However, if older stars are physically separated from the younger stars that are currently being formed, this implies that the bulk of the ONC is not necessarily co-located with molecular gas, with stars sitting in front of the reservoir of gas. Indeed, a 3D model of ionized gas of the nebula assumes that even the central cluster is located in front of the cloud (O’dell 2001).

This, combined with gravitational infall (which has also been previously observed by Getman et al. 2019), can in part explain the peculiar RV structure in the ONC towards NGC 1977 (δ ∼ −4.9°). The RVs towards it are not only more redshifted than what is found in the south of the cluster, the stellar RVs (of the infalling stars) are also somewhat more redshifted than the gas. Various scenarios for its motions have been considered by Getman et al. (2019). However, this is likely a signature of gravitational infall, sources in NGC 1977 being attracted to the younger stars in the central cluster that are located further away. This is similar to the scenario considered by Tobin et al. (2009) and Proszkow et al. (2009).

Similarly, Hacar et al. (2017) have observed a blueshifted ‘wedge’ in the N2H + gas in the inner 1 pc within the central cluster that is less apparent in the more diffuse CO gas. The fresh dense molecular gas that is sitting behind the cluster can none the less feel gravitational attraction to the cluster, and it is also infalling as it is pulled towards it.

4.2 Ejected stars and their effect on the observed boundness of the ONC

4.2.1 Framework

The Orion Nebula is the singularly most massive cluster that is found within more than 1 kpc of the Sun. It has long since been considered that the ONC is a cluster that may be similar to the progenitor of the Pleiades (Kroupa, Aarseth & Hurley 2001; Moraux, Kroupa & Bouvier 2004).

The comprehensive membership catalogue from McBride & Kounkel (2019) indicates that the maximum projected central stellar density is ∼2 × 105 stars pc−2, well in excess of ∼50 stars pc−2 in Pleiades (Gaia Collaboration 2018). This does not change significantly even estimating the conversion from surface density to the spatial density. While the entire ONC does extend ∼10 pc along the line of sight, the bulk of the sources towards the central cluster are young, and are concentrated at a much smaller spread of distances (hence the shape of the ONC in the works of Stutz et al. 2018 and Getman et al. 2019). The Pleiades also contains ∼3 times fewer stars than the ONC in a similar volume. The overall density of the ONC increases further if one were to consider a sizable mass of the molecular gas that is still forming stars. Over the course of its evolution, Pleiades has lost a substantial fraction of its stars (Dinnbier & Kroupa 2020), in part due to dynamical heating by binaries (Converse & Stahler 2010).

Given its mass, if ONC is not bound, no other cluster associated with a star-forming regions is expected to be – not just in the Orion Complex, but across the solar neighbourhood as a whole, within at least 1 kpc, as no other young cluster can match it in terms of its mass or its density within that volume. Since approximately 16 per cent of stars appear to form in bound clusters (Anders et al. 2021), this seems statistically unlikely.

However, while boundness of the ONC seems to be apparent conceptually, not just in comparison to other clusters, but also based on its morphology, and the effect its presence has on the surrounding stars (e.g. it appears to have substantially large self-gravity to attract other parts of the filament), there has been a long-standing difficulty in definitively proving it. Virial parameter (Bertoldi & McKee 1992) is often used to test boundness of a cluster, but while the ONC is a young cluster that is closest to being considered in a virial equilibrium (Kuhn et al. 2019), depending on the precise assumptions used to infer the mass distribution, and given the uncertainty in the total mass and radius of the cluster, it may fall just short of being considered unequivocally bound (Da Rio, Tan & Jaehnig 2014).

We propose an method that could help to mitigate this long-standing issue. The σv of the ONC may be inflated due to ejections of stars in unstable N-body interactions, i.e. incidents where three or more stars are in close proximity exchanging energy and kicking one of the stars outwards. This can occur in primordial unstable triple systems, or a close approach of a binary to another star in a dense cluster. Such ejection events can explain the apparent excess of stars moving on outward trajectories as well as a significantly wider σv of the outflowing versus infalling or rotating stars. We note that not every star that is moving on the outward trajectory would necessarily have to have been ejected; the excess can be explained by as few ∼5 per cent of stars in the total sample within 0.4° radius of the central clusters being runaways.

It may be difficult to fully disentangle which specific stars have been ejected, and which stars have proper motions that only coincidentally appear to point away from the cluster centre. However, stars that are outside of the innermost cluster core, stars that are found in the parameter space not balanced by an equal number of infalling stars, stars that project back directly to the centre without any angular offset are particularly likely candidates of bona fide orphans from a past ejection.

4.2.2 Likelihood of ejection events

To test the likelihood of the ONC having a suitable number of unstable interactions between its members, we model a simple cluster in GalPy (Bovy 2015) with AMUSE framework (Portegies Zwart et al. 2019). While this model is not complex enough to accurately reproduce such unstable interactions, in part due to a lack of a prescription for binary stars, it allows to examine how often close encounters are expected to occur, and compare their proximity to the typical orbital separations.

The model cluster consists of 2000 stars within a radius of 2 pc, overall potential of 4500 M, and velocity dispersion of ∼2.2 km s−1. The cluster evolved for 5 Myr, recording position of the stars every 100 yr. In this simulation, there are 2700 unique encounters closer than 1000 au (i.e. on average each star has experienced 1.35 such encounters). Of these, 517 encounters (26 per cent) are closer than 200 au, and even 32 encounters (1.6 per cent) closer than 30 au.

Approximately 50 per cent of the field stars are found in multiple systems. In the field, ∼50 per cent of systems have companions with separations >30 au, 33 per cent have separations >200 au, and ∼20 per cent of systems have separations >1000 au (Raghavan et al. 2010). That is to say, ∼25 per cent, 16 per cent, and 10 per cent of all stars are expected to have companions with such separations, respectively. Considering the above encounter rate, there is 10 per cent chance that there will be an encounter <1000 au in which at least one of a star has a companion >1000 au, 3 per cent chance for an encounter <200 au with >200 au companion, and 0.3 per cent chance for an encounter of an encounter <30 au with >30 au companion. Integrating across a full range of separations and considering the sizable number of stars within the ONC, this amounts to hundreds of binary stars that could have been involved in a dynamical altercation with another star in a dense cluster.

Given the frequency of encounters with other stars approaching closer than <1000 au, from this simulation it is expected that most of binaries with separations >1000 au to be largely unstable. This is supported by observations: throughout the ONC the wide binary fraction of sources with companions with separations >1000 au is only ∼5 per cent (Jerabkova et al. 2019), meaning that the orbits of at least half of wide binaries would have been processed in some degree (either hardened or became unbound). This fraction would grow higher if systems with separations <1000 au are also considered. Similarly, the estimates above include only a chance encounter between a binary and a single star; primordial triple systems would increase likelihood of unstable interactions to a degree that is difficult to quantify, but it should be noted that Class 0 YSOs have a much higher multiplicity fraction than the field stars, and that 50 per cent of multiples among them found in a high-order (3+) system in comparison to only ∼20 per cent among the field stars (Chen et al. 2013).

In some configurations, the energy exchange in the encounter would lead to unbinding of the binary. Often, one of the stars would form a binary system with an interloper/wider triple star while ejecting the original/closer companion. In other configurations, the binary would grow tighter, giving energy to an interloper/wider triple. As such, even stars for which a wide companion is detected, this does not mean that they were not participating in a dynamical interaction.

4.2.3 Properties of runaways

In dynamical simulations of triple systems, the mean ejection velocity of the unstable N-body interactions is 2.8 km s−1, which later decays to 1.1 km s−1 through intracluster interactions (Reipurth et al. 2010). This is similar to the typical velocity of the stars that are expanding away that we observe. Comparatively, true high-velocity walkaway and runaway stars are rare, though there are several dozen that are currently known to be associated with the ONC (McBride & Kounkel 2019; Farias, Tan & Eyer 2020; Schoettler et al. 2020). As such, hundreds of stars with lower ejection velocity are expected. As we have deliberately cut the amplitude of the proper motions to exclude high-velocity stars that could be considered as runaways, the bulk of these expanding stars may still remain gravitationally bound to the cluster, and they would be forced to turn around as they climb out the potential well. The typical free-fall time in the cluster is ∼0.5–1 Myr, comparable to the time it would take for a star ejected with a speed of 2 km s−1 to fully decelerate. Stars with ejection speed of >5 km s−1 are likely to be unbound. The primary reason why so many could be detected can be attributed to the youth of these stars and the recency of their ejection.

Almost all of the expanding stars (∼94 per cent) have traceback ages that are smaller than the age assigned to a star based on their log g. This is to be expected, as it would be impossible to eject them from the central cluster otherwise. Interestingly, the candidate ejected stars in the younger age bins appear to have a larger high-velocity (>4 km s−1) tail in comparison to their older counterparts (Fig. 2, third column). As the catalogue on which this analysis is performed consists of sources that have been targeted for spectroscopic observations, the older stars moving with comparable speeds are likely to be out of bounds of the targeted area, assuming they were ejected shortly after their formation.

4.2.4 Implications

The virial parameter is a somewhat imperfect metric of a cluster boundness, as it is difficult to measure directly, without any assumptions pertaining to e.g. cluster morphology and mass distribution within it. Moreover, it can be overestimated due to the binary stars. For example, σv can be inflated by the presence of spectroscopic binaries (Kouwenhoven & de Grijs 2008; Gieles, Sana & Portegies Zwart 2010). In the more evolved clusters boundness of which is not in question, mass of a cluster estimated from the virial parameter may be inflated by as much as 8–42 per cent (Rastello, Carraro & Capuzzo-Dolcetta 2020).

Here, we also introduce another consideration which can inflate the velocity dispersion, and thus affect the virial parameter, namely unstable interactions with multiple systems exchanging orbital energy, resulting in ejection of the individual stars. These individual stars may or may not become unbound depending on their velocity kick they receive relative to the escape velocity of the cluster as a whole. Such sources are able to contribute to the wider wings of σv (Fig. 8), but they would not affect boundness of the stars that did not experience such interactions, since, at the moment, the energy injection is applied to the individual stars, not the entire cluster as a whole.

Eventually the compounding effect of such ejections will propagate to the other stars. If an ejected star is unbound, it would lead to an incremental mass-loss for a cluster. If a star would remain bound, it would eventually be able to equilibrate its velocity with other stars, incrementally heating up the cluster. The final σv of the cluster will be higher than the initial σv that were not involved in the aforementioned dynamical interactions, but as all of the excess energy would be shared equally, it would be lower than σv in the sample in which the ejected stars are dominating the wings of the velocity dispersion.

But, after some time, a cluster would be able to ‘settle’: all the wide binaries that could be disrupted from their primordial separations would already do so, the hard binaries would increasingly become harder. Thus, eventually a close encounter with another star would be less likely to result in an ejection event in more evolved clusters compared to those that are still near their infancy.

Such a state of equilibrium has not yet been achieved in ONC, as can be inferred by the presence of a large number of runaways with an ejection age <1 Myr ago (McBride & Kounkel 2019), including a rather spectacular event mere ∼540 yr ago that has led to the formation of the BN/KL nebula and single-handedly ejected at least four stars (Gómez et al. 2008; Luhman et al. 2017; Rodríguez et al. 2020). Arguably, such a state of equilibrium cannot be achieved in a young cluster that is still actively forming stars, and thus is still actively producing binary stars at the full range of their orbital separations.

Orbital evolution of binary systems have often been considered to play a vital role in characterizing the dynamics in young populations as a whole, and young dense clusters such as the ONC in particular (Kroupa & Burkert 2001; Clarke, Mathieu & Reid 2015). However, there is still much regarding young multiples that is not well understood – in part because the observational data towards them are still sparse, in part because multiple systems add a significant degree of complexity that is often disregarded in theoretical framework (e.g. even in simulations that specifically aim to characterize ejected runaways from young clusters Schoettler et al. 2019). This also holds true for the studies that have so far applied the virial parameter to test boundness of young clusters, including the ONC (e.g. Da Rio et al. 2014; Kuhn et al. 2014). More detailed numerical simulations of ONC-like clusters are needed to constrain the role of multiples in general and ejected stars specifically to test the manner in which they challenge the typical approximations of virial equilibrium.

4.3 Combined model of the star-forming history of the ONC

It has been previously noted by Kounkel, Covey & Stassun (2020) and Großschedl et al. (2021) that 6 Myr ago a supernova (or several supernovae) have triggered the global expansion of the Orion Complex. Given that the ONC is being pushed into the Orion A molecular cloud (Großschedl et al. 2018) from the direction of the centre of the expansion (several pc away from the current position of the ONC), it is likely that the shockwave from a supernova has swept along the gas through the filamentary cloud, compressing the gas, and jump-starting the formation of the cluster in a ‘conveyor belt’ like fashion (Krumholz & McKee 2020).

This shock-driven star formation scenario can account for the variable distance to the ONC as a function of age. If the shockwave is propagating into the cloud at a slightly faster rate (faster by ∼2 km s−1) than the typical radial velocity of the gas and the stars, this would displace the star-forming front. Such propagation velocity is comparable to the sound speed in the ONC (Goicoechea et al. 2015).

Early on, star formation has occurred all throughout the length of the filament at an equal rate, but, eventually, as the cluster grew more and more massive, self-gravity became more important. The molecular gas was pulled increasingly more towards the middle, forming the central cluster, even though the bulk of the stars in the central cluster were formed further away than the initial burst of star formation in the region. As the star-forming front of the ONC was being pushed into the Orion A molecular cloud (Fig. 10) and continuing to access fresh gas, this provided sufficient fuel to form multiple generations of stars in the central cluster, increasing its mass and gravitational pull. On the other hand, NGC 1980 to the south and NGC 1977 to the north of it could not sustain star formation beyond the initial burst. To various degree, the molecular gas that was originally there was (a) consumed in forming stars, (b) attracted towards the central cluster, and (c) dissipated through stellar feedback. They could not replenish their gas content with Trapezium hoarding all of the new material.

A conceptual model showing a side view of the Orion A molecular cloud (corresponding to distance versus δ projection), and the relation of ONC relative to it. The cluster shows the age gradient with distance, with the younger stars located closer to the cloud; direction of star formation propagation is indicated by a yellow arrow. Meanwhile, the northern and southern parts of the cluster are contracting towards the central cluster, indicated by the black arrows.
Figure 10.

A conceptual model showing a side view of the Orion A molecular cloud (corresponding to distance versus δ projection), and the relation of ONC relative to it. The cluster shows the age gradient with distance, with the younger stars located closer to the cloud; direction of star formation propagation is indicated by a yellow arrow. Meanwhile, the northern and southern parts of the cluster are contracting towards the central cluster, indicated by the black arrows.

Getman et al. (2019) have noted that the more distant stars in the ONC are receding from the observer at a slower rate than the closer parts of the cluster. This does appear to be the case of the entire population, particularly due to NGC 1977. However, this distance–velocity relation is not immediately apparent in the RV distributions of stars segregated into different age bins along a given line of sight. None the less, this may be a signature of either the shockwave slowing down as it is encountering more and more gas of the filament, or the self-gravity of the ONC is protesting against the sheering of the spatially differential star formation and is attempting to bring the older nearby stars and younger distant stars closer together.

4.4 Order of magnitude mass estimate of the ONC

The sources that are infalling may offer a possible way estimate the dynamical mass of both the ONC and the Orion Complex; however, as has been previously stated, such a calculation is difficult as we do not fully know the initial velocity field that was present in the region beforehand. Similarly, as acceleration changes velocity incrementally over time, it is important to consider the time-scales over which the force of gravity affects the surroundings, which is non-trivial to estimate, as the conditions in a star-forming cloud can and do drastically change over time.

However, it is possible to make a rough order of magnitude inference assuming simple conditions to compare the derived mass with what is typically assumed for the population in question.

If we treat ONC as an isolated system and if we assume typical infall velocities of 2 km s−1 at 5 pc, which is a rough order of magnitude of the observed velocities, starting from rest, over ∼3 Myr, it would require a cluster mass of 2500 M to achieve such an acceleration. This is well comparable to the estimate of mass of the ONC from Hillenbrand & Hartmann (1998) of 4500 M.

5 CONCLUSIONS

We analyse the dynamics of young stars within the Orion Nebula and the solar neighbourhood in the vicinity of the Orion Complex as a whole. We find that

  • Examining the orientation of the proper motions of Orion members, we detect a significant (>4σ) excess of sources whose proper motions are consistent with infall towards the ONC. These infall signatures are most prominent to the north and south of the ONC, within the integral filament where stars are least likely to have been dynamically processed by a passage through the ONC. We interpret this signature as evidence that the stars in and around the ONC are contracting lengthwise along the integral filament, likely due to the self-gravity of the Trapezium.

  • We also detect a more modest (∼3σ) excess of sources with proper motions consistent with expansion from the ONC; the significance of this excess is largest for sources offset in RA, rather than Dec., from the ONC, placing them outside the integral filament, and thus more likely to be on orbits that have been dynamically processed via interactions in the ONC. We interpret this signature as evidence of dynamical processing of the ONC’s stellar population, not due to the cluster itself being unbound. However, dynamical ejections could inflate the estimate of the virial parameter of a young population, which is a common metric that tests for boundness.

  • Among the youngest stars in the ONC that have purely tangential motion, there is a preferred direction of rotation around the central cluster, there is no such preference among the older stars. As such, the angular momentum of the cluster has either evolved to develop organized rotation after the cluster has accreted enough mass, or early signatures of organized rotation in older stars has since been washed out through the dynamical evolution.

  • The distance to the ONC depends on the age of the sample; it varies from ∼385 pc for the oldest stars, to ∼395 pc for the younger stars; the star formation is continuously propagating into the Orion A molecular cloud at a faster rate than the typical velocity of the gas and the stars, consuming the outer layers of the gas in the process.

  • These observations reveal that young star-forming clusters are highly dynamic entities, and the structure of these clusters as well as velocities of stars within them are affected by various processes occurring across different scales. The ONC is not homogeneous. Its dynamical evolution is influenced in different measure by the gravitational potential of its surroundings, gravitational interactions between the stars, and the activity within the rest of the Complex. All these effects leave different kinematical signatures that are difficult to separate when examining the entire cluster as a whole, but they become more apparent through identifying appropriate subsets of stars.

DATA AVAILABILITY

This work is based on publicly available data. APOGEE Net reduction of APOGEE spectra are available at http://vizier.u-strasbg.fr/viz-bin/VizieR?-source = J/AJ/163/152, cross-matched with the catalogue of members of the ONC http://vizier.u-strasbg.fr/viz-bin/VizieR?-source = J/ApJ/884/6 and with Gaia DR3 photometry and astrometry. No additional data products are computed in this work.

Footnotes

1

LSR velocity correction is from (Schönrich, Binney & Dehnen 2010).

REFERENCES

Anders
F.
,
Cantat-Gaudin
T.
,
Quadrino-Lodoso
I.
,
Gieles
M.
,
Jordi
C.
,
Castro-Ginard
A.
,
Balaguer-Núñez
L.
,
2021
,
A&A
,
645
,
L2

Bally
J.
,
Langer
W. D.
,
Stark
A. A.
,
Wilson
R. W.
,
1987
,
ApJ
,
312
,
L45

Beccari
G.
et al. ,
2017
,
A&A
,
604
,
A22

Bertoldi
F.
,
McKee
C. F.
,
1992
,
ApJ
,
395
,
140

Bovy
J.
,
2015
,
ApJS
,
216
,
29

Chen
X.
et al. ,
2013
,
ApJ
,
768
,
110

Clarke
C. J.
,
Mathieu
R. D.
,
Reid
I.
,
2015
,
Dynamics of Young Star Clusters and Associations, 1 edn. Saas-Fee Advanced Course
.
Springer
,
Berlin

Converse
J. M.
,
Stahler
S. W.
,
2010
,
MNRAS
,
405
,
666

Da Rio
N.
et al. ,
2017
,
ApJ
,
845
,
105

Da Rio
N.
,
Tan
J. C.
,
Jaehnig
K.
,
2014
,
ApJ
,
795
,
55

David
T. J.
,
Hillenbrand
L. A.
,
Gillen
E.
,
Cody
A. M.
,
Howell
S. B.
,
Isaacson
H. T.
,
Livingston
J. H.
,
2019
,
ApJ
,
872
,
161

Dinnbier
F.
,
Kroupa
P.
,
2020
,
A&A
,
640
,
A85

Dzib
S. A.
et al. ,
2017
,
ApJ
,
834
,
139

Farias
J. P.
,
Tan
J. C.
,
Eyer
L.
,
2020
,
ApJ
,
900
,
14

Fűrész
G.
,
Hartmann
L. W.
,
Megeath
S. T.
,
Szentgyorgyi
A. H.
,
Hamden
E. T.
,
2008
,
ApJ
,
676
,
1109

Gaia Collaboration
,
2018
,
A&A
,
616
,
A10

Gaia Collaboration
,
2021
,
A&A
,
649
,
A1

Getman
K. V.
,
Feigelson
E. D.
,
Kuhn
M. A.
,
Garmire
G. P.
,
2019
,
MNRAS
,
487
,
2977

Gieles
M.
,
Sana
H.
,
Portegies Zwart
S. F.
,
2010
,
MNRAS
,
402
,
1750

Goicoechea
J. R.
et al. ,
2015
,
ApJ
,
812
,
75

Gómez
L.
,
Rodríguez
L. F.
,
Loinard
L.
,
Lizano
S.
,
Allen
C.
,
Poveda
A.
,
Menten
K. M.
,
2008
,
ApJ
,
685
,
333

Großschedl
J. E.
et al. ,
2018
,
A&A
,
619
,
A106

Großschedl
J. E.
,
Alves
J.
,
Meingast
S.
,
Herbst-Kiss
G.
,
2021
,
A&A
,
647
,
A91

Hacar
A.
,
Alves
J.
,
Tafalla
M.
,
Goicoechea
J. R.
,
2017
,
A&A
,
602
,
L2

Hartmann
L.
,
Burkert
A.
,
2007
,
ApJ
,
654
,
988

Hillenbrand
L. A.
,
Hartmann
L. W.
,
1998
,
ApJ
,
492
,
540

Hsieh
C.-H.
,
Arce
H. G.
,
Mardones
D.
,
Kong
S.
,
Plunkett
A.
,
2021
,
ApJ
,
908
,
92

Jerabkova
T.
,
Beccari
G.
,
Boffin
H. M. J.
,
Petr-Gotzens
M. G.
,
Manara
C. F.
,
Prada Moroni
P. G.
,
Tognelli
E.
,
Degl’Innocenti
S.
,
2019
,
A&A
,
627
,
A57

Jones
B. F.
,
Walker
M. F.
,
1988
,
AJ
,
95
,
1755

Kounkel
M.
et al. ,
2018
,
AJ
,
156
,
84

Kounkel
M.
et al. ,
2021
,
AJ
,
162
,
184

Kounkel
M.
,
Hartmann
L.
,
Tobin
J. J.
,
Mateo
M.
,
Bailey
J. I. III
,
Spencer
M.
,
2016
,
ApJ
,
821
,
8

Kounkel
M.
,
Covey
K.
,
Stassun
K. G.
,
2020
,
AJ
,
160
,
279

Kouwenhoven
M. B. N.
,
de Grijs
R.
,
2008
,
A&A
,
480
,
103

Kroupa
P.
,
Burkert
A.
,
2001
,
ApJ
,
555
,
945

Kroupa
P.
,
Aarseth
S.
,
Hurley
J.
,
2001
,
MNRAS
,
321
,
699

Krumholz
M. R.
,
McKee
C. F.
,
2020
,
MNRAS
,
494
,
624

Kuhn
M. A.
et al. ,
2014
,
ApJ
,
787
,
107

Kuhn
M. A.
,
Hillenbrand
L. A.
,
Sills
A.
,
Feigelson
E. D.
,
Getman
K. V.
,
2019
,
ApJ
,
870
,
32

Kuznetsova
A.
,
Hartmann
L.
,
Ballesteros-Paredes
J.
,
2015
,
ApJ
,
815
,
27

Luhman
K. L.
,
Robberto
M.
,
Tan
J. C.
,
Andersen
M.
,
Giulia Ubeira Gabellini
M.
,
Manara
C. F.
,
Platais
I.
,
Ubeda
L.
,
2017
,
ApJ
,
838
,
L3

Marigo
P.
et al. ,
2017
,
ApJ
,
835
,
77

McBride
A.
,
Kounkel
M.
,
2019
,
ApJ
,
884
,
6

McBride
A.
,
Lingg
R.
,
Kounkel
M.
,
Covey
K.
,
Hutchinson
B.
,
2021
,
AJ
,
162
,
282

Megeath
S. T.
et al. ,
2016
,
AJ
,
151
,
5

Moraux
E.
,
Kroupa
P.
,
Bouvier
J.
,
2004
,
A&A
,
426
,
75

O’dell
C. R.
,
2001
,
ARA&A
,
39
,
99

Olney
R.
et al. ,
2020
,
AJ
,
159
,
182

Portegies Zwart
S.
et al. ,
2019
,
AMUSE: the Astrophysical Multipurpose Software Environment

Proszkow
E.-M.
,
Adams
F. C.
,
Hartmann
L. W.
,
Tobin
J. J.
,
2009
,
ApJ
,
697
,
1020

Raghavan
D.
et al. ,
2010
,
ApJS
,
190
,
1

Rastello
S.
,
Carraro
G.
,
Capuzzo-Dolcetta
R.
,
2020
,
ApJ
,
896
,
152

Reipurth
B.
,
Mikkola
S.
,
Connelley
M.
,
Valtonen
M.
,
2010
,
ApJ
,
725
,
L56

Rodríguez
L. F.
,
Dzib
S. A.
,
Zapata
L.
,
Lizano
S.
,
Loinard
L.
,
Menten
K. M.
,
Gómez
L.
,
2020
,
ApJ
,
892
,
82

Schoettler
C.
,
Parker
R. J.
,
Arnold
B.
,
Grimmett
L. P.
,
de Bruijne
J.
,
Wright
N. J.
,
2019
,
MNRAS
,
487
,
4615

Schoettler
C.
,
de Bruijne
J.
,
Vaher
E.
,
Parker
R. J.
,
2020
,
MNRAS
,
495
,
3104

Schönrich
R.
,
Binney
J.
,
Dehnen
W.
,
2010
,
MNRAS
,
403
,
1829

Sprague
D.
et al. ,
2022
,
AJ
,
163
,
152

Stutz
A. M.
,
Gonzalez-Lobos
V.
,
Gould
A.
,
2018
,
preprint (arXiv:1807.11496)

Theissen
C. A.
,
Konopacky
Q. M.
,
Lu
J. R.
,
Kim
D.
,
Zhang
S. Y.
,
Hsu
C.-C.
,
Chu
L.
,
Wei
L.
,
2022
,
ApJ
,
926
,
141

Tobin
J. J.
,
Hartmann
L.
,
Furesz
G.
,
Mateo
M.
,
Megeath
S. T.
,
2009
,
ApJ
,
697
,
1103

This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://dbpia.nl.go.kr/journals/pages/open_access/funder_policies/chorus/standard_publication_model)