Our current understanding of the stellar initial mass function and massive star evolution suggests that young globular clusters may have formed hundreds to thousands of stellar-mass black holes, the remnants of stars with initial masses from $\sim 20 - 100\, M_\odot$. Birth kicks from supernova explosions may eject some black holes from their birth clusters, but most should be retained. Using a Monte Carlo method we investigate the long-term dynamical evolution of globular clusters containing large numbers of stellar black holes. We describe numerical results for 42 models, covering a range of realistic initial conditions, including up to $1.6\times10^6$ stars. In almost all models we find that significant numbers of black holes (up to $\sim10^3$) are retained all the way to the present. This is in contrast to previous theoretical expectations that most black holes should be ejected dynamically within a few Gyr. The main reason for this difference is that core collapse driven by black holes (through the Spitzer "mass segregation instability") is easily reverted through three-body processes, and involves only a small number of the most massive black holes, while lower-mass black holes remain well-mixed with ordinary stars far from the central cusp. Thus the rapid segregation of stellar black holes does not lead to a long-term physical separation of most black holes into a dynamically decoupled inner core, as often assumed previously. Combined with the recent detections of several black hole X-ray binary candidates in Galactic globular clusters, our results suggest that stellar black holes could still be present in large numbers in many globular clusters today, and that they may play a significant role in shaping the long-term dynamical evolution and the present-day dynamical structure of many clusters.

Globular clusters should be born with significant numbers of stellar-mass black holes (BHs). It has been thought for two decades that very few of these BHs could be retained through the cluster lifetime. With masses ~10 MSun, BHs are ~20 times more massive than an average cluster star. They segregate into the cluster core, where they may eventually decouple from the remainder of the cluster. The small-N core then evaporates on a short timescale. This is the so-called Spitzer instability. Here we present the results of a full dynamical simulation of a globular cluster containing many stellar-mass BHs with a realistic mass spectrum. Our Monte Carlo simulation code includes detailed treatments of all relevant stellar evolution and dynamical processes. Our main finding is that old globular clusters could still contain many BHs at present. In our simulation, we find no evidence for the Spitzer instability. Instead, most of the BHs remain well-mixed with the rest of the cluster, with only the innermost few tens of BHs segregating significantly. Over the 12 Gyr evolution, fewer than half of the BHs are dynamically ejected through strong binary interactions in the cluster core. The presence of BHs leads to long-term heating of the cluster, ultimately producing a core radius on the high end of the distribution for Milky Way globular clusters (and those of other galaxies). A crude extrapolation from our model suggests that the BH--BH merger rate from globular clusters could be comparable to the rate in the field.

We study the dynamical evolution of globular clusters using our H\'enon-type Monte Carlo code for stellar dynamics including all relevant physics such as two-body relaxation, single and binary stellar evolution, Galactic tidal stripping, and strong interactions such as physical collisions and binary mediated scattering. We compute a large database of several hundred models starting from broad ranges of initial conditions guided by observations of young and massive star clusters. We show that these initial conditions very naturally lead to present day clusters with properties including the central density, core radius, half-light radius, half-mass relaxation time, and cluster mass, that match well with those of the old Galactic globular clusters. In particular, we can naturally reproduce the bimodal distribution in observed core radii separating the "core-collapsed" vs the "non core-collapsed" clusters. We see that the core-collapsed clusters are those that have reached or are about to reach the equilibrium "binary burning" phase. The non core-collapsed clusters are still undergoing gravo-thermal contraction.

