Preparing the next gravitational million-body simulations: evolution of single and binary stars in nbody6++gpu , mocca , and mcluster

We present the implementation of updated stellar evolution recipes in the codes NBODY 6 ++ GPU , MOCCA , and MCLUSTER . We test them through numerical simulations of star clusters containing 1.1 × 10 5 stars (with 2.0 × 10 4 in primordial hard binaries) performing high-resolution direct N -body ( NBODY 6 ++ GPU ) and Monte Carlo ( MOCCA ) simulations to an age of 10 Gyr. We compare models implementing either delayed or core-collapse supernovae mechanisms, a different mass ratio distribution for binaries, and white dwarf (WD) natal kicks enabled/disabled. Compared to NBODY 6 ++ GPU , the MOCCA models appear to be denser, with a larger scatter in the remnant masses, and a lower binary fraction on average. The MOCCA models produce more black holes (BHs) and helium WDs, while NBODY 6 ++ GPU models are characterized by a much larger amount of WD–WD binaries. The remnant kick velocity and escape speed distributions are similar for the BHs and neutron stars (NSs), and some NSs formed via electron-capture supernovae, accretion-induced collapse, or merger-induced collapse escape the cluster in all simulations. The escape speed distributions for the WDs, on the other hand, are very dissimilar. We categorize the stellar ev olution recipes a vailable in NBODY 6 ++ GPU into four levels: the one implemented in previous NBODY 6 ++ GPU and MOCCA versions ( level A ), state-of-the-art prescriptions ( level B ), some in a testing phase ( level C ), and those that will be added in future versions of our codes.