Hostname: page-component-74656f9884-npj9w Total loading time: 0 Render date: 2026-11-06T18:54:52.192Z Has data issue: false hasContentIssue false
  • English
  • Français

DRAGON-III simulations: modelling million-body globular and nuclear star clusters over cosmic time

Published online by Cambridge University Press:  15 October 2026

Kai Wu Open the ORCID record for Kai Wu [Opens in a new window] ,
Philip Cho Open the ORCID record for Philip Cho [Opens in a new window] ,
Rainer Spurzem Open the ORCID record for Rainer Spurzem [Opens in a new window] ,
Long Wang ,
Francesco Flammini Dotti  and
Vahid Amiri
Kai Wu*
Affiliation:
Astronomisches Rechen-Institut, Zentrum für Astronomie der Universität Heidelberg, Mönchhofstr. 12-14, D-69120 Heidelberg, Germany
Philip Cho
Affiliation:
Astronomisches Rechen-Institut, Zentrum für Astronomie der Universität Heidelberg, Mönchhofstr. 12-14, D-69120 Heidelberg, Germany
Rainer Spurzem*
Affiliation:
Astronomisches Rechen-Institut, Zentrum für Astronomie der Universität Heidelberg, Mönchhofstr. 12-14, D-69120 Heidelberg, Germany National Astronomical Observatories, Chinese Academy of Sciences, 20A Datun Rd., Chaoyang District, 100101, Beijing, China Kavli Institute for Astronomy and Astrophysics, Peking University, Yiheyuan Lu 5, Haidian Qu, 100871, Beijing, China
Long Wang
Affiliation:
School of Physics and Astronomy, Sun Yat-sen University (Zhuhai Campus), Zhuhai 519082, China
Francesco Flammini Dotti
Affiliation:
Department of Physics, New York University Abu Dhabi, PO Box 129188 Abu Dhabi, UAE Center for Astrophysics and Space Science (CASS), New York University Abu Dhabi, PO Box 129188, Abu Dhabi, UAE
Vahid Amiri
Affiliation:
Astronomisches Rechen-Institut, Zentrum für Astronomie der Universität Heidelberg, Mönchhofstr. 12-14, D-69120 Heidelberg, Germany
*
email: kaiwu.astro@gmail.com
email: kaiwu.astro@gmail.com
Rights & Permissions [Opens in a new window]

Abstract

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).

Keywords

methods: numericalstars: generalstars: black holesgalaxies: star clusters: general

Information

Type
Contributed Paper
Information
Proceedings of the International Astronomical Union , Volume 21 , Symposium S398: Compact Objects and Binaries in Dense Star Clusters , September 2025 , pp. 184 - 188
Check for updates
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
© The Author(s), 2026. Published by Cambridge University Press on behalf of The International Astronomical Union

1. Introduction

1.1 Multi-messenger astrophysics of globular and nuclear star clusters

Globular clusters (GCs), abundant in galactic disks and spheroids, serve as ideal laboratories for studying stellar evolution alongside Newtonian and relativistic dynamics. Nuclear star clusters (NSCs), which often host supermassive black holes in massive elliptical galaxies, generate tidal disruption events as stars or compact objects interact with the central black hole.

From a comprehensive set of computer models, multi-messenger observational signatures of GCs and NSCs can be deduced, including (1) star-by-star kinematic information for nearby star clusters (Gaia, LSST, JWST), (2) numbers of X-ray binaries in local globular clusters and in our own Galactic Center, (3) gravitational wave (GW) signals (and possible electromagnetic counterparts) of merging binary black holes (BHs) and neutron stars (NSs) in a local universe (LIGO/Virgo/KAGRA and future Einstein Telescope Branchesi et al. Reference Branchesi, Maggiore and Alonso2023 Cosmic Explorer Ricci et al. Reference Ricci, Carpenter and Fu2015), and (4) gravitational inspiral events produced by objects falling into BHs to be discovered years before their coalescence (LISA Amaro-Seoane et al. Reference Amaro-Seoane, Audley and Babak2017, TianQin Luo et al. Reference Luo, Chen and Duan2016, and Taiji Chen et al. Reference Chen, Yan and Lu2021; Ruan et al. Reference Ruan, Liu and Guo2021).

