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 \texttt{Nbody6++GPU, MOCCA} and \texttt{McLuster}. We test them through numerical simulations of star clusters containing $1.1\times 10^5$ stars (with $2.0\times 10^4$ in primordial hard binaries) performing high-resolution direct $N$-body (\texttt{Nbody6++GPU}) and Monte-Carlo (\texttt{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 natal kicks enabled/disabled. Compared to \texttt{Nbody6++GPU}, the \texttt{MOCCA} models appear to be denser, with a larger scatter in the remnant masses, and a lower binary fraction on average. The \texttt{MOCCA} models produce more black holes (BHs) and helium white dwarfs (WDs), whilst \texttt{Nbody6++GPU} models are characterised 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 categorise the stellar evolution recipes available in \texttt{Nbody6++GPU}, \texttt{MOCCA} and \texttt{Mcluster} into four levels: the one implemented in previous \texttt{Nbody6++GPU} and \texttt{MOCCA} versions (\texttt{level A}), state-of-the-art prescriptions (\texttt{level B}), some in a testing phase (\texttt{level C}), and those that will be added in future versions of our codes.