Direct N-body Simulations Research Articles

On the survivability of planets in young massive clusters and its implication of planet orbital architectures in globular clusters

ABSTRACT As of 2019 August, among the more than 4000 confirmed exoplanets, only one has been detected in a globular cluster (GC) M4. The scarce of exoplanet detections motivates us to employ direct N-body simulations to investigate the dynamical stability of planets in young massive clusters (YMC), which are potentially the progenitors of GCs. In an N = 128 k cluster of virial radius 1.7 pc (comparable to Westerlund-1), our simulations show that most wide-orbit planets (a ≥ 20 au) will be ejected within a time-scale of 10 Myr. Interestingly, more than $70{{\ \rm per\ cent}}$ of planets with a < 5 au survive in the 100 Myr simulations. Ignoring planet–planet scattering and tidal damping, the survivability at t Myr as a function of initial semimajor axis a0 in au in such a YMC can be described as fsurv(a0, t) = −0.33log10(a0)(1 − e−0.0482t) + 1. Upon ejection, about $28.8{{\ \rm per\ cent}}$ of free-floating planets (FFPs) have sufficient speeds to escape from the host cluster at a crossing time-scale. The other FFPs will remain bound to the cluster potential, but the subsequent dynamical evolution of the stellar system can result in the delayed ejection of FFPs from the host cluster. Although a full investigation of planet population in GCs requires extending the simulations to multiGyr, our results suggest that wide-orbit planets and free-floating planets are unlikely to be found in GCs.

N-body simulations of terrestrial planet growth with resonant dynamical friction

ABSTRACTWe investigate planetesimal accretion via a direct N-body simulation of an annulus at 1 au orbiting a 1 M$\odot$ star. The planetesimal ring, which initially contains N = 106 bodies is evolved into the oligarchic growth phase. Unlike previous lower resolution studies, we find that the mass distribution of planetesimals develops a bump at intermediate mass after the oligarchs form. This feature marks a boundary between growth modes. The smallest planetesimals are packed tightly enough together to populate mean motion resonances with the oligarchs, which heats the small bodies, enhancing their growth. If we depopulate most of the resonances by decreasing the width of the annulus, this effect becomes weaker. To clearly demonstrate the dynamics driving these growth modes, we also examine the evolution of a planetary embryo embedded in an annulus of collisionless planetesimals. In this case, we find that the resonances push planetesimals away from the embryo, decreasing the surface density of the bodies adjacent to the embryo. This effect only occurs when the annulus is wide enough and the mass resolution of the planetesimals is fine enough to populate the resonances. The bump we observe in the mass distribution resembles the 100 km power-law break seen in the size distribution of asteroid belt objects. Although the bump produced in our simulations occurs at a size larger than 100 km, we show that the bump location is sensitive to the initial planetesimal mass, which implies that this feature is potentially useful for constraining planetesimal formation models.

Modelling the Effects of Dark Matter Substructure on Globular Cluster Evolution with the Tidal Approximation

Abstract We present direct N-body simulations of tidally filling 30,000 M⊙ star clusters orbiting between 10 kpc and 100 kpc in galaxies with a range of dark matter substructure properties. The time-dependent tidal force is determined based on the combined tidal tensor of the galaxy’s smooth and clumpy dark matter components, the latter of which causes fluctuations in the tidal field that can heat clusters. The strength and duration of these fluctuations are sensitive to the local dark matter density, substructure fraction, sub-halo mass function, and the sub-halo mass-size relation. Based on the cold dark matter framework, we initially assume sub-halos are Hernquist spheres following a power-law mass function between 105 and 1011M⊙ and find that tidal fluctuations are too weak and too short to affect star cluster evolution. Treating sub-halos as point masses, to explore how denser sub-halos affect clusters, we find that only sub-halos with masses greater than 106M⊙ will cause cluster dissolution times to decrease. These interactions can also decrease the size of a cluster while increasing the velocity dispersion and tangential anisotropy in the outer regions via tidal heating. Hence increased fluctuations in the tidal tensor, especially fluctuations that are due to low-mass halos, do not necessarily translate into mass loss. We further conclude that the tidal approximation can be used to model cluster evolution in the tidal fields of cosmological simulations with a minimum cold dark matter sub-halo mass of 106M⊙, as the effect of lower-mass sub-halos on star clusters is negligible.