1.2 The DRAGON simulation family

The DRAGON simulations aim to model the largest star clusters to the most realistic degree, using accurate N-body methods with up-to-date treatment for stellar evolution recipes.

DRAGON-I (Wang et al. Reference Wang, Spurzem and Aarseth2016) is the first N-body simulation to simulate million-body GCs for Giga-years while DRAGON-II (Arca Sedda et al. Reference Arca Sedda, Kamlah and Spurzem2023; Arca Sedda et al. Reference Arca Sedda, Kamlah and Spurzem2024a; Arca Sedda et al. Reference Arca Sedda, Kamlah and Spurzem2024b) significantly improves stellar evolution models involving compact objects, and predicts gravitational wave signals from compact object mergers that would be measured by LIGO/Virgo/KAGRA.

The aim of DRAGON-III is to get step by step closer to realistic particle numbers for massive globular and nuclear star clusters. Furthermore, the inclusion of a novel feature, a large number of initially wide binaries, significantly changes the entire binary population during the evolution. Moreover, we follow all escaping stars to investigate how they build tidal tails with a more recent galaxy model.

2. Numerical methods

2.1 Code

NBODY6++GPU (Wang et al. Reference Wang, Spurzem and Aarseth2015) is one of the state-of-the-art codes capable of modelling star clusters to the most realistic degree possible. It has a long history as an offspring of Sverre’s N-body code (Aarseth Reference Aarseth1999), and is being continuously maintained by the Silkroad Team Footnote † and made publicly available Footnote ‡, with a major update for stellar evolution by Kamlah et al. (Reference Kamlah, Leveque and Spurzem2022). We refer readers to the review article Spurzem and Kamlah (Reference Spurzem and Kamlah2023) for its physics, algorithms, features, and achievements. Beyond the descriptions therein, we recently ported the Galactic potential MWPotential2014 from Galpy (Bovy Reference Bovy2015), for better comparison with results from other codes. In recent years, the MWPotential2014 has been adopted by other star cluster codes (e.g., MOCCA by Giersz Reference Giersz1998; PeTar by Wang et al. Reference Wang, Iwasawa and Nitadori2020). We have ported the CUDA code in NBODY6++GPU to HIP Footnote § software stacks, enabling our program to run large-scale simulations on next-generation non-CUDA architectures, including DCU-accelerated computing clusters from China and AMD-GPU-based LUMI supercomputer in Finland.

PeTar Footnote is one of the most popular codes for star cluster N-body simulations. It was first developed by Wang et al. (Reference Wang, Iwasawa and Nitadori2020), and then enriched with functions such as BSE (Wang et al. Reference Wang, Fujii and Tanikawa2021), BSEEMP (Wang et al. Reference Wang, Tanikawa and Fujii2022), and MWPotential2014 (Wang and Jerabkova Reference Wang and Jerabkova2021). We refer interested readers to the article of Wang (Reference Wang2026) in the same proceeding volume.

2.2 Initial conditions

Table 1 shows the initial conditions of the ongoing and planned DRAGON-III simulations. We follow the model with 1 million stars in DRAGON-II for unlisted initial parameters, and refer readers to (Arca Sedda et al. Reference Arca Sedda, Kamlah and Spurzem2024a, Section 2.3) for them. Inspired by Offner et al. (Reference Offner, Moe, Kratter, Inutsuka, Aikawa and Muto2023), we explore three different fractions of initial soft (wide) binaries. We do not remove stars from simulations when they are far enough from the cluster centre and usually considered as escaping. Instead, we keep evolving them in the galactic tidal field. This allows us to investigate possible stellar streams formed by massive clusters.

Table 1.

Initial conditions for DRAGON-III simulations.

A table comparing initial conditions for DRAGON-III simulations.
Table 1 Long description