We present a new parallel code for computing the dynamical evolution of collisional N-body systems with up to N~10^7 particles. Our code is based on the the Henon Monte Carlo method for solving the Fokker-Planck equation, and makes assumptions of spherical symmetry and dynamical equilibrium. The principal algorithmic developments involve optimizing data structures, and the introduction of a parallel random number generation scheme, as well as a parallel sorting algorithm, required to find nearest neighbors for interactions and to compute the gravitational potential. The new algorithms we introduce along with our choice of decomposition scheme minimize communication costs and ensure optimal distribution of data and workload among the processing units. The implementation uses the Message Passing Interface (MPI) library for communication, which makes it portable to many different supercomputing architectures. We validate the code by calculating the evolution of clusters with initial Plummer distribution functions up to core collapse with the number of stars, N, spanning three orders of magnitude, from 10^5 to 10^7. We find that our results are in good agreement with self-similar core-collapse solutions, and the core collapse times generally agree with expectations from the literature. Also, we observe good total energy conservation, within less than 0.04% throughout all simulations. We analyze the performance of the code, and demonstrate near-linear scaling of the runtime with the number of processors up to 64 processors for N=10^5, 128 for N=10^6 and 256 for N=10^7. The runtime reaches a saturation with the addition of more processors beyond these limits which is a characteristic of the parallel sorting algorithm. The resulting maximum speedups we achieve are approximately 60x, 100x, and 220x, respectively.

Decades after the first predictions of intermediate-mass black holes (IMBHs) in globular clusters (GCs) there is still no unambiguous observational evidence for their existence. The most promising signatures for IMBHs are found in the cores of GCs, where the evidence now comes from the stellar velocity distribution, the surface density profile, and, for very deep observations, the mass-segregation profile near the cluster center. However, interpretation of the data, and, in particular, constraints on central IMBH masses, require the use of detailed cluster dynamical models. Here we present results from Monte Carlo cluster simulations of GCs that harbor IMBHs. As an example of application, we compare velocity dispersion, surface brightness and mass-segregation profiles with observations of the GC M10, and constrain the mass of a possible central IMBH in this cluster. We find that, although M10 does not seem to possess a cuspy surface density profile, the presence of an IMBH with a mass up to 0.75% of the total cluster mass, corresponding to about 600 Msun, cannot be excluded. This is also in agreement with the surface brightness profile, although we find it to be less constraining, as it is dominated by the light of giants, causing it to fluctuate significantly. We also find that the mass-segregation profile cannot be used to discriminate between models with and without IMBH. The reason is that M10 is not yet dynamically evolved enough for the quenching of mass segregation to take effect. Finally, detecting a velocity dispersion cusp in clusters with central densities as low as in M10 is extremely challenging, and has to rely on only 20-40 bright stars. It is only when stars with masses down to 0.3 Msun are included that the velocity cusp is sampled close enough to the IMBH for a significant increase above the core velocity dispersion to become detectable.

A promising mechanism to form intermediate-mass black holes (IMBHs) is the runaway merger in dense star clusters, where main-sequence stars collide and form a very massive star (VMS), which then collapses to a black hole. In this paper we study the effects of primordial mass segregation and the importance of the stellar initial mass function (IMF) on the runaway growth of VMSs using a dynamical Monte Carlo code for N-body systems with N as high as 10^6 stars. Our code now includes an explicit treatment of all stellar collisions. We place special emphasis on the possibility of top-heavy IMFs, as observed in some very young massive clusters. We find that both primordial mass segregation and the shape of the IMF affect the rate of core collapse of star clusters and thus the time of the runaway. When we include primordial mass segregation we generally see a decrease in core collapse time (tcc). Moreover, primordial mass segregation increases the average mass in the core, thus reducing the central relaxation time, which also decreases tcc. The final mass of the VMS formed is always close to \sim 10^-3 of the total cluster mass, in agreement with the previous studies and is reminiscent of the observed correlation between the central black hole mass and the bulge mass of the galaxies. As the degree of primordial mass segregation is increased, the mass of the VMS increases at most by a factor of 3. Flatter IMFs generally increase the average mass in the whole cluster, which increases tcc. For the range of IMFs investigated in this paper, this increase in tcc is to some degree balanced by stellar collisions, which accelerate core collapse. Thus there is no significant change in tcc for the somewhat flatter global IMFs observed in very young massive clusters.