Black hole and neutron star mergers in galactic nuclei: the role of triples

ABSTRACT Nuclear star clusters that surround supermassive black holes (SMBHs) in galactic nuclei are thought to contain large numbers of black holes (BHs) and neutron stars (NSs), a fraction of which form binaries and could merge by Kozai–Lidov oscillations (KL). Triple compact objects are likely to be present, given what is known about the multiplicity of massive stars, whose life ends either as an NS or a BH. In this paper, we present a new possible scenario for merging BHs and NSs in galactic nuclei. We study the evolution of a triple black hole (BH) or neutron star (NS) system orbiting an SMBH in a galactic nucleus by means of direct high-precision N-body simulations, including post-Newtonian terms. We find that the four-body dynamical interactions can increase the KL angle window for mergers compared to the binary case and make BH and NS binaries merge on shorter time-scales. We show that the merger fraction can be up to ∼5–8 times higher for triples than for binaries. Therefore, even if the triple fraction is only ∼10–$20\rm{\,per\,cent}$ of the binary fraction, they could contribute to the merger events observed by LIGO/VIRGO in comparable numbers.

StePS: A multi-GPU cosmological N-body Code for compactified simulations

Abstract We present the multi-GPU realization of the StePS (Stereographically Projected Cosmological Simulations) algorithm with MPI-OpenMP-CUDA hybrid parallelization and nearly ideal scale-out to multiple compute nodes. Our new zoom-in cosmological direct N-body simulation method simulates the infinite universe with unprecedented dynamic range for a given amount of memory and, in contrast to traditional periodic simulations, its fundamental geometry and topology match observations. By using a spherical geometry instead of periodic boundary conditions, and gradually decreasing the mass resolution with radius, our code is capable of running simulations with a few gigaparsecs in diameter and with a mass resolution of ∼ 1 0 9 M ⊙ in the center in four days on three compute nodes with four GTX 1080 Ti GPUs in each. The code can also be used to run extremely fast simulations with reasonable resolution for fitting cosmological parameters. These simulations are useful for prediction needs of large surveys. The StePS code is publicly available for the research community.

Resonant Relaxation in Globular Clusters

Abstract Resonant relaxation has been discussed as an efficient process that changes the angular momenta of stars orbiting around a central supermassive black hole due to the fluctuating gravitational field of the stellar cluster. Other spherical stellar systems, such as globular clusters, exhibit a restricted form of this effect where enhanced relaxation rate only occurs in the directions of the angular momentum vectors, but not in their magnitudes; this is called vector resonant relaxation (VRR). To explore this effect, we performed a large set of direct N-body simulations, with up to 512k particles (where k =1024) and ∼500 dynamical times. Contrasting these simulations, which naturally include the collective effects, with Spitzer-style Monte Carlo simulations, which by design only exhibit two-body relaxation, we show that the temporal behavior of the angular momentum vectors in N-body simulations cannot be explained by two-body relaxation alone. VRR operates efficiently in globular clusters with N > 104. The fact that VRR operates in globular clusters may open a way to use powerful tools in statistical physics for their description. In particular, since the distribution of orbital planes relaxes much more rapidly than the distribution of the magnitude of angular momentum and the radial action, the relaxation process reaches an internal statistical equilibrium in the corresponding part of phase space while the whole cluster is generally out of equilibrium, in a state of quenched disorder. We point out the need to include effects of VRR in Monte Carlo simulations of globular clusters.

Spectroscopic and dynamical properties of comet C/2018 F4, likely a true average former member of the Oort cloud