A table comparing initial conditions for DRAGON-III simulations. The table has two columns: DRAGON-III GC and DRAGON-III NSC with accreting CMBH. It has six rows with the following headers: Number of stars, Integration Time, Orbit in the Galaxy, Tidal field, Binary, and Stellar evolution packages. Row 1: Number of stars, DRAGON-III GC: 1M, 4M, 10M; DRAGON-III NSC with accreting CMBH: 1M. Row 2: Integration Time, DRAGON-III GC: 1 – 10 Gyr; DRAGON-III NSC with accreting CMBH: ~ Gyr. Row 3: Orbit in the Galaxy, DRAGON-III GC: a = 10 kpc; e = 0.5; inc = 60°; DRAGON-III NSC with accreting CMBH: In the galaxy centre. Row 4: Tidal field, DRAGON-III GC: Galpy MWPotential2014; DRAGON-III NSC with accreting CMBH: (not specified). Row 5: Binary, DRAGON-III GC: hard: 5% & soft: 0% / 20% / 60%; DRAGON-III NSC with accreting CMBH: (not specified). Row 6: Stellar evolution packages, DRAGON-III GC: Kamlah et al. (2022); Banerjee et al. (2020); Banerjee (2021), IMF & profile: Kroupa (2001) 0.08 – 150 M☉ & King (1966); DRAGON-III NSC with accreting CMBH: (not specified).

Notes: [1] The accreting central supermassive black hole is added in NSC simulations. [2] The orbit of the GC is taken typical observational values from Vasiliev and Baumgardt (Reference Vasiliev and Baumgardt2021); Bajkova and Bobylev (Reference Bajkova and Bobylev2021) [3] In the initial setup, hard binaries means those with semi-major axis of 0.1 normal a normal u less than or equals a less than or equals 50 normal a normal u $0.1~\mathrm{au} \leq a \leq 50~\mathrm{au}$ , while soft binary means 50 normal a normal u less than or equals a less than or equals 5000 normal a normal u $50~\mathrm{au} \le a \leq 5000~\mathrm{au}$ .

3. Selected preliminary results

Fig. 1 shows binaries related to compact objects in early time of one simulation with NBODY6++GPU. Many interesting objects appear after evolving for 84.16 Myr, especially:

  • • An IMBH is present in a wide BH–BH binary, with primary mass of 127 Me and secondary mass of 89 Me. This binary is escaping from the star cluster, and will likely be dynamically unbound to the cluster.

  • • An BH–MS binary ( a equals 1.29 $a=1.29$ au, upper M Subscript normal upper B normal upper H Baseline equals 5.82 $M_\mathrm{BH}=5.82$ Me) is located next to the recently observed BH–giant wide binary in the mass gap (Wang et al. Reference Wang, Zhao and Feng2024) ( 2.1 less than a divided by normal a normal u less than 5.1 $2.1 \lt a/\mathrm{au} \lt 5.1$ , 3.1 less than upper M Subscript normal upper B normal upper H Baseline divided by normal upper M Subscript circled dot Baseline less than 4.1 $3.1 \lt M_\mathrm{BH}/\mathrm{{\mathrm{M}_\odot}{}} \lt 4.1$ ). With further evolution, we will see whether this simulated binary will resemble the observed binary.

Figure 1.

Semi-major axis as a function of primary mass for binaries with at least one compact object component, in the simulation with 20% initial soft binary after evolving for 84.16 Myr, which equals 3749 initial crossing time, and 0.5 initial half-mass relaxation time. In the bottom left corner is the count of binaries by stellar types of their components. Meaning of the stellar type abbreviation in this table and in the legend can be found in Hurley et al. (Reference Hurley, Pols and Tout2000). Meaning of marker style: square - black hole as primary star, octagon - neutron star as primary star, pentagon - white dwarf as primary star.

Scatter plot showing the relationship between semi-major axis and primary mass for binaries with compact objects.

As shown in Table 2, after evolving the cluster for 100 Myr, there are 41 pulsar sources, 191 X-ray binaries, 17 gravitational wave sources, and one black hole-black hole merger due to losing orbital energy in the form of gravitational wave emission. These results show that young globular clusters produce numerous observational candidates for multi-messenger astrophysics. Some of these binaries will escape from the cluster and scatter into the galactic field. Their distribution and subsequent evolution will be continuously tracked in the running DRAGON-III simulations.

