After the third LIGO–Virgo–KAGRA observing run, the number of detected binary black hole (BBH) mergers became sufficient to identify statistical features of the population. We explore how different prescriptions for the final fate of massive stars and key binary-evolution processes shape isolated binaries and their remnants. Using sevn, we evolved $10^{7}$ binaries across 15 metallicities, 3 core-collapse supernova models, 4 PPISN models, and 6 common-envelope (CE) prescriptions, for a total of 990 runs ($9.9 \times 10^{9}$ systems). Both single- and binary-star physics shape the BH mass distribution: single-star processes control the high-mass tail ($M_{\rm BH} \geq 45\,M_{\odot}$), while binary evolution produces pile-ups in specific intervals. In particular, the bump at $\sim 35\,M_{\odot}$, often attributed to PPISNe, also emerges from binaries evolving only through stable mass transfer, without CE. Finally, we test a top-heavy IMF, finding it boosts merger numbers and alters the abundance of systems with a given primary BH mass.

We present a review of recent advances in exploring and interpreting the full six-dimensional phase space of globular clusters, as enabled by Gaia, Hubble Space Telescope, and large spectroscopic surveys. Observations reveal widespread internal rotation, velocity anisotropy, and kinematic differences among multiple stellar populations, requiring a departure from classical dynamical models. State-of-the-art million-body N-body simulations reproduce these features and demonstrate the impact of anisotropy, rotation, and tides on star cluster evolution, including core-collapse timescales. We outline three key scientific goals: distinguishing primordial from evolutionary signatures, bridging the gap between idealised initial conditions and realistic formation scenarios, and probing the presence of black holes and dark matter within star clusters. Taken together, these developments are confirming globular clusters as essential proving grounds for stellar dynamics, with significant implications for galaxy assembly, gravitational wave astrophysics, and even fundamental physics.

Omega Centauri ($\omega$ Cen) is the most massive globular cluster in the Galaxy, providing the ideal environment for a large number of potential stellar mass black holes to be initially formed, as well as being retained due to the large escape velocity and extended relaxation time. However, $\omega$ Cen has been the focus of a long debate over whether there is a sub-population of stellar mass black holes or an intermediate mass black hole. In this work, we use roughly 14 years of archival Hubble Space Telescope observations to search for potential microlensing events caused by stellar mass black holes. We create synthetic data from simulations of $\omega$ Cen, which allows us to optimize our regression model and provides an idea of what real microlensing event statistics will reflect. We also report our detection completeness of the regression model and detail ongoing candidate follow-up.

Millisecond pulsars (MSPs) are neutron stars with spin periods as short as a few milliseconds, formed through mass accretion from companion stars. In the dense environments of globular clusters (GCs), MSPs are likely to originate through dynamically assembled interacting binaries. Over 300 MSPs have been detected in GCs to date, more than half of the known Galactic MSP population. In this work, we model MSP populations in intermediate-mass GCs using the direct N-body code NBODY6++GPU. We update the code by implementing pulsar spin-down due to magnetic braking and spin-up through accretion, and use this framework to model the pulsar population in the globular cluster M71 to investigate the pulsar population within and the associated gravitational wave transients.

Understanding the origin of hypervelocity stars (HVSs) in the Galactic Center is a critical question, particularly regarding the dynamic processes that have shaped this region. In this study, we assert that HVSs likely formed as counterparts to S stars, specifically through the Hills mechanism originating in the young stellar disk. By employing N-body simulations, we thoroughly analyze the statistical likelihood of HVS production in our models and robustly compare these results with observational data of a recently detected HVS, S5-HVS1, with an origin in the Galactic Center. Our findings provide a consistent model that enables the production of HVSs like S5-HVS1 and captures the known dynamics of the Galactic Center.

The discovery of wide, eccentric Gaia neutron stars (NSs) in binaries with still evolving (likely main sequence) companions offers a new probe of mass transfer and pulsar recycling beyond the compact-binary regime. We model their origins and fates using population synthesis with POSYDON and detailed binary evolution with MESA, contrasting two limiting prescriptions at Roche-lobe overflow (RLOF): enforced circularization versus explicitly eccentric mass transfer. Our MESA setups include updated treatments of eccentric, non-conservative transfer, magnetic-braking torques for cool stars, and neutron-star spin evolution with accretion and dipole spindown. Under optimistic assumptions, isolated evolution yields Gaia-like systems at only $\leq$1.5% relative rates of NS-evolving companion binaries, yet absolute numbers remain consistent with detections for continuous star formation. Synthetic populations indicate that many canonical millisecond pulsar-white dwarf (WD) binaries arise from unstable mass transfer and common envelope recycling, whereas Gaia systems typically avoid common envelope and only undergo stable mass transfer. In the case of capping accretion onto the NS at the Eddington rate, circular RLOF keeps the donor’s mass-loss rate hovering around the Eddington limit and sustained over long timescales. Eccentric mass transfer instead produces briefer burst-y signature where the donor’s mass-loss rate can climb up to a thousand times higher than in the circular case. The eccentric channel then leaves wide, eccentric NS–helium WD binaries with only mild recycling, whereas the circular channel yields long-lived transfer, circular NS–WD binaries (helium or carbon-oxygen core), and fully recycled millisecond pulsars.