Context. The population of comets hosted by the Oort cloud is heterogeneous. Most studies in this area have focused on highly active objects, those with small perihelion distances or examples of objects with peculiar physical properties and/or unusual chemical compositions. This may have produced a biased sample of Oort cloud comets in which the most common objects may be rare, particularly those with perihelia well beyond the orbit of the Earth. Within this context, the known Oort cloud comets may not be representative of the full sample meaning that our current knowledge of the appearance of the average Oort cloud comet may not be accurate. Comet C/2018 F4 (PANSTARRS) is an object of interest in this regard. Aims. Here, we study the spectral properties in the visible region and the cometary activity of C/2018 F4, and we also explore its orbital evolution with the aim of understanding its origin within the context of known minor bodies moving along nearly parabolic or hyperbolic paths. Methods. We present observations obtained with the 10.4 m Gran Telescopio Canarias (GTC) that we use to derive the spectral class and visible slope of C/2018 F4 as well as to characterise its level of cometary activity. Direct N-body simulations are carried out to explore its orbital evolution. Results. The absolute magnitude of C/2018 F4 is Hr > 13.62 ± 0.04 which puts a strong limit on its diameter, D < 10.4 km, assuming a pV = 0.04 cometary-like value of the albedo. The object presents a conspicuous coma, with a level of activity comparable to those of other comets observed at similar heliocentric distances. Comet C/2018 F4 has a visible spectrum consistent with that of an X-type asteroid, and has a spectral slope S′ = 4.0 ± 1.0%/1000 Å and no evidence of hydration. The spectrum matches those of well-studied primitive asteroids and comets. The analysis of its dynamical evolution prior to discovery suggests that C/2018 F4 is not of extrasolar origin. Conclusions. Although the present-day heliocentric orbit of C/2018 F4 is slightly hyperbolic, both its observational properties and past orbital evolution are consistent with those of a typical dynamically old comet with an origin in the Oort cloud.

Globular cluster tidal interactions and mergers in the Galactic disc

AbstractThe Galactic globular cluster system went and is still going through dynamical processes that require to be explored in detail. Here we illustrate how primordial massive globular clusters born in the Milky Way’s disc evolved by stripping material from each other or even merging very early during their lives. These processes might explain the puzzling presence of star-by-star spreads in iron content observed in massive globular clusters and should be taken into account when studying globular cluster stellar populations. In this context, we show how the direct comparison between the predictions provided by our direct N-body simulations and observations can shed light on the origin and chemo-dynamical evolution of globular clusters.

Brownian motion of supermassive black holes in galaxy cores

AbstractWe investigate the dynamics of supermassive black holes (SMBHs) in galactic cores by means of a semi-analytic model based on the Langevin equation, including dynamical friction and stochastic noise accounting for the gravitational interactions with stars. The model is validated against direct N-body simulations of intermediate-mass black holes in stellar clusters where a realistic number of particles is accessible. For the galactic case, we find that the SMBH experiences a Brownian-like motion with a typical displacement from the geometric center of the Galaxy of a few parsecs, for system parameters compatible with M87.

Black Hole dynamics in Young Star Clusters

AbstractYoung star clusters are a promising environment for forming binary black holes. Such binaries may form dynamically or via binary star evolution or through the interplay of these two channels. To study these formation pathways, we have performed high precision direct N-body simulations of low-mass (M < 1000 M⊙) young star clusters. The simulations were carried out with the code Nbody6++GPU coupled with the population synthesis code MOBSE. Our results highlight the importance of dynamics to form massive black hole binaries even in low-mass young star clusters.

Unexpectedly high formation rate of merging binary black holes in open clusters