Table 2.

Binaries with at least one compact object component, categorized by stellar types of stellar components, in the simulation with 20% initial soft binary after evolving for 106.8 Myr. Due to page limit, we refer readers to Hurley et al. (Reference Hurley, Pols and Tout2000), Section 4 for the meaning of abbreviations.

A table categorizing binaries with compact object components by stellar types.
Table 2 Long description

A table categorizing binaries with compact object components by stellar types. The table has three rows and six columns. Column headers are Stellar Type, ONeWD-MS, NS-MS, COWD-MS, BH-MS, BH-BH, and ONeWD-CHeB. Row 1: Stellar Type, ONeWD-MS, NS-MS, COWD-MS, BH-MS, BH-BH, ONeWD-CHeB. Row 2: count, 46, 26, 26, 16, 16, 14. Row 3: Stellar Type, BH-LMS, ONeWD-ONeWD, NS-ONeWD, NS-CHeB, ONeWD-HEMS, COWD-LMS. Row 4: count, 9, 8, 7, 6, 2, 2. Row 5: Stellar Type, ONeWD-LMS, COWD-CHeB, NS-HG, BH-CHeB, BH-NS, NS-HEMS. Row 6: count, 2, 2, 1, 1, 1, 1.

The inclusion of initial soft binaries produced not only more wide binaries, but also more close binaries with unexpected properties as discussed above. An upcoming article (Wu et al. Reference Wu, Spurzem and Flammini Dotti2026) will highlight how different initial soft binaries affect the simulation. The DRAGON-III simulations are ongoing. The first DRAGON-III data-release will be publicly available soon.

References

