On the true shape of the upper end of the stellar initial mass function
Key Words.:Galaxies: star clusters: general – Methods: numerical – Open clusters and associations: individual (R136) – Stars: kinematics and dynamics – Stars: luminosity function, mass function – Stars: massive
Context:The shape of the stellar initial mass function (IMF) of a star cluster near its upper mass limit is a focal topic of investigation as it determines the high mass stellar content and hence the dynamics of the cluster at its embedded phase as well as during its young gas-free phase. The massive stellar content of a young cluster, however, can be substantially modified due to the dynamical ejections of the massive stars so that the present-day high-mass stellar mass function (hereafter MF) can be different than that with which the cluster is born.
Aims:In the present study, we provide a preliminary estimate of this evolution of the high-mass IMF of a young cluster due to early ejections of massive stars, using the Large Magellanic Cloud massive, young cluster R136 as an example.
Methods:To that end, we utilize the results of the state-of-the-art calculations by Banerjee, Kroupa & Oh (2012) comprising direct N-body computations of realistic, binary-rich, mass-segregated models of R136. In particular, these calculations provide the ejection fraction of stars as a function of stellar mass.
Results:We find that if the measured IMF of R136 is granted to be canonical, as observations indicate, then the “true” high-mass IMF of R136 at its birth must be at least moderately top-heavy when corrected for the dynamical escape of massive stars.
Conclusions:The top-heaviness of the true high-mass IMF over the observationally determined one is a general feature of massive, young clusters where the dynamical ejection of massive stars is efficient. We discuss its implications and possible improvements over our current estimate.
It is widely believed that stars form following an initial mass function (IMF) that has a universal form (Bastian et al., 2010; Kroupa et al., 2012). In course of time however, the mass function (MF) of a stellar population evolves as the stellar members are removed from the system. For bound stellar clusters, the low mass stars are naturally vulnerable to removal since they evaporate out of the system as the cluster relaxes (Spitzer, 1987; Vesperini & Heggie, 1997; Baumgardt & Makino, 2003). The depletion of the low-mass stars from the bound system can be greatly augmented due to the expulsion of the residual gas that did not form stars, in the early phase of the cluster evolution (Baumgardt & Kroupa, 2007; Marks & Kroupa, 2010). A marked positive correlation of the global low-mass MF index () with the cluster’s concentration has been pointed out for the first time by De Marchi et al. (2007). Baumgardt et al. (2008) showed that such a correlation can be explained by clusters born with a canonical IMF which are primordially mass-segregated but are filling their tidal radii. Marks & Kroupa (2010), however, demonstrated that the De Marchi trend can as well be explained by residual gas expulsion from initially compact and mass-segregated clusters without introducing a variable low-mass IMF. This result is consistent with the compactness of freshly formed clusters (Marks & Kroupa, 2012), whereas the tidal-radius-filling ansatz would require the clusters to condense from a cloud with large radius matching the tidal radius.
On the other hand, the distribution of the most massive stars bound to a cluster is influenced by the ejection of the massive members via dynamical encounters (Clarke & Pringle, 1992). Massive stars in a young cluster, which are usually in binaries, segregate in the cluster’s center, where their density can be greatly augmented. The massive stars therefore frequently interact through super-elastic binary-single and binary-binary encounters to eject each other from the cluster, that way depleting the massive stellar population as the cluster evolves. In general, more massive stars are more centrally concentrated and hence are more severely depleted (Pflamm-Altenburg & Kroupa, 2006). This is found to be true in theoretical studies of massive runaway stars from star clusters by several authors (Fujii & Portegies Zwart, 2011; Banerjee et al., 2012a).
The increasing depletion of (bound) massive stars with increasing stellar mass in star clusters has remarkable implications on its IMF. If a young cluster is found with the massive end of its MF in the canonical form (Kroupa, 2001) at its present age, then the upper end of its IMF has to be top-heavy to achieve the canonical form at the current age (taking into account the mass loss due to stellar evolution; see Sec. 2). In general, given a measured MF of a young cluster near the massive limit, the corresponding stellar mass distribution at birth, i.e., the actual IMF, has to be “top-heavier” depending on the age of the cluster due to stellar-dynamical ejections.
In the present work we provide a preliminary quantification of the top-heaviness of the IMF of a star cluster as inferred from the above considerations. To that end we utilize the ejected fraction of massive stars, , as a function of the mass of the ejected stars from model clusters as computed by Banerjee et al. (2012a; hereafter Paper I). The above quantity is defined as the ratio of the number of stars, within a bin around mass , moving away from the cluster beyond pc from the cluster’s center to the total number of stars in the whole system, i.e., including all the bound and the ejected members, within the same mass bin. Note that was denoted by in Paper I but we stick to this new symbol in this paper for better clarity. The above authors studied the runaway O-stars from the Large Megallanic Cloud (LMC) cluster R136 by computing the dynamical evolution of realistically modelled R136-type clusters using direct N-body integration. We point out the essentials of these computed models of Paper I in Sec. 3.
2 Dynamical evolution of the bound stellar mass function of a star cluster
In reality, the present-day observationally determined bound stellar mass function of a young star cluster is corrected for unresolved multiple companions (Kroupa, Gilmore & Tout, 1991; Weidner, Kroupa, & Maschberger, 2009), and for the wind mass loss of the stars to determine the so-called cluster IMF, . The correction for unresolved multiple components is negligible for massive stars (Weidner, Kroupa, & Maschberger, 2009). As pointed out in Sec. 1, needs to be further corrected for the stars removed from the bound cluster to determine the true IMF of the cluster, , which is the bound stellar mass function of the cluster at its birth. To correct over the entire mass range, one should, in principle, reassemble the cluster including all the stars that escaped from the cluster. If the cluster is formed mass-segregated (as considered in Paper I and inferred to be true for several young clusters; see, e.g., Bonnell & Davies 1998), then the massive stars reside in the cluster’s center from the very beginning of its secular evolution and they can only be ejected by dynamical encounters (i.e., they cannot be removed by the residual gas expulsion or the external galactic field). Therefore, can be obtained from by only correcting for the runaway stars, for the massive end.
In Paper I, the ejection fraction is found to increase significantly with stellar mass for at the cluster age of Myr (c.f. Fig. 4 of Paper I). The number of ejected stars over the mass range to is,
Here, the decrement of the zero-age mass, , is not included since in evaluating , the stellar masses at age are used for both the total population between to and the corresponding ejected population so that the same mass shift () is applicable to all the three functions , and . For brevity, we also suppress the explicit dependence in Eqn. (1).
The number of stars in the mass range to remaining bound to the cluster is therefore,
where we have used Eqn. (1). The right-hand-side of Eqn. (2) is the “observed” (high-mass) IMF of the cluster without including the ejected stars (but corrected for stellar evolution). In other words,
Eqn. (4) connects the observed cluster IMF, , with its true IMF, .
3 Previous computations
State-of-the-art calculations have been performed in Paper I to evolve model star clusters that mimic R136 in terms of observable parameters. The primary objective of these calculations, which are so far the largest sized ones with fully mass-segregated hard massive binaries (see below), was to study the ejected OB-stars from R136, in particular, whether massive runaway stars like VFTS 682 and 30 Dor 016 can indeed be ejected dynamically from R136 as suspected (Bestenlehner et al., 2011; Evans et al., 2010). These computations demonstrated good agreement of the kinematics of the massive ejecta with those of these noted runaways. Furthermore, they suggest that the tight massive binaries that drive these runaways can merge due to the dynamical interactions leading to the formation of “super-canonical” (Kroupa et al., 2012) stars, i.e., stars with masses exceeding the canonical upper limit of (see Paper I and Banerjee et al. 2012b), in agreement with the observed stellar population in R136 (Crowther et al., 2010). In this study, we continue to utilize the results of these computations. The initial conditions and the method of these calculations are accounted in detail in Sec. 2 of Paper I which we briefly recapitulate here.
The above computations comprise direct N-body integrations of Plummer spheres (Heggie & Hut, 2003; Kroupa, 2008) with parameters conforming with those of R136. The initial mass of the Plummer spheres is which is an upper limit for R136 (Crowther et al., 2010) and the half mass radii are taken to be pc. The clusters are made of stars with ZAMS masses drawn from a canonical IMF (Kroupa, 2001) over the range and are of metallicity . With the mean stellar mass of , each cluster comprises stars.
As for the primordial binary population, stars with are all in binaries while all lighter stars are kept single. This termination of the binary population is due to computational ease; binaries bottleneck the calculation speed of direct N-body integration significantly so that adapting a full spectrum of primordial binaries (i.e., 100% initial binary fraction; Marks & Kroupa 2012) in models of the size that we compute becomes prohibitive. Disregarding the binary population for of course does not affect the ejection of massive stars and mergers of massive binaries significantly (see Paper I for details). Following the observed period distribution of O-star binaries (Sana & Evans, 2011), the orbital periods of the binaries having primary masses are chosen from a uniform distribution between , where is the orbital period in days. For primary masses , the “birth period distribution” of Kroupa (1995) is adapted. The eccentricity distribution is chosen to be thermal (Spitzer, 1987; Kroupa, 2008). These distributions, however, are not modified via “eigenevolution” (Kroupa, 1995) as the interactions between massive pre-main-sequence stars are not yet quantified. The initial configurations are fully mass-segregated (Baumgardt et al., 2008) so that the massive binaries concentrate within the clusters’ central region mimicking primordial mass segregation which is inferred for several Galactic globular clusters and young star clusters.
To perform the direct N-body computations, the state-of-the-art “NBODY6” integrator (Aarseth, 2003) has been used. In addition to integrating the particle orbits using the highly accurate fourth-order Hermite scheme and dealing with the diverging gravitational forces in close encounters through regularizations, NBODY6 also employs the analytical stellar and binary evolution recipes of Hurley et al. (2000, 2002). The models are computed for Myr four times, each with a different random number seed initialization.
4 Ejected fraction of massive stars : a top-heavy IMF in R136
Figure 4 of Paper I shows the ejected fraction of stars as a function of stellar mass , where a significant and monotonic increase of can be seen beyond . As explained in Paper I, this trend is a direct consequence of mass segregation; more massive systems are more concentrated towards the cluster’s center and hence interact more efficiently (see below). Fig. 1 re-plots the at the cluster age Myr, as obtained from the computations of Paper I, for where we use the same symbols as in Figure 4 of Paper I to differentiate among the data from different computed models.
A possible form of the fitting function to the data in Fig. 1 would be asymptotic to 1 with increasing as by definition but monotonically increasing as the trend of the data suggests. Due to the scatter of the points in Fig. 1 over a large range and to ensure reliability, we do a least-squares fit to the logarithm of the data in Fig. 1. We adopt the exponential functional form,
to best-fit with the logarithm of the data-coordinates of Fig. 1. This functional form asymptotically increases to 0, i.e., approaches 1 when transformed into the linear scale. We use all the 28 data-points of the above figure to obtain the best fit values and . The central solid line in Fig. 1 represents the best-fit line to the data and the two dashed lines represent the upper and the lower bounds of the fit corresponding to the errors (1- or the standard errors) in and .
As such, the data-points in Fig. 1 can be found to follow an overall power-law index of if one fits a power-law to them. A basic physical argument that supports such an increment is as follows. The massive single stars get ejected by super-elastically (Heggie & Hut, 2003) encountering with binaries within the cluster (as noted in Paper I, most of the ejected entities are single stars, given the initial distribution of the binaries’ orbital periods adapted here111 Note that ejected single stars may in reality be binaries if the models would incorporate a large fraction of O-stars being in very tight ( few au) binaries within wider triple or quadruple systems with outer orbital periods as chosen here (see Chini et al. 2012).). These single stars are liberated from binaries either by binary-binary encounters or by direct disruption of binaries and they can, in general, be expected to remain mass-segregated within the cluster. A given segregated group of stars of mass remain confined within a radius of the cluster (see Chapter 16 of Heggie & Hut 2003) so that their density . This density is also proportional to the single-stars’ encounter rate with the binaries and hence to the probability of the formers’ ejection from the cluster. Hence, one can expect that the ejection fraction which is indeed close to the best-fitted power-law index as obtained here.
The dotted line in Fig. 1 represents the best-fit power law (truncated at ). By doing a linear least-squares fit to the logarithm of the data-points, we obtain the best-fit values and . Unfortunately, due to the sparsity of data-points near the massive end it is presently difficult to determine the true functional behavior of over the entire data of Fig. 1. A larger number of such computations would help to settle this issue.
The true IMFs, , corresponding to the lines in Fig. 1, as obtained from Eqn. (4), is shown in Fig. 2 where the observed IMF is taken to be the canonical one, i.e., following the Salpeter form for (shown by the dash-dotted line). The upper dashed, central solid, lower dashed and the dotted lines in Fig. 2 correspond to those in Fig. 1 respectively. Noticeably, the IMF of R136 is indeed observed to be very close to canonical over both low-mass (Andersen et al. 2009; ) and high-mass (Massey & Hunter 1998; ) ranges.
It can be seen in Fig. 2 that the true IMF, , is top-heavy near its massive end even for . In other words, the computed ejection fraction of massive stars from R136-type model clusters (in Paper I) indicates from a moderate to a substantial top-heaviness in the high-mass IMF of R136. This point becomes more vivid in Fig. 3 where we plot the ratio of the true IMF, , to the Salpeter law .
Notably, the fitted bounds to the data-points in Fig. 1 (the dashed lines) represent the effect of statistical variation of the ejection-fraction, , between the computed cluster models. Such cluster to cluster variation affects the top-heaviness of the inferred true IMF as can be seen in the corresponding bounds of in Fig. 2. The inherent statistical fluctuations can result in the inference of a marginally top-heavy to a moderately top-heavy true IMF (c.f. Fig. 3). In other words, if the present-day observed IMF in young clusters is canonical, then there has to be a cluster to cluster stochastic variation of the top-heaviness of the true IMF with which they are born. Further detailed study (c.f. Sec. 5) is needed to quantify this better. Dib et al. (2010) has studied the growth of massive stellar cores by gas accretion and the subsequent formation of stars including the effects of stellar wind feedback on the protocluster which indicates cluster-to-cluster variation of the stellar IMF (see also Dib et al. 2007).
5 Discussions and outlook
The above inference of top-heaviness near the upper end of the “true” IMF is not limited to R136 but is true for any massive cluster that is observed to have a normal (i.e., canonical) IMF. The extent of the top-heaviness of the true IMF for a particular cluster, of course, depends on its ejection fraction . It is to be noted that the above analysis is only of the lowest order and is essentially a qualitative inference. This is because had we started the N-body computations with the inferred top-heavy true IMF as corrected for the runaways, would have been larger than the presently used one which is obtained from a canonical IMF. Therefore, in general, the (top-heavier) true massive end of the IMF can be determined quantitatively from an observed IMF by iteratively applying the correction in Sec. 4 for the dynamical ejection of the massive stars which would tend to augment the top-heaviness. Such a study would reveal relationships between the top-heaviness of a cluster’s IMF with its initial parameters such as the total mass, the degree of primordial mass segregation and the concentration subject to the constraint that a nearly canonical mass distribution is reached at the present day which is true for many Galactic young clusters. Such relationships would, in turn, have fundamental implications of massive star formation in a clustered environment (Bonnell et al., 1998; Zinnecker & Yorke, 2007).
It is important to note that the top-heaviness that has been inferred above for the massive (true) IMF is influenced by the adopted initial complete mass-segregation (see Sec. 3) that provides the maximum efficiency in ejecting massive stars at young ages. However, such an assumption of primordial mass-segregation is in accordance with the observations of young star clusters; see Paper I and references therein (and, e.g., Bonnell & Davies 1998). As pointed out above, an interesting outlook would be to study the degree of top-heaviness of the true IMF as a function of the initial degree of segregation. The timescale of segregation of the massive stars/binaries is also shortened by the cluster’s concentration. Therefore, it would also be worthwhile to study the dependence of the effect on the cluster’s initial concentration. A physical basis of the true IMF being highly top-heavy can be found in Dib et al. (2010). In this study of the appearance of the stellar IMF via accretion onto proto-stellar cores including wind feedback, these authors do find a substantial flattening over the massive end of the IMF.
In Fig. 1, the best-fit plot extends beyond the canonical upper limit of . As explained in Paper I and Banerjee et al. (2012b), this is to take into account the “super-canonical” stars with zero-age masses between that form via mergers of massive binaries due to dynamical encounters. Although the corrected IMF, , is plotted in Fig. 2 for consistently with Fig. 1, only for is relevant as the cluster’s IMF.
In conclusion, our above analysis with the example of R136 indicates that taking into account the dynamical ejection of massive stars, a young cluster with a seemingly normal or canonical IMF can actually have formed with an IMF which is at least moderately top-heavy near its high-mass end. This may well constitute the first ever detection of the flattening of the IMF as predicted by Dib et al. (2010). For completeness, we mention that by tracing the present-day constitution of globular clusters to their initial configurations, a trend with birth density and the power-law index of the IMF has emerged such that the IMF is increasingly top-heavy with increasing star-forming density above about on a cluster-forming spatial scale (Marks et al., 2012). The present results thus add to the increasingly solid evidence that the IMF becomes top-heavy in star bursts.
- Aarseth (2003) Aarseth, S.J. 2003, “Gravitational N-Body Simulations”, Cambridge University Press.
- Andersen et al. (2009) Andersen, M. et al. 2009, ApJ, 707, 1347.
- Banerjee et al. (2012a) Banerjee, S., Kroupa, P. & Oh, S. 2012a, ApJ, 746, 15.
- Banerjee et al. (2012b) Banerjee, S., Kroupa, P. & Oh, S. 2012b, MNRAS, arXiv:1208.0826 (preprint).
- Bastian et al. (2010) Bastian, N., Covey, K.R. & Meyer, M.R. 2010, ARA&A, 48, 339.
- Baumgardt & Makino (2003) Baumgardt, H. & Makino, J. 2003, MNRAS, 340, 227.
- Baumgardt & Kroupa (2007) Baumgardt, H. & Kroupa, P. 2007, MNRAS, 380, 1589.
- Baumgardt et al. (2008) Baumgardt, H., De Marchi, G. & Kroupa, P. 2008, ApJ, 685, 247.
- Bestenlehner et al. (2011) Bestenlehner, J.M., Vink, J.S., Gräfener, G., et al. 2011, A&A, 530, L14.
- Bonnell et al. (1998) Bonnell, I.A., Bate, M.R. & Zinnecker, H. 1998, MNRAS, 298, 93.
- Bonnell & Davies (1998) Bonnell, I.A. & Davies, M.B. 1998, MNRAS, 295, 691.
- Chini et al. (2012) Chini, R., Hoffmeister, V., Nasseri, A., et al. 2012, MNRAS, 424, 1925.
- Clarke & Pringle (1992) Clarke, C.J. & Pringle, J.E. 1992, MNRAS, 255, 423.
- Crowther et al. (2010) Crowther, P.A., Schnurr, O., Hirschi, R., et al. 2010, MNRAS, 408, 731.
- De Marchi et al. (2007) De Marchi, G., Paresce, F. & Pulone, L. 2007, ApJ, 656, L65.
- Dib et al. (2007) Dib, S., Kim, J. & Shadmehri, M. 2007, MNRAS, 381, L40.
- Dib et al. (2010) Dib, S., Shadmehri, M., Padoan, P., et al. 2010, MNRAS, 405, 401.
- Evans et al. (2010) Evans, C.J., et al. 2010, ApJ, 715, L74.
- Fujii & Portegies Zwart (2011) Fujii, M.S. & Portegies Zwart, S. 2011, Science, 334, 1380.
- Heggie & Hut (2003) Heggie, D.C., Hut, P., 2003, “The Gravitational Millon-Body Problem: A Multidisciplinary Approach to Star Cluster Dynamics”, Cambridge University Press, Cambridge, UK.
- Hurley et al. (2000) Hurley, J.R., Pols, O.R. & Tout, C.A. 2000, MNRAS, 315, 543.
- Hurley et al. (2002) Hurley, J.R., Tout, C.A. & Pols, O.R. 2002, MNRAS, 329, 897.
- Kroupa, Gilmore & Tout (1991) Kroupa, P., Gilmore, G. & Tout, C.A. 1991, MNRAS, 251, 293.
- Kroupa (1995) Kroupa, P. 1995, MNRAS, 277, 1507.
- Kroupa (2001) Kroupa, P. 2001, MNRAS, 332, 231.
- Kroupa (2008) Kroupa, P. 2008, in Aarseth, S.J., Tout, C.A., Mardling, R.A., eds, Lecture Notes in Physics Vol. 760, Initial Conditions for Star Clusters. Springer-Verlag, Berlin, p. 181.
- Kroupa et al. (2012) Kroupa, P., Weidner, C., Pflamm-Altenburg, J., et al. 2012, submitted to Stellar Systems and Galactic Structure, arXiv:1112.3340 (preprint).
- Marks & Kroupa (2010) Marks, M. & Kroupa, P. 2010, MNRAS, 406, 2000.
- Marks et al. (2012) Marks, M., Kroupa, P., Dabringhausen, J. & Pawlowski, M.S. 2012, MNRAS, 422, 2246.
- Marks & Kroupa (2012) Marks, M. & Kroupa, P. 2012, arXiv:1205.1508 (preprint).
- Massey & Hunter (1998) Massey, P. & Hunter, D.A. 1998, ApJ, 493, 180.
- Pflamm-Altenburg & Kroupa (2006) Pflamm-Altenburg, J. & Kroupa, P. 2006, MNRAS, 373, 295.
- Sana & Evans (2011) Sana, H. & Evans, C.J. 2011, IAUS, 272, 474.
- Spitzer (1987) Spitzer, L.Jr. 1987, “Dynamical Evolution of Globular Clusters”, Princeton University Press.
- Vesperini & Heggie (1997) Vesperini, E. & Heggie, D.C. 1997, MNRAS, 289, 898.
- Weidner, Kroupa, & Maschberger (2009) Weidner, C., Kroupa, P. & Maschberger, T. 2009, MNRAS, 393, 663.
- Zinnecker & Yorke (2007) Zinnecker, H. & Yorke, H.W. 2007, ARA&A, 45, 481.