AbstractGravitational wave direct detections suggest that 30 M⊙ binary black holes (BBHs) commonly exist in the universe. One possible formation scenario of such BBHs is dynamical three-body encounters in dense star clusters. We performed a series of direct N-body simulations with a mass of 2500 and 10000 M⊙ and found a new channel for the formation of BBHs which is dominant in open clusters. In open clusters, the core-collapse time is shorter than in globular clusters, and therefore massive main-sequence (MS) binaries can form before they evolve to BHs. These MS binaries experience common envelope evolution and evolve to hard BBHs, which can merge within the Hubble time. The number of BBH mergers per unit mass obtained from our simulations reached 20–50 % of that for globular clusters, assuming an initial cluster mass function. Thus, open clusters can be a dominant formation site of hard BBHs.

N-body chaos, phase-space transport and relaxation in numerical simulations

AbstractUsing direct N-body simulations of self-gravitating systems we study the dependence of dynamical chaos on the system size N. We find that the N-body chaos quantified in terms of the largest Lyapunov exponent Λmax decreases with N. The values of its inverse (the so-called Lyapunov time tλ) are found to be smaller than the two-body collisional relaxation time but larger than the typical violent relaxation time, thus suggesting the existence of another collective time scale connected to many-body chaos.

Searching for the GD-1 stream progenitor inGaiaDR2 with directN-body simulations

We perform a large suite of direct N-body simulations aimed at revealing the location of the progenitor, or its remnant, of the GD-1 stream. Data from \gaia\ DR2 reveals the GD-1 stream extends over $\approx 100^\circ$, allowing us to determine the stream's leading and trailing ends. Our models suggest the length of the stream is consistent with a dynamical age of between 2-3 Gyr and the exact length, width and location of the GD-1 stream correspond to the stream's progenitor being located between $-30^\circ < \phi_{1,\mathrm{pro}} < -45^{\circ}$ in the standard GD-1 coordinate system. The model stream density profiles reveal that intact progenitors leave a strong over-density, recently-dissolved progenitors appear as gaps in the stream as escaped stars continue to move away from the remnant progenitor's location, and long-dissolved progenitors leave no observational signature on the remaining stream. Comparing our models to the GD-1 stream yields two possible scenarios for its progenitor's history: a) the progenitor reached dissolution approximately 500 Myr ago during the cluster's previous perigalactic pass and is both located at and responsible for the observed gap at $\phi_1=-40^{\circ}$ or b) the progenitor reached dissolution over 2.5 Gyr ago, the fully-dissolved remnant is at $-30^\circ < \phi_1 < -45^{\circ}$, and an observational signature of its location no longer exists. That the dissolved progenitor is in the range $-30^\circ < \phi_1 < -45^{\circ}$ implies that density fluctuations outside of this range, e.g., a deep gap at $\phi_1 \approx -20^\circ$, are likely produced by compact baryonic or dark-matter perturbers.

Correlation between mass segregation and structural concentration in relaxed stellar clusters

Abstract The level of mass segregation in the core of globular clusters has been previously proposed as a potential indicator of the dynamical constituents of the system, such as presence of a significant population of stellar-mass black holes (BHs), or even a central intermediate-mass black hole (IMBH). However, its measurement is limited to clusters with high-quality Hubble Space Telescope data. Thanks to a set of state-of-the-art direct N-body simulations with up to 200k particles inclusive of stellar evolution, primordial binaries, and varying BH/neutron stars, we highlight for the first time the existence of a clear and tight linear relation between the degree of mass segregation and the cluster structural concentration index. The latter is defined as the ratio of the radii containing 5 per cent and 50 per cent of the integrated light (R5/R50), making it robustly measurable without the need to individually resolve low-mass stars. Our simulations indicate that given R5/R50, the mass segregation Δm (defined as the difference in main-sequence median mass between centre and half-light radius) is expressed as Δm/M⊙ = −1.166R5/R50 + 0.3246, with a root-mean-square error of 0.0148. In addition, we can explain its physical origin and the values of the fitted parameters through basic analytical modelling. Such correlation is remarkably robust against a variety of initial conditions (including presence of primordial binaries and IMBHs) and cluster ages, with a slight dependence in best-fitting parameters on the prescriptions used to measure the quantities involved. Therefore, this study highlights the potential to develop a new observational tool to gain insight on the dynamical status of globular clusters and on its dark remnants.