Aarseth, S. J. 1999, PASP, 111(765), 1333–1346.10.1086/316455CrossRefGoogle Scholar
Amaro-Seoane, P., Audley, H., Babak, S., et al. 2017, arXiv e-prints, arXiv:1702.00786.Google Scholar
Arca Sedda, M., Kamlah, A. W. H., Spurzem, R., et al. 2024a, MNRAS, 528(3), 5119–5139.10.1093/mnras/stad3952CrossRefGoogle Scholar
Arca Sedda, M., Kamlah, A. W. H., Spurzem, R., et al. 2024b, MNRAS, 528(3), 5140–5159.10.1093/mnras/stad3951CrossRefGoogle Scholar
Arca Sedda, M., Kamlah, A. W. H., Spurzem, R., et al. 2023, MNRAS, 526(1), 429–442.10.1093/mnras/stad2292CrossRefGoogle Scholar
Bajkova, A. T. & Bobylev, V. V. 2021, Research in Astronomy and Astrophysics, 21(7), 173.10.1088/1674-4527/21/7/173CrossRefGoogle Scholar
Banerjee, S. 2021, MNRAS, 500(3), 3002–3026.10.1093/mnras/staa2392CrossRefGoogle Scholar
Banerjee, S., Belczynski, K., Fryer, C. L., et al. 2020, A&A, 639, A41.Google Scholar
Bovy, J. 2015, ApJS, 216(2), 29.10.1088/0067-0049/216/2/29CrossRefGoogle Scholar
Branchesi, M., Maggiore, M., Alonso, D., et al. 2023, J. Cosmology Astropart. Phys., 2023(7), 068.10.1088/1475-7516/2023/07/068CrossRefGoogle Scholar
Chen, J., Yan, C.-S., Lu, Y.-J., et al. 2021, Research in Astronomy and Astrophysics, 21(11), 285.10.1088/1674-4527/21/11/285CrossRefGoogle Scholar
Giersz, M. 1998, MNRAS, 298(4), 1239–1248.10.1111/j.1365-8711.1998.01734.xCrossRefGoogle Scholar
Hurley, J. R., Pols, O. R., & Tout, C. A. 2000, MNRAS, 315(3), 543–569.10.1046/j.1365-8711.2000.03426.xCrossRefGoogle Scholar
Kamlah, A. W. H., Leveque, A., Spurzem, R., et al. 2022, MNRAS, 511(3), 4060–4089.10.1093/mnras/stab3748CrossRefGoogle Scholar
King, I. R. 1966, AJ, 71, 64.10.1086/109857CrossRefGoogle Scholar
Kroupa, P. 2001, MNRAS, 322(2), 231–246.10.1046/j.1365-8711.2001.04022.xCrossRefGoogle Scholar
Luo, J., Chen, L.-S., Duan, H.-Z., et al. 2016, Classical and Quantum Gravity, 33(3), 035010.10.1088/0264-9381/33/3/035010CrossRefGoogle Scholar
Offner, S. S. R., Moe, M., Kratter, K. M., et al. In Inutsuka, S., Aikawa, Y., Muto, T., et al., editors, Protostars and Planets VII 2023, volume 534 of Astronomical Society of the Pacific Conference Series, 275.10.26624/LSUL6799CrossRefGoogle Scholar
Ricci, L., Carpenter, J. M., Fu, B., et al. 2015, ApJ, 798(2), 124.10.1088/0004-637X/798/2/124CrossRefGoogle Scholar
Ruan, W.-H., Liu, C., Guo, Z.-K., et al. 2021, Research, 2021, 6014164.Google Scholar
Spurzem, R. & Kamlah, A. 2023, Living Reviews in Computational Astrophysics, 9(1), 3.10.1007/s41115-023-00018-wCrossRefGoogle Scholar
Vasiliev, E. & Baumgardt, H. 2021, MNRAS, 505(4), 5978–6002.10.1093/mnras/stab1475CrossRefGoogle Scholar
Wang, L. 2026, Proceedings of the International Astronomical Union.Google Scholar
Wang, L., Fujii, M. S., & Tanikawa, A. 2021, MNRAS, 504(4), 5778–5787.10.1093/mnras/stab1157CrossRefGoogle Scholar
Wang, L., Iwasawa, M., Nitadori, K., et al. 2020, MNRAS, 497(1), 536–555.10.1093/mnras/staa1915CrossRefGoogle Scholar
Wang, L. & Jerabkova, T. 2021, A&A, 655, A71.Google Scholar
Wang, L., Spurzem, R., Aarseth, S., et al. 2016, MNRAS, 458(2), 1450–1465.10.1093/mnras/stw274CrossRefGoogle Scholar
Wang, L., Spurzem, R., Aarseth, S., et al. 2015, MNRAS, 450(4), 4070–4080.10.1093/mnras/stv817CrossRefGoogle Scholar
Wang, L., Tanikawa, A., & Fujii, M. 2022, MNRAS, 515(4), 5106–5120.10.1093/mnras/stac2043CrossRefGoogle Scholar
Wang, S., Zhao, X., Feng, F., et al. 2024, Nature Astronomy, 8, 1583–1591.10.1038/s41550-024-02359-9CrossRefGoogle Scholar
Wu, K., Spurzem, R., Flammini Dotti, F., et al. 2026, working in progress.Google Scholar
Figure 0

Table 1. Initial conditions for DRAGON-III simulations.Table 1 long description.

Figure 1

Figure 1. Semi-major axis as a function of primary mass for binaries with at least one compact object component, in the simulation with 20% initial soft binary after evolving for 84.16 Myr, which equals 3749 initial crossing time, and 0.5 initial half-mass relaxation time. In the bottom left corner is the count of binaries by stellar types of their components. Meaning of the stellar type abbreviation in this table and in the legend can be found in Hurley et al. (2000). Meaning of marker style: square - black hole as primary star, octagon - neutron star as primary star, pentagon - white dwarf as primary star.

Figure 2

Table 2. Binaries with at least one compact object component, categorized by stellar types of stellar components, in the simulation with 20% initial soft binary after evolving for 106.8 Myr. Due to page limit, we refer readers to Hurley et al. (2000), Section 4 for the meaning of abbreviations.Table 2 long description.