We present a new technique to quantify cluster-to-cluster variations in the observed present-day stellar mass functions of a large sample of star clusters. Our method quantifies these differences as a function of both the stellar mass and the total cluster mass, and offers the advantage that it is insensitive to the precise functional form of the mass function. We applied our technique to data taken from the ACS Survey for Globular Clusters, from which we obtained completeness-corrected stellar mass functions in the mass range 0.25-0.75 M$_{\odot}$ for a sample of 27 clusters. The results of our observational analysis were then compared to Monte Carlo simulations for globular cluster evolution spanning a range of initial mass functions, total numbers of stars, concentrations, and virial radii. We show that the present-day mass functions of the clusters in our sample can be reproduced by assuming an universal initial mass function for all clusters, and that the cluster-to-cluster differences are consistent with what is expected from two-body relaxation. A more complete exploration of the initial cluster conditions will be needed in future studies to better constrain the precise functional form of the initial mass function. This study is a first step toward using our technique to constrain the dynamical histories of a large sample of old Galactic star clusters and, by extension, star formation in the early Universe.

We present results of a series of Monte Carlo simulations investigating the imprint of a central intermediate-mass black hole (IMBH) on the structure of a globular cluster. We investigate the three-dimensional and projected density profiles, and stellar disruption rates for idealized as well as realistic cluster models, taking into account a stellar mass spectrum and stellar evolution, and allowing for a larger, more realistic, number of stars than was previously possible with direct N-body methods. We compare our results to other N-body and Fokker-Planck simulations published previously. We find, in general, very good agreement for the overall cluster structure and dynamical evolution between direct N-body simulations and our Monte Carlo simulations. Significant differences exist in the number of stars that are tidally disrupted by the IMBH, which is most likely an effect of the wandering motion of the IMBH, not included in the Monte Carlo scheme. These differences, however, are negligible for the final IMBH masses in realistic cluster models as the disruption rates are generally much lower than for single-mass clusters. As a direct comparison to observations we construct a detailed model for the cluster NGC 5694, which is known to possess a central surface brightness cusp consistent with the presence of an IMBH. We find that not only the inner slope but also the outer part of the surface brightness profile agree well with observations. However, there is only a slight preference for models harboring an IMBH compared to models without.

We study the dynamical evolution of globular clusters containing primordial binaries, including full single and binary stellar evolution using our Monte Carlo cluster evolution code updated with an adaptation of the single and binary stellar evolution codes SSE/BSE from Hurley et. al (2000, 2002). We describe the modifications we have made to the code. We present several test calculations and comparisons with existing studies to illustrate the validity of the code. We show that our code finds very good agreement with direct N-body simulations including primordial binaries and stellar evolution. We find significant differences in the evolution of the global properties of the simulated clusters using stellar evolution compared to simulations without any stellar evolution. In particular, we find that the mass loss from stellar evolution acts as a significant energy production channel simply by reducing the total gravitational binding energy and can significantly prolong the initial core contraction phase before reaching the binary-burning quasi steady state of the cluster evolution as noticed in Paper IV. We simulate a large grid of clusters varying the initial cluster mass, binary fraction, and concentration and compare properties of the simulated clusters with those of the observed Galactic globular clusters (GGCs). We find that our simulated cluster properties agree well with the observed GGC properties. We explore in some detail qualitatively different clusters in different phases of their evolution, and construct synthetic Hertzprung-Russell diagrams for these clusters.

We study the dynamical evolution of the young star cluster Arches and its dependence on the assumed initial stellar mass function (IMF). We perform many direct $N$-body simulations with various initial conditions and two different choices of IMFs. One is a standard Kroupa IMF without any mass segregation. The other is a radially dependent IMF, as presently observed in the Arches. We find that it is unlikely for the Arches to have attained the observed degree of mass segregation at its current age starting from a standard non-segregated Kroupa IMF. We also study the possibility of a collisional runaway developing in the first $\sim 2-3 \rm{Myr}$ of dynamical evolution. We find that the evolution of this cluster is dramatically different depending on the choice of IMF: if a primordially mass segregated IMF is chosen, a collisional runaway should always occur between $2-3 \rm{Myr}$ for a broad range of initial concentrations. In contrast, for a standard Kroupa IMF no collisional runaway is predicted. We argue that if Arches was created with a mass segregated IMF similar to what is observed today then at the current cluster age a very unusual, high-mass star should be created. However, whether a collisional runaway leads to the formation of an intermediate-mass black hole (IMBH) depends strongly on the mass loss rate via winds from massive stars. Growth of stellar mass through collisions can be quenched by strong wind mass loss. In that case, the inter-cluster as well as intra-cluster medium are expected to have a significant Helium enrichment which may be observed via Helium recombination lines. The excess amount of gas lost in winds may also be observed via X-ray observations as diffused X-ray sources.