Mergers, tidal interactions, and mass exchange in a population of disc globular clusters

Globular clusters (GCs), the oldest stellar systems observed in the Milky Way, have long been considered single stellar populations. As such, they provided an ideal laboratory to understand stellar dynamics and primordial star formation processes. However, during the last two decades, observations have unveiled their true, complex nature. Beside their pristine stars, GCs host one or more helium enriched and possibly younger stellar populations whose formation mechanism is still unknown. Even more puzzling is the existence of GCs showing star-by-star iron spreads. Using detailed N-body simulations we explore the hypothesis that these anomalies in metallicity could be the result of mutual stripping and mergers between a primordial population of disc GCs. In the first paper of this series we proved, both with analytical arguments and short-term N-body simulations, that disc GCs have larger fly-by and close-encounter rates with respect to halo clusters. These interactions lead to mass exchange and even mergers that form new GCs, possibly showing metallicity spreads. Here, by means of long-term direct N-body simulations, we provide predictions on the dynamical properties of GCs that have undergone these processes. The comparison of our predictions with currently available and future observational data could provide insights into the origin of GCs and the Milky Way build-up history as a whole.

Hydrodynamic Simulations of Moonlet-induced Propellers in Saturn’s Rings: Application to Blériot

Abstract One of the biggest successes of the Cassini mission is the detection of small moons (moonlets) embedded in Saturns rings that cause S-shaped density structures in their close vicinity, called propellers. Here, we present isothermal hydrodynamic simulations of moonlet-induced propellers in Saturn’s A ring that denote a further development of the original model. We find excellent agreement between these new hydrodynamic and corresponding N-body simulations. Furthermore, the hydrodynamic simulations confirm the predicted scaling laws and the analytical solution for the density in the propeller gaps. Finally, this mean field approach allows us to simulate the pattern of the giant propeller Blériot, which is too large to be modeled by direct N-body simulations. Our results are compared to two stellar occultation observations by the Cassini Ultraviolet Imaging Spectrometer (UVIS), which intersect the propeller Blériot. Best fits to the UVIS optical depth profiles are achieved for a Hill radius of 590 m, which implies a moonlet diameter of about 860 m. Furthermore, the model favors a kinematic shear viscosity of the surrounding ring material of ν 0 = 340 cm2 s−1, a dispersion velocity in the range of 0.3 cm s−1 &lt; c 0 &lt; 1.5 cm s−1, and a fairly high bulk viscosity 7 &lt; ξ 0/ν 0 &lt; 17. These large transport values might be overestimated by our isothermal ring model and should be reviewed by an extended model including thermal fluctuations.

Obliquity Tides May Drive WASP-12b’s Rapid Orbital Decay

Abstract Recent analyses have revealed a mystery. The orbital period of the highly inflated hot Jupiter, WASP-12b, is decreasing rapidly. The rate of inspiral, however, is too fast to be explained by either eccentricity tides or equilibrium stellar tides. While dynamical stellar tides are possible, they require a subgiant structure for the star, whereas stellar models point toward a main-sequence host. Here, we show that these hitherto irreconcilable observations might be explained by planetary obliquity tides if planet b’s spin vector is trapped in a high-obliquity state maintained by a secular spin–orbit resonance with an unseen exterior perturbing planet. We derive constraints on the obliquity ( ), reduced tidal quality factor ( ), and perturbing planet parameters ( , ) required to generate the observed orbital decay. Direct N-body simulations that include tidal and spin dynamics reinforce the plausibility of the scenario. Furthermore, we show that the resonance could have been captured when planet b’s obliquity was small, making the proposed sequence of events easy to explain. The hypothetical perturbing planet is within the limits of current radial velocity constraints on the system, yet it is also detectable. If it exists, it could provide evidence in favor of the in situ formation hypothesis for hot Jupiters.