We present a study using the (BPASS) code suite that predicts the Galactic population of binaries that contain a black hole or neutron star. By incorporating the stellar evolution models from BPASS with a Milky Way analogue galaxy from the Feedback in Realistic Environment (FIRE) simulation, we can generate a theoretical population of quiescent compact remnant binary systems as well as x-ray faint binaries within our Galaxy. The study of how compact object masses evolve over time plays a key role for time domain astronomy and gravitational wave progenitors. We find that the agreement between our synthetic and observed systems is reasonable, but there are strong indications that new physical processes should be explored within BPASS to be able to accurately reproduce the observed compact remnant distributions.

The recent detection of high-redshift supermassive black holes with JWST has renewed interest in the processes driving black hole growth. At the same time, both simulations and observations point to a widespread population of off-center intermediate-mass black holes in dwarf galaxies. Their ability to merge outside galactic centers may play a key role in shaping black hole mass assembly. Here, we investigate the dynamics of off-center black holes in dwarf galaxies hosting cored dark matter haloes, where long dynamical friction timescales and core stalling naturally arise. By embedding off-center black holes into an idealized galactic potential and scanning a wide range of orbital configurations, we assess the likelihood of close interactions through Jacobi capture. We find that captures occur in about 13% of cases. Such captures, possibly sustained within compact stellar systems like stripped nuclei or globular clusters, represent a crucial first step toward assembling massive black hole binaries beyond galactic centers.

The hierarchical nature of galaxy formation in the $\Lambda$CDM framework often leads to multiple supermassive black holes (SMBHs) in the galactic nuclei. The timescale over which galaxies merge, plays a crucial role in shaping the dynamical evolution and the merger dynamics of their central SMBHs. While binary SMBH evolution is well studied, the long-term dynamics of triple SMBH systems, particularly in non-spherical potentials, remain less understood. We investigate the role of triaxiality in the evolution and dynamics of triple SMBHs with initial conditions drawn from the ROMULUS25 cosmological simulation, using high-resolution gravitodynamical $\textit{N-body}$ simulations. We explore different orbital configurations and host shapes, tracking the evolution from galactic inspiral to hard binary formation at sub-parsec scales. In all cases, the two most massive SMBHs form a rapidly hardening binary that coalesces within a fraction of a Hubble time, while the third forms a stable hierarchical triple system with the heavier binary, or remains on a wide orbit.

M67 (NGC 2682) is a dynamically evolved old (4 Gyr) open cluster in the Milky Way. Due to its nearness and relatively low levels of dust in the line-of-sight, it has been extensively observed both in surveys and targeted observations. We perform new N-body simulations of M67 using the updated nbody6++gpu and also compare these models to observational data from the Gaia survey. While overall cluster parameters such as mass and main-sequence turnoff mass match observations, further simulations are needed to better match the blue straggler population. Further work is needed to investigate the binary transient systems in the models.

Most intermediate-mass galaxies, including the Milky Way, host both a supermassive black hole (SMBH) and a nuclear stellar cluster (NSC). We simulate the 10 Gyr evolution of lower-mass binaries near the hard-soft boundary, including three-body encounters and von Zeipel-Lidov-Kozai oscillations. Encounter-driven kicks cause both inward and outward migration for hard binaries on timescales shorter than classical binary mass segregation, determining whether the systems survive, merge, or evaporate. The process is, however, too slow to bring binaries into the inner arcsecond within the ages of the S-stars or G-objects, but suggests that i) the S-cluster is contaminated by merger remnants from earlier star formation events and more generally that ii) a primordial population of blue stragglers ($\lesssim 2\mathrm{M}_\odot$) should be scattered towards the outer parts of the inner parsec – providing new constraints on NSC stellar populations and dynamical histories.

Populations of stellar-mass black holes (BHs) in globular clusters (GCs) strongly influence their dynamical evolution and lifetimes. Recently, we used multimass equilibrium models of Milky Way GCs to explore the present-day BH populations of a large sample of clusters, based on several observables, including velocity dispersions, density profiles and stellar mass functions. We have now combined these equilibrium models with the new rapid cluster evolution model “clusterBH”, which simulates the bulk properties of tidally-limited GCs and their BH subsystems over time. These coupled models allow us to place constraints on the initial conditions of real GCs, based only on their observable present-day conditions. These models also provide a framework for probing the highly uncertain physics surrounding the formation of BHs, such as their natal kicks, by experimenting with flexible prescriptions and common assumptions and analyzing the impacts on the formation, evolution and present-day structure of Milky Way GCs. We describe here these new coupled models and present preliminary applications to both mock observations of N-body models and a large sample of real Milky Way GCs.