Recent Spitzer observations of the globular cluster M15 detected dust associated with its intracluster medium. Surprisingly, these observations imply that the dust must be very short-lived compared to the time since the last Galactic plane crossing of the cluster.Here we propose a simple mechanism to explain this short lifetime. We argue that the kinetic energy of the material ejected during a stellar collision may be sufficient to remove the gas and dust entirely from a cluster, or to remove the gas as a wind, in addition to partially destroying the dust. By calculating the rate of stellar collisions using an N-body model for the cluster, we find remarkable agreement between the average time between collisions and the inferred dust lifetime in this cluster, suggesting a possible close relation between the two phenomena. Furthermore, we also obtain the birthrate of blue stragglers formed through collisions in M15. By comparing with the observed number of blue stragglers, we derive an upper limit for their average lifetime which turns out to be consistent with recent model calculations, thereby lending further support to our model.

We present the first results of a series of Monte-Carlo simulations investigating the imprint of a central black hole on the core structure of a globular cluster. We investigate the three-dimensional and the projected density profile of the inner regions of idealized as well as more realistic globular cluster models, taking into account a stellar mass spectrum, stellar evolution and allowing for a larger, more realistic, number of stars than was previously possible with direct N-body methods. We compare our results to other N-body simulations published previously in the literature.

Simulations of the collapse and fragmentation of turbulent molecular clouds and dense young clusters show that encounters between disc-surrounded stars are relatively common events which should significantly influence the resulting disc structure. In turn this should alter the accretion rate of disc matter onto the star and the conditions under which planet formation occurs. Although the effects of star-disc encounters have been previously investigated, very little is known about encounters where both stars are surrounded by discs. In this paper encounters of such disc-disc systems are studied quantitatively. It is found that for low-mass discs ($M_D$= 0.01 $M_\sun$) the results from star-disc encounters can be straightforwardly generalized to disc-disc encounters as long as there is no mass transport between the discs. Differences to star-disc encounters occur naturally where significant amounts of matter are transported between the discs. In this case it is found that although the mass distribution does not change significantly, matter caught onto highly eccentric orbits is transported surprisingly far inside the disc. The captured mass partly replenishes the disc, but has a much lower angular momentum. This can lead to a reduction of the angular momentum in the entire disc and thus considerably increased accretion shortly after the encounter as well as in the long term.

We investigate the dynamical decay of non-hierarchical accreting triple systems and its implications on the ejection model as Brown Dwarf formation scenario. A modified chain-regularization scheme is used to integrate the equations of motion, that also allows for mass changes over time as well as for momentum transfer from the accreted gas mass onto the bodies. We integrate an ensemble of triple systems within a certain volume with different accretion rates, assuming several prescriptions of how momentum is transferred onto the bodies. We follow their evolution until the systems have decayed. We analyze the end states and decay times of these systems and determine the fraction of Brown Dwarfs formed, their escape speeds as well as the semi-major axis distribution of the formed Brown Dwarf binaries. We find that the formation probability of Brown Dwarfs depends strongly on the assumed momentum transfer which is related to the motion of the gas. Due to ongoing accretion and consequent shrinkage of the systems, the median escape velocity is increased by a factor of 2 and the binary separations are decreased by a factor of 5 compared with non-accreting systems. Furthermore, the obtained semi-major axis distribution drops off sharply to either side of the median, which is also supported by observations. We conclude that accretion and momentum transfer of accreted gas during the dynamical decay of triple systems is able to produce the observed distribution of close binary Brown Dwarfs, making the ejection model a viable option as Brown Dwarf formation scenario.

Observations of brown dwarfs provide important feedback on theories of atmospheres and inner structure of substellar objects. Brown dwarf binary systems furthermore offer the unique opportunity to determine the mass of individual brown dwarfs, which is one of the fundamental astrophysical quantities.