A new hybrid technique for modeling dense star clusters

The “gravitational million-body problem,” to model the dynamical evolution of a self-gravitating, collisional N-body system with ∼106 particles over many relaxation times, remains a major challenge in computational astrophysics. Unfortunately, current techniques to model such systems suffer from severe limitations. A direct N-body simulation with more than 105 particles can require months or even years to complete, while an orbit-sampling Monte Carlo approach cannot adequately model the dynamics in a dense cluster core, particularly in the presence of many black holes. We have developed a new technique combining the precision of a direct N-body integration with the speed of a Monte Carlo approach. Our Rapid And Precisely Integrated Dynamics code, the RAPID code, statistically models interactions between neighboring stars and stellar binaries while integrating directly the orbits of stars or black holes in the cluster core. This allows us to accurately simulate the dynamics of the black holes in a realistic globular cluster environment without the burdensome N^{2} scaling of a full N-body integration. We compare RAPID models of idealized globular clusters to identical models from the direct N-body and Monte Carlo methods. Our tests show that RAPID can reproduce the half-mass radii, core radii, black hole ejection rates, and binary properties of the direct N-body models far more accurately than a standard Monte Carlo integration while remaining significantly faster than a full N-body integration. With this technique, it will be possible to create more realistic models of Milky Way globular clusters with sufficient rapidity to explore the full parameter space of dense stellar clusters.

Blue Straggler Bimodality: A Brownian Motion Model

Abstract The shape of the radial distribution of blue straggler stars (BSS), when normalized to a reference population of horizontal branch (HB) stars, has been found to be a powerful indicator of the dynamical evolution reached by a globular cluster (GC). In particular, observations suggest that the BSS distribution bimodality is modulated by the dynamical age of the host GC; dynamically unrelaxed GCs show a flat BSS distribution, and more relaxed GCs show a minimum at a radius that increases for increasing dynamical age, resulting in a natural “dynamical clock.” While direct N-body simulations are able to reproduce the general trend, thus supporting its dynamical origin, the migration of the minimum of the distribution appears to be noisy and not well defined. Here we show that a simple unidimensional model based on dynamical friction (drift) and Brownian motion (diffusion) correctly reproduces the qualitative motion of the minimum, without adjustable parameters except for the BSS to HB stars mass ratio. Differential dynamical friction effects combine with diffusion to create a bimodality in the BSS distribution and to determine its evolution, driving the migration of the minimum to larger radii over time. The diffusion coefficient is strongly constrained by the need to reproduce the migratory behavior of the minimum, and the radial dependence of diffusion set by fundamental physical arguments automatically satisfies this constraint. Therefore, our model appears to capture the fluctuation–dissipation dynamics that underpins the dynamical clock.

Complete ejection of OB stars from very young star clusters and the formation of multiple populations

Recently, three stellar sequences separated in age by about 1 Myr were discovered in the Orion Nebula Cluster (ONC; Beccari et al. 2017). Kroupa et al. (2018) suggest that such small dense subpopulations eject all their OB stars via the decay of unstable few-body systems such that the gas can recombine and form new stars. This explains the multisequence phenomenon without introducing an extra mechanism into star formation theory. In this work, we apply the recently updated primordial binary distribution model (Belloni et al. 2017; implemented here in a new version of MCLUSTER) and perform a large set of direct N-body simulations to investigate the feasibility of this dynamical scenario. Our results suggest that if 3-4 OB stars in the ONC formed primordially mass-segregated in the cluster centre with a maximum separation of about 0.003 pc, all OB stars have a high chance (50-70 per cent) to escape from the centre and do not come back within 1 Myr and the dynamical ejection scenario is a viable channel to form short-age-interval multipopulation sequences as observed in the ONC. This is also consistent with self-regulated star formation.