We perform high-resolution simulations of a $10^{10}{\,\rm M_\odot}$ dwarf galaxy hosting a $10^4{\,\rm M_\odot}$ intermediate-mass black hole (IMBH), reaching $\lesssim0.01$ pc resolution. Our model resolves a multiphase ISM with star formation, radiative feedback, and individual supernovae, while BH accretion is captured using an $\alpha$-disc model that includes Lense-Thirring torques and hence follows spin evolution self-consistently. We find that a nuclear star cluster (NSC) efficiently captures inflowing gas and promotes the formation of a dense circumnuclear disc (CND) on $\sim6$ pc scales. Gravitational torques from the NSC drain angular momentum, circularizing the gas onto the accretion disc and enabling BH growth together with spin-up. The CND undergoes in-situ star formation, which regulates BH accretion through feedback. Our results highlight the interplay between BH accretion, star formation, and feedback in the multiphase ISM, emphasize the role of NSCs in mediating inflows, and point to future cosmological applications. Details are presented in Shin et al. (2025).

As a continuation of DRAGON-II, we present the DRAGON-III project, which focuses on the simulations of million-body globular clusters and nuclear clusters over 10 Gyr. We report on its preliminary results on globular clusters. The first 100 Myr of the simulations have produced 41 pulsars, 191 X-ray binaries, 17 gravitational wave sources, and one black hole-black hole merger due to the loss of orbital energy in the form of gravitational wave emission. The inclusion of initial soft binaries brings surprisingly interesting results, including one IMBH in a binary black hole, and compact object binaries resembling the Gaia-BH1 and the wide black hole-giant binary reported in Wang et al. (2024).

We performed a series of N-body simulations of a star cluster model in our Galactic Center (GC) using realistic scale and numbers of particles. We employed two kinds of models: with and without an Intermediate Mass Black Hole (IMBH). We found that: (1) Resonant Relaxation (RR) takes place. (2) If the change in inclination of S stars driven by RR is observed, it can be a clue to estimate of the enclosed mass and object composition within apocenter of S stars. (3) If there is an IMBH, the magnitude of the change in inclination is similar to that of the RR, but may differ by a few factors. This makes it difficult to accurately estimate enclosed mass of a GC using the RR.

We study tidal disruption events (TDEs) and compact object inspirals in nuclear star clusters (NSCs) hosting a central supermassive black hole (SMBH), focusing on their role in SMBH growth. Using the STARDISK version of the direct N-body code NBODY6++GPU, we perform pilot simulations with two improved models: one for mass fallback from TDEs and another for compact object plunges based on orbital decay timescales. Our results show that mass accretion via TDEs peaks within the first 2 Myr and decreases more rapidly for higher initial SMBH masses, with roughly half the disrupted stellar debris being accreted. Compact object accretion is confined mostly to orbits with pericenters between 4 and 27 Schwarzschild radii and is suppressed by an order of magnitude when inspiral criteria are applied.

Recent observations have demonstrated that globular clusters are systems of greater kinematic and morphological complexity than can be described using spherical non-rotating dynamical models. Here, a new family of self-consistent equilibrium models for quasi-relaxed stellar systems is presented, treating for the first time the case of asynchronous rotation within an external tidal field. Their construction as a perturbative expansion of the classical King models entails a free boundary problem, which is solved using matched asymptotics. The resulting equilibria are found to possess a triaxial structure whose aspherical distortion increases with increased asynchronicity. It is hoped that these models will provide a more realistic equilibrium description for use in investigations of the internal structure and kinematics of star clusters.

. The dense and gaseous environments of active galactic nuclei (AGNs) can catalyze repeated mergers of stellar-mass black holes (BHs), potentially explaining the high-mass tail of binary black hole (BBH) mergers observed by LIGO-Virgo-KAGRA (LVK). We present a semi-analytical population synthesis framework that captures key physical processes in AGN disks, including gas capture, migration, binary pair-up, gas hardening, and dynamical binary-single interactions. Our simulations show that AGN disks can produce hierarchical mergers, especially near migration traps, and may contribute to the high-mass, high-spin BBH population. This work opens prospects for constraining both the AGN channel and supermassive black hole growth with future gravitational-wave detections.

This paper, the second in a series, examines the evolutionary processes shaping pulsars and millisecond pulsars (MSPs). It explores mechanisms such as accretion-driven spin-up, magnetic dipole braking, and binary interactions, presenting seven distinct evolutionary pathways. These include the formation of isolated pulsars, MSPs, and black holes through asymmetric collapse, recycling in low-mass X-ray binaries, dynamical interactions in globular clusters (GCs), and neutron star-neutron star mergers. Insights from observations by FAST, CHIME, and the Parkes Telescope, supported by initial simulation efforts, highlight the role of environmental factors, such as GC density compared to the Galactic field. Ongoing challenges in modeling magnetic field decay and binary dynamics suggest further investigation is warranted, through the accurate simulator, NBODY6++GPU.