BSE versus StarTrack: implementations of new wind, remnant-formation, and natal-kick schemes in NBODY7 and their astrophysical consequences
Sambaran Banerjee, Krzysztof Belczynski, Christopher L. Fryer, Peter, Berczik, Jarrod R. Hurley, Rainer Spurzem, Long Wang
TL;DR
This paper introduces new stellar-remnant formation and natal-kick schemes in NBODY7, demonstrating their astrophysical implications and consistency with population synthesis models, especially regarding black hole and neutron star masses and retention in clusters.
Contribution
It develops and tests new implementations of wind, remnant formation, and natal-kick schemes in NBODY7, aligning results with StarTrack and exploring effects on black hole retention.
Findings
New schemes produce NS and BH masses matching StarTrack.
Different natal kicks influence BH retention in clusters.
Primordial binary mergers affect BH mass distribution and spins.
Abstract
The masses of stellar-remnant black holes (BH), as a result of their formation via massive single- and binary-stellar evolution, is of high interest in this era of gravitational-wave detection from binary black hole (BBH) and binary neutron star (BNS) mergers. Here we present new developments in the N-body evolution program NBODY7 in regards to its stellar-remnant formation and related schemes. We demonstrate that the newly-implemented stellar-wind and remnant-formation schemes in the NBODY7 code's BSE sector, such as the 'rapid' and the 'delayed' supernova (SN) schemes along with an implementation of pulsational-pair-instability and pair-instability supernova (PPSN/PSN), now produces neutron star (NS) and BH masses that agree nearly perfectly, over large ranges of zero-age-main sequence (ZAMS) mass and metallicity, with those from the StarTrack population-synthesis program. We also…
Click any figure to enlarge with its caption.
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5
Figure 6
Figure 7
Figure 8
Figure 9
Figure 10
Figure 11
Figure 12
Figure 13
Figure 14
Figure 15
Figure 16
Figure 17
Figure 18
Figure 19
Figure 20
Figure 21
Figure 22
Figure 23
Figure 24
Figure 25
Figure 26
Figure 27
Figure 28
Figure 29
Figure 30
Figure 31
Figure 32
Figure 33
Figure 34
Figure 35
Figure 36
Figure 37
Figure 38
Figure 39
Figure 40| Metallicity | Remnant scheme | ||||
|---|---|---|---|---|---|
| 0.0001 | F12-rapid+B16-PPSN/PSN | 6.130 | 0.771 | 1.833 | 0.656 |
| 0.0001 | F12-rapid | 6.562 | 0.786 | 1.844 | 0.658 |
| 0.0001 | F12-delayed | 6.973 | 0.601 | 2.580 | 0.353 |
| 0.0001 | B08 | 7.493 | 0.835 | 2.055 | 0.682 |
| 0.002 | F12-rapid+B16-PPSN/PSN | 4.999 | 0.777 | 1.751 | 0.633 |
| 0.002 | F12-rapid | 5.018 | 0.777 | 1.751 | 0.633 |
| 0.002 | F12-delayed | 5.334 | 0.566 | 2.533 | 0.318 |
| 0.002 | B08 | 5.604 | 0.808 | 1.938 | 0.675 |
| 0.02 | F12-rapid | 2.117 | 0.594 | 1.564 | 0.530 |
| 0.02 | F12-delayed | 1.944 | 0.185 | 2.055 | 0.074 |
| 0.02 | B08 | 2.220 | 0.538 | 1.716 | 0.374 |
| 0.0001 | F12-rapid+B16-PPSN/PSN (+MB) | 6.152 | 0.793 | 1.167 | 0.700 |
| 0.0001 | F12-rapid (+MB) | 7.418 | 0.813 | 1.249 | 0.710 |
| 0.0001 | F12-delayed (+MB) | 7.563 | 0.736 | 1.447 | 0.564 |
| 0.0001 | B08 (+MB) | 7.992 | 0.881 | 1.296 | 0.793 |
| StarTrack evolution | BSE evolution | ||||
| 118.5 | 25.5 | 73.6 | 0.73 | ||
| 138.5 | 20.0 | 3806.2 | 0.28 | PSN (no remnant) | |
| 82.3 | 80.0 | 38.8 | 0.58 | ||
| 48.8 | 38.9 | 48.7 | 0.00 | ||
| 39.4 | 17.5 | 26.5 | 0.46 | ||
| 27.2 | 17.0 | 288.1 | 0.57 | ||
| 60.6 | 28.0 | 63.0 | 0.60 | , , () | |
| 33.9 | 24.5 | 270.7 | 0.00 | , , () | |
| BBH merger () | |||||
| 68.7 | 23.1 | 54.2 | 0.78 | ||
| 96.6 | 21.3 | 4126.7 | 0.21 | ||
| 32.4 | 19.6 | 67.1 | 0.00 | ||
| 84.7 | 57.3 | 56.1 | 0.33 | , , () | |
| BBH merger () |
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Code & Models
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
11institutetext: Helmholtz-Instituts für Strahlen- und Kernphysik, Nussallee 14-16, D-53115 Bonn, Germany 22institutetext: Argelander-Institut für Astronomie, Auf dem Hügel 71, D-53121, Bonn, Germany 33institutetext: Nicolaus Copernicus Astronomical Centre of the Polish Academy of Sciences, ul. Bartycka 18, 00-716 Warszawa, Poland 44institutetext: Center for Theoretical Astrophysics, Los Alamos National Laboratory, P.O. Box 1663, Los Alamos, NM 87545, U.S.A. 55institutetext: National Astronomical Observatories and Key Laboratory of Computational Astrophysics, Chinese Academy of Sciences, 20A Datun Rd., Chaoyang District, Beijing 100101, China 66institutetext: Astronomisches Rechen-Institut am Zentrum für Astronomie der Universität Heidelberg, Mönchhofstraße 12-14, 69120 Heidelberg, Germany 77institutetext: Main Astronomical Observatory, National Academy of Sciences of Ukraine, 27 Akademika Zabolotnoho St., 03143 Kyiv, Ukraine 88institutetext: Swinburne University of Technology, Hawthorn VIC 3122, Australia 99institutetext: Kavli Institute for Astronomy and Astrophysics, Peking University, 5 Yi He Yuan Road, Haidian District, Beijing 100871, China 1010institutetext: Department of Astronomy, School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo, 113-0033, Japan 1111institutetext: RIKEN Center for Computational Science, 7-1-26 Minatojima-minami-machi, Chuo-ku, Kobe, Hyogo 650-0047, Japan
BSE versus StarTrack: implementations of new wind, remnant-formation,
and natal-kick schemes in NBODY7 and their astrophysical consequences
S. Banerjee Corresponding author. E-mail: [email protected]
K. Belczynski 33
C. L. Fryer 44
P. Berczik 556677
J. R. Hurley 88
R. Spurzem Research Fellow of the Frankfurt Institute for Advanced Study556699
L. Wang 10101111
Abstract
*Context. *The masses of stellar-remnant black holes (BH), as a result of their formation via massive single- and binary-stellar evolution, is of high interest in this era of gravitational-wave detection from binary black hole (BBH) and binary neutron star (BNS) merger events.
*Aims. *In this work, we present new developments in the state-of-the-art N-body evolution program NBODY7 in regards to its stellar-remnant formation and related schemes. We demonstrate that the newly-implemented stellar-wind and remnant-formation schemes in the NBODY7 code’s stellar-evolutionary sector or BSE , such as the “rapid” and the “delayed” supernova (SN) schemes along with an implementation of pulsational-pair-instability and pair-instability supernova (PPSN/PSN), now produces neutron star (NS) and BH masses that agree nearly perfectly, over large ranges of zero-age-main-sequence (ZAMS) mass and metallicity, with those from the widely-recognized StarTrack population-synthesis program. We also demonstrate the new, recipe-based implementations of various, widely-debated mechanisms of natal kicks on NSs and BHs such as the “convection-asymmetry-driven”, “collapse-asymmetry-driven”, and “neutrino-emission-driven” kicks, in addition to a fully consistent implementation of the standard, fallback-dependent, momentum-conserving natal kick.
*Methods. *All the above newly-implemented schemes are also shared with the standalone versions of SSE and BSE . All these demonstrations are performed with both the updated standalone BSE and the updated .
*Results. *When convolved with stellar and primordial-binary populations as observed in young massive clusters, such remnant-formation and natal-kick mechanisms crucially determine the accumulated number, mass, and mass distribution of the BHs retaining in young massive, open, and globular clusters, which would become available for long-term dynamical processing.
*Conclusions. *Among other interesting conclusions, we find that although the newer delayed-SN remnant formation model gives birth to the largest number (mass) of BHs, the older remnant-formation schemes cause the largest number (mass) of BHs to survive in clusters, when SN material fallback on to the BHs is incorporated. The SN material fallback also causes the convection-asymmetry-driven SN kick to effectively retain similar number and mass of BHs in clusters as for the standard, momentum-conserving kick. The collapse-asymmetry-driven SN kick would cause nearly all BHs to retain in clusters irrespective of its mass, the remnant formation model, and the metallicity, whereas the inference of a large population of BHs in GCs would potentially rule out the neutrino-driven SN kick mechanism. Pre-SN mergers of massive primordial binaries would potentially cause BH masses to deviate from the theoretical, single-star ZAMS mass-remnant mass relation unless a substantial, up to %, of the total merging stellar mass is lost during a merger process. In particular, such mergers, at low metallicities, have the potential to produce low-spinning BHs within the PSN mass gap that can be retained in a stellar cluster and be available for subsequent dynamical interactions. The new remnant-formation modelling reassures that young massive and open clusters would potentially contribute to the dynamical BBH merger detection rate to a similar extent as their more massive globular-cluster counterparts, as recent studies indicate.
Key Words.:
Stars: black holes — Stars: massive — Stars: mass-loss — Stars: kinematics and dynamics — supernovae: general — Methods: numerical
1 Introduction
The formation mechanisms of double compact object binaries undergoing general-relativistic (hereafter GR) inspiral and their occurrence rates has always been a topic of diverse interest (e.g., Abadie et al., 2010) and, naturally, it has become a focal concern since the first direct detection of gravitational waves (hereafter GW) from a stellar-origin binary-black hole (BBH) inspiral and merger by the Laser Interferometer Gravitational-wave Observatory (LIGO; Abbott et al., 2016b). Among the most important ingredients in theoretical modelling and population-synthesis of BBH mergers are the natal masses and kicks of stellar-remnant black holes (hereafter BH) — this, as well, applies to neutron stars (hereafter NS) while modelling NS-containing binaries. In fact, the importance of realistic BH mass estimates became clear right after the very first LIGO detection of BHs (Abbott et al., 2016a; Belczynski et al., 2016b).
In fact, as long as the proto-remnant mass, , in a supernova (hereafter SN) explosion and the amount of material fallback onto it during the SN depend on the carbon-oxygen (CO) core mass and on the pre-supernova stellar mass (Fryer, 1999; Hurley et al., 2000; Fryer & Kalogera, 2001; Fryer et al., 2012), the final remnant mass would depend on the entire life cycle of the progenitor star. In such cases, starting from a given zero-age main sequence (ZAMS) mass, the sequentially-dependent evolutions of the star’s radius, luminosity, mass loss via wind, and total mass would determine the final remnant mass. This is typically the case for the formation of relatively heavier () NSs and BHs. For example, by applying “alternative”, theoretical and semi-empirical wind recipes in the main sequence and evolved stages of O- and B-type progenitors, Belczynski et al. (2010) have demonstrated that BHs, as observed by the LIGO-Virgo (Abbott et al., 2016b, 2017a, 2017b), can easily form from progenitors having dwarf galaxy- and globular cluster (hereafter GC)-like metallicities, which masses are impossible to achieve with the “standard” Hurley et al. (2000) wind prescription even at a very low metallicity.
The formation of compact binaries through both field binary evolution (e.g., Belczynski et al., 2002; Voss & Tauris, 2003; Belczynski et al., 2010, 2016b; Stevenson et al., 2017; Mapelli et al., 2017; Giacobbo et al., 2018) and dynamical interactions in stellar clusters (e.g., Banerjee et al., 2010; Ziosi et al., 2014; Morscher et al., 2015; Rodriguez et al., 2016; Mapelli, 2016; Wang et al., 2016; Askar et al., 2017, 2018; Banerjee, 2017, 2018b, 2018a; Park et al., 2017; Rodriguez et al., 2018) are now widely studied. In both formation channels, stellar wind, BH masses, and natal kicks play a role. A few widely-visible stellar and binary population synthesis programs such as StarTrack (Belczynski et al., 2008, hereafter B08) and MOBSE (Giacobbo et al., 2018; Di Carlo et al., 2019) have now adopted what can be called the state-of-the-art in stellar-wind and remnant-formation prescriptions (Belczynski et al., 2010; Fryer et al., 2012; Belczynski et al., 2016a, 2017b). On the other hand, much simpler and crude prescriptions are often used in works that explore the dynamical-formation channel, especially, those using the publicly-available and widely used direct N-body codes (Aarseth, 2003; Nitadori & Aarseth, 2012), (Aarseth, 2012), and (Wang et al., 2015). Note that such studies of the dynamical-formation channel, nevertheless, incorporate adequate or near-adequate treatment of stellar content, Newtonian dynamical interactions, and their GR counterparts (given that both the field-population-synthesis and the dynamical codes presently resort to similar, parametrized recipes for treating the tidally-interacting and the mass-transferring binaries’ internal evolution; Hurley et al. 2002). This, in turn, often makes model ingredients such as BH retention and mass distribution, that critically influence a cluster’s evolution and compact-binary formation, questionable in an otherwise realistic simulation. Note, in particular, that the BH population that is retained in a cluster following the BHs’ birth and which then becomes available for dynamical processing, influences the formation of both BH- and NS-containing compact binaries (Banerjee, 2018b), alongside governing the cluster’s long-term evolution (Chatterjee et al., 2017; Banerjee, 2017; Arca Sedda et al., 2018).
This work, for the first time, introduces three major upgrades to the widely-used direct N-body evolutionary program NBODY7 (Aarseth, 2012), a direct descendant of NBODY6 (Aarseth, 2003; Nitadori & Aarseth, 2012): (a) the semi-empirical stellar wind prescriptions as in Belczynski et al. (2010, hereafter B10), (b) remnant formation and material fallback according to the “rapid” and the “delayed” SN models of Fryer et al. (2012, hereafter F12), incorporating the occurrences of pair-instability supernova (PSN) and pulsation pair-instability supernova (PPSN) according to the conditions of Belczynski et al. (2016a, hereafter B16), and (c) an explicit modulation of the BHs’ and the NSs’ natal kicks based on the fallback fraction during their formation. The stellar- and binary-evolution drivers of NBODY6/7 (as well of the Monte Carlo-based star cluster-evolution codes MOCCA ; Hypki & Giersz 2013; Giersz et al. 2013 and CMC ; Joshi et al. 2000) are the well-known, fitting formulae and recipe based (as also for, e.g., StarTrack ) SSE (Hurley et al., 2000) and BSE (Hurley et al., 2002) programs 111Since SSE is a part of BSE , BSE will imply both of them hereafter in this paper. Similarly NBODY7 , that is utilized in this work, will, hereafter, imply both NBODY6 and NBODY7 since both adopt SSE/BSE .. After the introduction of the above updates (a), (b), and (c), BSE and NBODY7 are now similarly rich and diverse in essentially all stellar- and binary-evolutionary aspects as modern population-synthesis codes such as StarTrack and MOBSE . This would not only enable the “NBODY6/7 community” to do more realistic computations but also allow for direct, cross-community comparisons of results.
Note that the above-mentioned upgrades (a) and (c) are partly available in the current, public version of NBODY7 but the implementations do not produce outcomes fully consistent with StarTrack (see, e.g., Banerjee, 2017), as it should since BSE and StarTrack utilize the same underlying fitting formulae for the stellar-structural parameters (Hurley et al., 2000). In particular, the “maximum BH mass effect” (Belczynski et al., 2010) is absent in that implementation.
In this work, we validate the above new updates of by comparing its outcomes, for wide ranges of ZAMS mass and metallicity, with those from . We note that detailed comparisons among several binary-population-synthesis codes, in the context of the formation of white dwarf (hereafter WD) containing binaries, have recently been conducted by Toonen et al. (2014). In this work, we will, however, focus mostly on NS and BH outcomes since the above upgrades mostly affect SN outcomes.
In Secs. 2, 2.1, and 2.2, we describe our new updates to . In Sec. 2.3, we make detailed comparisons of the remnant outcomes of the updated with those from . In Secs. 3, 3.1, and 3.2, we describe how the updated is adapted to the direct N-body integration code , the new implementation of the explicit use of SN fallback fraction in determining SN natal kicks, and the new implementations of alternative natal kick models. In Sec. 4, we discuss how the single-stellar initial mass-final mass relations get affected by the presence of an observationally-motivated population of massive primordial binaries in a dynamically-active environment as in young massive and open clusters. In Sec. 4.1, we make a preliminary comparison between (isolated) massive binary evolutions from the updated and . In Sec. 5, we summarize our results and indicate forthcoming studies along this line.
2 The updated BSE : stellar-wind and remnant-formation schemes
In the following, the new implementations in the standalone BSE are described and the outcomes are compared with those of StarTrack .
2.1 The stellar-wind prescription
The recipe for stellar wind follows the semi-empirical prescriptions of B10, which is currently considered as the state-of-the-art. In the subroutine MLWIND, these prescriptions are implemented in such a way that they are activated individually in the following sequence of priority when the star is evolved until the second asymptotic giant branch (AGB) (BSE stellar type ; Hurley et al., 2000): (1) metallicity()-independent luminous-blue-variable (LBV) wind (Eqn. 8 of B10; Humphreys & Davidson, 1994), (2) -dependent Vink et al. (2001) winds (Eqns. 6 & 7 of B10), (3) -dependent Nieuwenhuijzen & de Jager (1990) wind (Eqn. 3 of B10), (4) the Hurley et al. (2000) treatments for lower-mass/colder stars (Kudritzki & Reimers 1978 wind for the giant branch and beyond and Vassiliadis & Wood 1993 wind for the AGB; Eqns. 1 & 2 of B10 respectively). For naked-Helium and more evolved stars (), only the -dependent Wolf-Rayet (WR) wind (Eqn. 9 of B10; Hamann & Koesterke, 1998; Vink & de Koter, 2005) is applied. The conditions for activating the various winds remain the same as in B10. The key difference between this new function and that of the original (or of the currently-public version of ) is that in the original (older) codes the Nieuwenhuijzen & de Jager (1990) and the Kudritzki & Reimers (1978) winds are applied for low-mass and high-mass stars. In new , for massive, hot stars, Vink et al. (2001) winds are applied instead. As it turns out, this makes a substantial difference in the resultant wind mass loss and hence in the remnant mass, especially for sub-solar metallicities.
2.2 The remnant-mass prescriptions
The remnant-mass model is a very important ingredient of the stellar-evolution modelling in and , that determines the masses of the NSs and the BHs formed as the final product of stellar (and binary; see Sec. 4) evolution. In other words, such a prescription, together with the underlying stellar-structure model and mass-loss recipes (see Hurley et al. 2000 and Secs. 1 & 2.1), determines the “initial-final” relation of stellar evolution. Although the final remnant mass can be expected to be an overall increasing function of the initial ZAMS mass, such initial-final curves often posses intricate excursions depending on the details of the remnant-mass prescription, the wind mass loss prescription, and the underlying stellar structure (see, e.g., B10, F12; below).
The subroutine HRDIAG in NBODY7 (which is directly ported also into the standalone BSE after omitting the limit from /star), that derives the remnant mass from the CO core mass in the AGB or the naked-Helium phases (in addition to dealing with the earlier evolutionary phases), adequately incorporates the older remnant-mass prescriptions provided by Belczynski et al. (2002), B08, and Eldridge & Tout (2004). The subroutine is now enhanced by incorporating the remnant-mass prescriptions of F12 for the rapid and the delayed SN (Eqns. 15-17 and Eqns. 18-20 of F12 respectively) 222Formula for in Eq. 16 (rapid explosions) of F12 should be: . Moreover, the PPSN and PSN conditions, as described in B16 (see references therein), are implemented in the HRDIAG routine. In the present version of the routine, any of the five remnant-mass schemes can be chosen with a switch (, 2, 3, 4, 5 for Belczynski et al. 2002, B08, F12-rapid, F12-delayed, and Eldridge & Tout 2004 recipes respectively) and the PPSN/PSN mass cutoff can be (de)activated independently ( or 1). The fallback amount and fraction onto a BH, and respectively, are calculated for each case according to their definitions in F12. As in previous studies, a 10% neutrino mass loss is assumed to obtain the gravitational mass of the remnant from its baryonic mass () in the case of a BH formation and Eqn. 13 of F12 (Lattimer & Yahil, 1989; Timmes et al., 1996) to account for the neutrino loss during an NS formation. The already-available electron-capture-supernova(ECS; Podsiadlowski et al., 2004)-NS formation scheme in , which is analogous to such scheme in B08 producing the characteristic NSs, is kept intact (activated with ). Hereafter in this text, ‘NS’ will imply a regular NS produced through core-collapse SN and an ECS product will be denoted by ‘ECS-NS’.
Both the new and routines, as described above, are now shared between the standalone BSE and NBODY7 .
All the amendments in the standalone code, to implement the new ingredients as described above and also in Secs. 3.1 and 3.2, are elaborated in Sec. A. Note that this updated , which is used in the rest of this work, retains the same binary-interaction physics (i.e., treatments of tidal evolution, mass transfer, common envelope evolution, and mergers) as described in Hurley et al. (2002) along with its subsequent amendments as available in the public version of .
2.3 Comparison with StarTrack
Fig. 1 demonstrates the agreement in ZAMS mass-remnant mass relation between the updated BSE and StarTrack for the B08333The ZAMS mass-remnant mass relations for the B08 case correspond to those of B10 which adopts the B08 remnant-mass model along with the newer wind recipes as described in Sec. 2.1 of this paper. This same wind model accompanies all the remnant-mass models considered in this work., F12-delayed, and F12-rapid remnant-mass models without PPSN/PSN and for , 0.006, and 0.02. Fig. 2 demonstrates similar agreements, for the F12-rapid case and the metallicities indicated therein, including the PPSN/PSN recipes of B16. The middle panel of Fig. 2 exhibits the PSN “mass gap” from BSE for and that beyond this gap ( along the ZAMS axis) the BH formation resumes.
The bottom panel of Fig. 2 shows such ZAMS mass-remnant mass relations (for the F12-rapid remnant mass model including B16-PPSN/PSN) up to very large ZAMS masses, beyond the PSN mass gap, for a wide range of (colour coding). It also demonstrates the excellent agreement between BSE (dashed lines) and StarTrack (solid lines) ZAMS mass-remnant mass curves for such large ZAMS masses. For each , the extent of the PSN mass gap (i.e., both the ZAMS mass at which PSN commences and the ZAMS mass, BH mass point at which BH formation resumes) is also in agreement. Note that with increasing , the PSN mass gap diminishes. The BH formation resumes whenever the mass of the He core, in the evolved phase (AGB and beyond) of the parent star, exceeds (see B16; Woosley 2017). However, at low , depending on the total wind mass loss, an amount of H-envelope retains (i.e., the star is of AGB type) so that the least massive pre-SN star, beyond the PSN mass gap, that would directly collapse to BH (see B16; for a pre-SN star of such a large mass, the direct collapse condition is satisfied irrespective of the remnant-mass model) may well exceed . For solar-like , the much stronger winds eliminate the H-envelope nearly completely (i.e., the pre-remnant star is of naked Helium type) so that the pre-SN (pre-collapse) star, just beyond the PSN mass gap, is of which, given the 10% neutrino mass loss (see B16; Sec. 2.2), results in an BH. Hence, the overall PSN BH mass gap, as obtained from , ranges between as indicated in Fig. 2 (bottom panel; black, horizontal lines) which limits closely agree with those of Woosley (2017).
Fig. 3 individually shows the close agreements between BSE and StarTrack ZAMS mass-remnant mass curves, until ZAMS mass, for the F12-rapid remnant mass model without (top panels) and with (bottom panels) PPSN/PSN and for , 0.002, and 0.0002. The “clipping” of the BH mass at , due to PPSN, is apparent in the column.
It is worth considering such comparisons between initial mass-final mass relations also for lower ZAMS masses which are WD progenitors. Fig. 4 demonstrates the agreement between the initial mass-final mass relations from and over the ZAMS mass range of . As evident, the outcomes from and agree nearly perfectly beginning from low-mass WDs to massive NSs.
2.3.1 The time-step parameters in BSE
Note that such near-perfect agreement of remnant masses between BSE and StarTrack depends on the choice of the stellar-evolutionary time step parameters , , and (Hurley et al., 2000) for the main sequence (MS) and the evolved phases respectively. These BSE input parameters (, , and respectively) are essentially fractions of stellar lifetimes in the main sequence, sub-giant, and more evolved phases that are taken as stellar-evolutionary time steps in the respective evolutionary stages. For a given underlying stellar structure (Hurley et al., 2000), the accuracy and hence the convergence of BSE ’s stellar-evolutionary calculation would be compromised if the time step parameters are chosen to be too large (on the other hand, choosing too small values would result in impractically long computing time and as well numerical difficulties).
This is demonstrated in Fig. 5, where the B08 remnant-formation scheme (that is already available in the to-date-public ) is applied. As can be seen, the spurious spikes in the BSE ZAMS mass-remnant mass curves (thin lines), that appear in the curves’ “maximum-mass” part (see B10) especially at low , decline with decreasing time-step parameters. In particular, , , , as defaulted in , still produce considerable spikes (see top-right panel of Fig. 5) for the lower s. The combination of , , is found to be optimal between speed and convergence over a wide range of (third row, right panel of Fig. 5) which is what is chosen in the BSE and NBODY7 computations in Figs. 1, 2, 3 and the rest of the figures in this work. The final panel of Fig. 5 demonstrates that despite taking very small time-step parameters, the currently-public version of produces ZAMS mass-remnant mass relations that overshoot their StarTrack counterparts significantly, especially for low and ZAMS masses. In particular, the “saturation” effect in BH mass is completely absent unlike the cases with the new MLWIND (with the same implementation of the B08 remnant scheme in HRDIAG ) and StarTrack whose outcomes agree nearly perfectly, as demonstrated above.
3 The adoption of the updated BSE in NBODY7 : remnant natal kicks
Although the new and routines can be readily shared between the standalone BSE and NBODY7 , the main stellar- and binary-evolution engines in the two codes are implemented in different ways, which may produce different remnant masses. Fig. 6 demonstrates as good agreements between the ZAMS mass-remnant mass relations, over a range of and remnant-mass prescriptions, when the new and are adopted in NBODY7 . The NBODY7 data in Fig. 6 are obtained whilst evolving model clusters of initial mass (K), initially-Plummer (1911) profile with half-mass radius pc, having the standard stellar initial mass function (IMF; Kroupa, 2001), and with all stars initially single. With the choices of the BSE time-step parameters suggested in Sec. 2.3.1 (which needs to be parametrized in the routine trdot in NBODY7 as opposed to specifying them as input parameters in the standalone BSE ), such runs practically do not slow down during the evolution of the most massive remnant progenitors (the first Myr), since the primary bottleneck on computing time (up to a given physical time) still comes from the direct N-body integration.
3.1 The retention of black holes in stellar clusters: standard,
fallback-controlled natal kicks
In the present context, the retention of stellar remnants, especially of BHs, in stellar clusters (young massive clusters, open clusters, and GCs) after their birth is widely debated. How many and what masses of BHs remain gravitationally bound to their parent cluster after their birth through a core-collapse SN, depending on their natal kicks, is instrumental in determining their long-term population evolution in the cluster, their impact on the structural and internal-kinematic evolution of the cluster (e.g.; Kremer et al., 2018b; Askar et al., 2018), and the nature of their dynamical pairing and GR merger (e.g.; Banerjee, 2017; Farr et al., 2017; Rodriguez et al., 2018; Samsing, 2018). The BH natal kicks are also very important for the formation of BBHs and BH-star systems and their coalescences through field binary evolution (e.g., Belczynski et al., 2016b; Stevenson et al., 2017). However, BH natal kicks are, to date, poorly constrained and understood from both observational and theoretical point of view (Willems et al., 2005; Fragos et al., 2009; Repetto et al., 2012; Repetto & Nelemans, 2015; Mandel, 2016; Belczynski et al., 2016c; Repetto et al., 2017). The single Galactic-field NSs (the vast majority of which would be the products of core-collapse SN), on the other hand, are observationally inferred to have large natal kick magnitudes, distributed according to a Maxwellian with 1-dimensional dispersion (average speed of ; Hobbs et al., 2005).
A commonly-used model (as in, e.g., Belczynski et al., 2008; Giacobbo et al., 2018) for core-collapse-SN natal kick magnitude, , is to assume NS-like kicks also for BHs (see, e.g., Repetto et al., 2012; Janka, 2013; Repetto & Nelemans, 2015) but which are scaled down linearly with increasing material-fallback fraction, (see Sec. 2.2), so that for (i.e., a failed SN or direct collapse) the natal kick is necessarily zero:
[TABLE]
where is chosen randomly from a Maxwellian of . For the model computations presented in this subsection, both BHs and NSs are treated with this kick scheme, except for the ECS-NSs which are given zero or natal kicks (Podsiadlowski et al., 2004; Gessner & Janka, 2018).
The , according to Eqn. 1, is applied in NBODY7 by transporting the , as computed in (see Sec. 2.2), to the standard version of the NBODY7 ’s subroutine via a dedicated common block. The routine already includes an elaborate algorithm for generating an isotropic distribution of velocity vectors with their magnitudes, s, chosen from a Maxwellian distribution of a specified dispersion. As in the original routine, the ECS-NSs are distinguished based on their , which are subjected to lowered-dispersion or zero s and are exempted from the above fallback treatment (the ECS engine is assumed to produce a full explosion). Keeping in mind that the newly-planted remnant-formation schemes in such as F12-rapid (Sec. 2.2) produce non-monotonic ZAMS mass-NS mass relations and, in particular, sub-Chandrasekhar NSs (see, e.g., Fig. 10), a very narrow mass window around is allowed for the ECS-NS treatment (unlike in the default where the ECS-NS treatment was invoked for , since the ECS-NSs were anyway the distinctly least massive NSs produced according to the Belczynski et al. 2002 and B08 prescriptions). Analogous updates are now implemented also in the standalone BSE ’s kick treatment which would then accordingly modify a binary’s response to an SN.
Fig. 7 shows the early (up to 20 Myr) evolution of the total mass, (top-left panel), and number, (bottom-left panel), of the BHs bound within the (initially all single stars) model cluster (see Sec. 3) for the various remnant-formation schemes (i.e., B08, F12-delayed, F12-rapid, F12-rapid+B16-PPSN/PSN; see Sec. 2.2), for this standard, fallback-controlled kick prescription, as given by Eqn. 1, when adopted in NBODY7 as above ( assumed) 444Hereafter, the subscript ‘bound’ indicates quantities measured within the cluster’s tidal radius, pc, as is customary. Most of the BHs within are indeed gravitationally bound to the cluster but, at a given time, a few of them may be on their way to escape the cluster. The number of the latter depends on the BHs’ and on the escape speed from the cluster ( here).. The differences in the curves between the remnant-formation cases arise due to the differences in which quantity determines the of the BH (or NS). The differences in the curves additionally arise due to the differences in the BHs’ masses (for a given ZAMS mass and ) between the various remnant-formation cases. Note that in all of the remnant-formation prescriptions considered here, all BHs and the heaviest NSs form with a nonzero . In the F12-rapid/delayed schemes, all NSs form with a small amount of fallback (fallback mass, ; see F12). However, for all NSs, is small () and the corresponding is typically large enough (but see Sec. 3.2) to eject them from a young massive, open, or globular cluster. The only NSs that would retain in clusters following their birth would be the ECS-NSs that receive small natal kicks (see above), in the present (and the following; Sec. 3.2) kick scheme(s).
As can be seen, in terms of retained BHs (but see below), the B08 remnant-mass model produces more BHs (in both mass and number) than any other remnant model, due to the wider “direct-collapse hump” in its ZAMS mass-remnant mass relations for lower ZAMS masses (; see Fig. 1 and B10), and hence expands a cluster more (see the bottom-right panel of Fig. 7; the expansion beyond Myr is driven mainly by the centrally-segregated BH subsystem). For the F12-rapid case, having PPSN/PSN or not makes a little difference in practice (as long as bulk properties are concerned, the PPSN truncation has important implications for interpreting BBH merger masses; see, e.g., B16, Mandel & Farmer 2017): due to the standard IMF in the model cluster, there are only a few stars undergoing PPSN/PSN. The dips in the and curves for the F12-rapid cases (top-left and bottom-left panels of Fig. 7) are due to a significant number of BHs receiving partial fallback, occurring in the dip between the two direct-collapse humps of the corresponding ZAMS-remnant curves (see Fig. 6), which escape marginally, taking time. For massive GCs and nuclear clusters, whose escape speeds, , are much higher (), this feature would vanish. The BH subsystems, for all cases, nevertheless sink and concentrate similarly for all the remnant-formation cases (top-right panel of Fig. 7).
Fig. 8 shows the BH mass distributions for the four remnant-formation scenarios as indicated in the panels’ legends, for . On each panel, both the mass distribution with which the BHs are born (steel blue histogram) and the mass distribution of the BHs retaining at Myr (blue histogram), in the cluster (standard IMF; initially single-only stars), are shown. At a time such as Myr, which is well after the completion of the BH formation (the last BH forms at Myr), all the unbound BHs (due to ) have escaped through the tidal radius but the remaining, bound BHs are still in the midst of segregating towards the cluster’s center (Fig. 7, top-right panel) and their dynamical processing is yet to begin. Therefore, the blue histograms fairly represent the mass distributions of the BHs which become available for long-term dynamical processing, for the various remnant-formation schemes (corresponding to with , standard IMF, ).
In fact, as seen in Fig. 8, the F12-delayed case inherently produces the largest number of BHs from stellar evolution (as opposed to the resulting retained BHs where B08 wins; see above), due the largest proto-remnant core + fallback masses (see F12) resulting in the lowest BH-(and NS-)formation threshold with respect to the core-collapsing progenitor mass, but it also produces the largest number of escapees due to its smallest direct-collapse hump, ultimately resulting in the smallest number of retainers (see Fig. 7). The F12-rapid cases produce less escapees (more direct- and near direct-collapse BHs) due to the “double hump” (see Fig. 6). The BH distribution is truncated at , as expected, when PPSN/PSN is included.
The discrepancy between the natal and the retained mass distribution, as in Fig. 8, gives the number and the mass retention fraction of BHs in a cluster, given the and the remnant-formation model. Table 1 provides such BH retention fractions, corresponding to the standard, fallback-controlled natal kicks (Eqn. 1), for the various remnant-formation models and s considered here.
3.2 The retention of black holes in stellar clusters: alternative natal-kick
prescriptions
While the fallback-modulated natal kick prescription has some observational and theoretical motivation (see Sec. 3.1 and references therein) and is often applied in N-body and population-synthesis studies, it is far from any robust constraint on observational or theoretical grounds, especially for BHs. Therefore, there exists ample room for alternatives of such prescriptions, in particular, for translating the NSs-like kicks to BHs that comprise a much wider remnant-mass range (see, e.g., Fig. 8). From the theoretical perspective, a key uncertainty lies in identifying the mechanism(s) that generate spatial and directional asymmetries in the SNs’ ejected (baryonic) material, which, in turn, momentum-propel the remnant. On the other hand, the asymmetry in the neutrino flux during an SN could alone be enough to kick an NS with the typical observed speed (Hobbs et al., 2005). If we assume one dominant SN-kick mechanism and take the observed NS natal kicks as the constraint, then we are already left with several possibilities as formulated below.
3.2.1 Convection-asymmetry-driven natal kick
The convection-asymmetry-driven kick mechanism assumes that the natal kicks are produced by asymmetries in the convection within the collapsing SN core (Scheck et al., 2004, 2008; Fryer & Young, 2007). This is in the line of the standard kick prescription but a more accurate prescription would be to boost the kick for systems with longer convection durations: a lower-mode convection produces stronger kicks and it develops with time, so the longer it takes to explode the stronger is the kick. In that case,
[TABLE]
Here, is a typical NS mass, is the remnant (NS or BH) mass, is the CO core mass, and is an efficiency factor (somewhere between 2-10) with the theory that we are more likely to get a low-mode convection in explosions that take longer to develop, i.e., those with larger .
3.2.2 Collapse-asymmetry-driven natal kick
Another alternative ejecta kick mechanism would be the one produced when asymmetries in the pre-collapse silicon and oxygen shells produce asymmetric explosions (Burrows & Hayes, 1996; Fryer, 2004; Meakin & Arnett, 2006, 2007). In contrast to the convection-asymmetry mechanism (see above), this mechanism would produce stronger kicks for shorter convective durations (smaller ) since the asymmetries would be washed out by long periods of convection:
[TABLE]
Here, is the suppression factor for larger s.
3.2.3 Neutrino-driven natal kick
The neutrino mechanism (Fuller et al., 2003; Fryer & Kusenko, 2006) would produce the kick through asymmetric neutrino emission. In that case, a formulation would be
[TABLE]
Here, () is the effective remnant mass beyond which the neutrino emission would not increase significantly with increasing mass of the SN core. Unlike the baryonic material, neutrinos do not fall back on to the remnant to return a part (or whole) of the ejecta momentum, despite the occurrence of such a fallback (or of a failed SN). Hence the omission of the term in Eqn. 4.
3.2.4 The dependence of BH retention on natal-kick mechanism
It would be interesting to study the impact of these different natal-kick prescriptions on the retention of BHs and NSs in a cluster and, potentially, formulate signatures to observationally distinguish between such cases. For this purpose, we continue to utilize the , pc model cluster (initially single-only stars; Sec. 3) here, which would facilitate comparisons. For definiteness, (Eqn. 2), (Eqn. 4), and are used.
Fig. 9 shows the early time evolution of , , and in the cluster for the different natal-kick prescriptions formulated in Secs. 3.1, 3.2.1, 3.2.2, and 3.2.3, as obtained from NBODY7 computations (). For implementing the alternative kick recipes in and , that can be selected with a switch, KMECH, in the routines, the corresponding is also transported from to via the dedicated common block (see Sec. 3.1).
As seen, the standard and the convection-asymmetry-driven cases collect very similar number and mass of BHs within the cluster while the collapse-asymmetry-driven case collects more BHs. This is expected since the latter kick model imparts much lower kicks (*c.f. *, Eqns. 1, 2, and 3) onto the BHs and also onto the relatively massive NSs. The neutrino mechanism does the worst in this regard with no BHs retained in the cluster: if the observed high NS kicks are indeed due to asymmetric neutrino emission alone then all BHs (provided the BH natal kick mechanism is the same as for NSs) would get substantially higher kick despite significant fallback (the emitted neutrinos would not return any momentum since they do not participate in fallback; see Sec. 3.2.3). If there is indeed, depending on dynamical age, a significant number () of stellar-mass BHs retaining and comprising a BH subsystem in several Galactic GCs (e.g., Askar et al., 2018; Kremer et al., 2018b, a), then, perhaps, the neutrino-driven kick mechanism is disfavoured.
The converges to a few after Myr after which ECS-NSs, due to their slow/zero dispersion (see Sec. 3.1), start to build up. However, there is a significant number of slow-escaping NSs’ buildup around 10-20 Myr age (see Fig. 9, bottom panel) which is unique for the case of the collapse-asymmetry-driven kick that may be interesting for radio observations in young massive clusters and for supporting or ruling-out such a kick scenario. Note that the model clusters utilized in these NBODY7 runs have which is 1/2 to 1/3 of typical s in massive GCs, so that one would potentially retain some BHs in massive, compact GCs in all the kick-mechanism cases anyway, although the overall trends would be similar as in Fig. 9. A survey of such GC models exploring these different natal-kick mechanisms, based on the Monte Carlo cluster evolution program MOCCA , is underway (Leveque, A., et al., in preparation).
Fig. 10 (left panel) compares, analogously to Fig. 8 (see Sec. 3.1), the mass spectra of the BHs retained in the parent cluster (at Myr), for the different natal-kick prescriptions considered here ( is assumed). As above, the standard and the convection-asymmetry-driven mechanisms produce similar retained BH mass spectra and nearly same BH retention in terms of number and mass. BHs escape the cluster at or shortly after birth. In contrast, the collapse-symmetry-driven mechanism produces low natal kicks for all BHs (see also Fig. 11) so that all of them retain in the cluster, following their natal mass spectrum (the dashed, magenta histogram of Fig. 10). Note that such complete retention of BHs for the collapse-asymmetry-driven kick mechanism is due to the low kicks for (Eqn. 3) which is below the BH-formation threshold. Therefore, this is a generic feature and it holds true for all remnant-mass schemes and for all metallicities. In particular, as seen in Fig. 10 (left panel), BHs would also be retained in GCs and open clusters; this particular feature is unique to the collapse-asymmetry-driven kick mechanism which can, therefore, be used to support or rule out such a kick scenario. Note, in this context, that the BH identified in the GC NGC 3201, through radial-velocity measurements, has a minimum mass of (Giesers et al., 2018).
The natal mass distribution of the NSs is also shown in Fig. 10 (right panel, black histogram). The prominent sub-Chandrasekhar peak is a feature of the F12-rapid scheme applied, convolved with the cluster’s standard IMF - these NSs are just the ( proto-remnant fallback) NSs (see F12). As mentioned before (Sec. 3.1), only the ECS-NSs retain in the cluster for longer term which is demonstrated in Fig. 10 (right panel, blue histogram). For massive GCs, some of the slow-escaping NSs would additionally retain in the collapse-asymmetry-driven-kick case (see above; Fig. 10, right panel, gray histogram).
Fig. 11 shows the s generated by NBODY7 , as a function of (left panel) and (right panel), for the different natal-kick recipes considered here (due to the logarithmic vertical axis, direct-collapse BHs with and are not shown in these panels). Comparing with the present model clusters’ typical (the solid, blue line), these panels clarify that for such clusters the standard and the convection-asymmetry-driven kick models would yield similar BH retention, the collapse-asymmetry-driven case would retain all BHs, whereas the neutrino-driven case would eject nearly all BHs. The ECS-NSs (see Sec. 3.1), which are given s from a Maxwellian distribution of much lowered dispersion (1-dimensional dispersion of ; see, e.g., Gessner & Janka 2018), make up the vertical stripe along the axis (Fig. 11, right panel). Most of these ECS-NSs retain in these model clusters due to their low natal kicks. From Fig. 11, it is also apparent that such conclusions concerning the retention of BHs and NSs in clusters remain nearly unaltered for , i.e., for young massive clusters, moderately-massive open clusters, and typical GCs. The mass gap between NSs and BHs (Fig. 11, right panel) is a characteristic of the F12-rapid remnant scheme adopted in this example.
4 The effect of massive primordial binaries
Until now, model clusters initially comprising only single stars are considered in this work. It would be worth looking at the case where a population of primordial binaries is present in the model, especially among the BH-progenitor stars. The merger among the members of a massive-stellar binary, especially during their advanced (non-remnant) stellar-evolutionary stages, can lead to a stellar structure (the relative masses of the CO-core, He-core, and H-envelope) that is not achievable through the evolution of a single star (even of the ZAMS mass that is equal to the sum of the ZAMS masses of the members; see below, Fig. 13) and hence can produce a (single) BH that deviates significantly from the corresponding single-star ZAMS mass-remnant mass relation (see, e.g., Marchant et al. 2016; Kruckow et al. 2018; Spera et al. 2019). Such merger events would, therefore, modify the natal and the retained BH mass spectra which would, in turn, modify the BH-driven long-term cluster evolution and the dynamical GR-merger events involving BHs. Note that in the dynamically-active environment of a dense stellar cluster, massive-stellar mergers can happen not only in the course of a massive-binary evolution (leading to Roche lobe overflow and common-envelope phases; see, e.g., Ivanova et al. 2013) but also during close binary-single and binary-binary interactions (e.g., in-orbit stellar collisions during Kozai oscillations in a dynamically-formed hierarchical triple or during a resonant triple-interaction; see, e.g., Banerjee et al. 2012; Leigh & Geller 2013). Even the outcomes of the massive-binary evolution channel can get significantly altered via close dynamical encounters that would alter the binaries’ orbital parameters. In principle, in such an environment, a given massive star can undergo multiple mergers until a BH is formed (e.g., Fujii & Portegies Zwart, 2013; Wang et al., 2020) or even a BH non-progenitor can grow in mass to become a BH progenitor. Finally, note that BH masses can get altered by not only such star-star mergers but also other stellar-hydrodynamical processes such as envelope ejection during a common-envelope (CE) phase and enhanced stellar wind in a tidally-interacting binary; such processes can also be either an outcome of massive-binary evolution or triggered dynamically (or both).
A detailed study of the influence of dynamical interactions on massive-binary evolution is beyond the scope of the present work. To have an idea of to what extent BH masses can get altered due to their progenitors being binary members, we evolve with , pc model clusters as before but with all stars (also sampled from a standard IMF) of ZAMS mass initially paired in binaries among themselves. A high fraction of massive-stellar (O-star) binaries (% initially in these models) is indeed consistent with observations of young massive clusters (Sana & Evans, 2011; Sana et al., 2013). Motivated by such observations, these O-type-stellar binaries are taken to initially follow the orbital-period distribution of Sana & Evans (2011) and a uniform mass-ratio distribution (an O-star is paired only with another O-star, as typically observed, and the pairing among the lower-mass stars is obtained separately; see below). The orbital periods of the non-O-star ( ZAMS) primordial binaries are taken to follow a Duquennoy & Mayor (1991) distribution that represents a dynamically-processed binary population (Kroupa, 1995) and their mass-ratio distribution is also taken to be uniform. The initial binary fraction among the non-O-type stars is taken to be small; %. The initial eccentricities of the O-type stellar binaries follow the Sana & Evans (2011) eccentricity distribution and those for the rest of the binaries are drawn from the thermal eccentricity distribution (Spitzer, 1987)555We note the recent work by Geller et al. (2019) questioning the validity of the assumption of an initial thermal eccentricity distribution of the primordial binaries. The sub-population of primordial binaries that directly influences the present results is the O-star-binary population for which the adopted (sub-thermal) initial Sana & Evans (2011) eccentricity distribution is motivated by observations.. As explained in Banerjee (2018b), such a scheme for including primordial binaries provides a reasonable compromise between the economy of computing and consistencies with observations.
The filled, black squares in Fig. 12 (top row) show the remnant mass versus ZAMS mass for the above model cluster with massive primordial binaries, when evolved with that includes the new (Sec. 2); for the BH progenitors that have undergone a star-star merger before the BH formation, the ZAMS mass of the primary (the more massive of the members participating in the star-star merger, at the time of the merger) is plotted in the abscissa. Due to BH progenitors undergoing mergers in their evolved phases (and/or due to other stellar-hydrodynamical processes) before the BH formation (see above) 666The BH masses are printed right at their formation by arranging a special output in the new routine (Secs. 3.1, 3.2), as also for the data points in Fig. 6. Hence, any merger(s) or stellar-hydrodynamical process(es), influencing an BH mass in Fig. 12, must have happened before the BH formation. Otherwise, the BH mass will necessarily lie on the ZAMS mass-remnant mass relation as in Fig. 6. In particular, the post-formation mass gain of a BH (due to, e.g., mass transfer, BH-star collision) is not captured in Fig. 12., the resulting BH masses deviate significantly from the ZAMS mass-remnant mass relations (which relations, by definition, are for single ZAMS stars) for the assumed remnant-formation scheme (F12-rapid with and without PPSN/PSN) and metallicity (). The majority of the BHs are over-massive, showing broader natal and retained BH mass distributions as contrasted with those from single stars only, for models with and without PPSN/PSN physics included (Fig. 12, second row; the standard, fallback-controlled natal kick is assumed). Also, due to the formation of more massive BHs, the total mass of the BHs formed/retained in the cluster is similar to that of the initially-single-star-only model clusters (Fig. 12, third row) despite the BHs forming/retaining in a smaller number (due to the mergers; Fig. 12, fourth row). These dependencies, as demonstrated in Fig. 12, suggest that BH retention in a cluster and the retained BHs’ mass spectrum both depend on the BHs’ progenitor stars’ binary properties. Table 1 demonstrates the somewhat increased mass and number retention fractions of BHs, at Z=0.0001, in the presence of a massive primordial-binary population as adopted here (*c.f. *the top and the bottom sections of the table). Note that although complex, dynamical-interaction-mediated channels leading to star-star mergers (see above) among members from different binaries can and do occur in these models, most of these mergers are between the members of the same O-star primordial binary occurring due to Roche lobe overflow (the overflow itself can, of course, be triggered via dynamical perturbation of the binary’s orbit) 777Some authors distinguish between “coalescence” occurring due to mass transfer in a binary and “merger” due to an in-orbit or a hyperbolic collision. In this paper, we will call the event a “merger”, irrespective of the channel, whenever two stars combine to become a single entity..
The BH masses in Fig. 12 should, however, be taken with caution. The BH masses depend on how the evolved stars are “mixed” at their merger to compose the stellar collision product. This is rather simplistic in where a complete mixing is assumed to compute the stellar-evolutionary age and hence the CO-core and He-core masses of the merged star. A simple recipe for the “dynamical” mass loss during the merger process is assumed which results in a mass loss of % of the total available mass budget per merger (e.g., Glebbeek et al., 2009) 888In the version of BSE adopted in NBODY7 , a 30% mass loss from the secondary (the less massive merger participant) is assumed when the kinetic energy of approach of the merging non-compact stars exceeds the internal binding energy of the secondary, provided the stars geometrically overlap at the merger and provided they are in the main sequence. Merger among more evolved stars and grazing encounters are even poorly understood/studied and no mass loss is associated with such mergers in . If one of the merging member is an NS or a BH, an unstable Thorne-Żytkow object is formed. In such a case, the entire secondary is assumed to be dissipated in the case of an NS core and half of the secondary mass is retained onto the BH in the case of a BH core.. In reality, a star-star merger is likely to be a much more complex process whose study is beyond the scope of the present work. Nevertheless, this study suggests that the natal/retained BH mass function in clusters can potentially be substantially modified and broadened, despite PPSN/PSN, provided the clusters contain a high fraction of massive primordial binaries at their young age, as observed in young massive and open clusters.
This is due to the “unusual” structure obtained in merger products that contain a much more massive H-envelope when the He-core is below the PPSN threshold of (that see B16), compared to that of a single-stellar-evolution product of the same He-core mass. The substantial H-envelope in the merger products results in BHs (formed via direct collapse) well exceeding despite the sub-PPSN He-core. Note that the BH mass surpasses the PPSN cutoff (which mass would result from the collapse of a He-core after subtracting for the 10% neutrino mass loss) for ZAMS mass approaching from below (resulting in the He-core approaching the PPSN threshold of ) also for pure single-star evolution, but only by a few (see Figs. 2 & 6). Owing to much more massive H-envelope for similar He-core masses, this surpassing happens to a much larger extent for the BHs obtained from merger products that result from the observationally-motivated massive binary population, as adopted here, in a dynamically-active environment. However, this still retains the notion that BBH mergers comprising one or both of the members having (the single-stellar PPSN cutoff limit) are necessarily of dynamical origin (see, e.g., Rodriguez et al., 2018). This is because after being formed from a star-star merger product, a BH can continue to interact with stellar members and remnants and participate in a GR merger only in a dynamically-active environment. However, the occurrence of BBH mergers within the “PSN mass gap” in this way would not require low natal spins of (first-generation) BHs. On the other hand, if stellar progenitors indeed produce very low spin BHs irrespective of their own spin, internal structure, and chemical composition (Spruit, 2002; Fuller et al., 2019), then, again, in contrast to the widely-conceived notion of second-generation BHs having high spins, such mass-gap BBH mergers would exhibit a low effective spin parameter, provided the first-generation, mass-gap BH merges with another first-generation BH. See also Spera et al. (2019) in this context.
On the other hand, when the He-core mass exceeds the PPSN threshold of and until the end of the PSN gap is reached (He-core mass of ; see Sec. 2.3), the outcome becomes independent of the mass of the H-envelope around it. Namely, when the He-core mass lies between (the PPSN range) a BH of is born via direct collapse as a result of PPSN and for He-core between no remnant is left due to PSN (the BH formation is resumed and the BH mass will again depend on the H-envelope mass if the He-core exceeds ; see Sec. 2.3 and Fig. 2). Note that in the present BSE scheme a stellar object is always treated in the same, single-star way (see Sec. 2.2) irrespective of it being a merger product or not. Hence, even for a merger product with unusually-massive H-envelope, a BH (no BH) is formed whenever the PPSN (PSN) condition is reached, such conditions and the corresponding outcomes being dependent solely on the He-core mass. This is why in Fig. 12 (top row, right panel, filled, black squares), there are only or no BHs beyond the start of the horizontal, single-star PPSN line (at ZAMS mass for ) even in the presence of massive primordial binaries: a merger, which typically leads to a similarly- or more-massive He-core (depending on the mass ratio and the epochs of the merger members), will exclusively satisfy either the PPSN or the PSN condition when a star of ZAMS mass is involved. In this example, all but one of the PPSN derived BHs and PSN points, beyond ZAMS, come from star-star merger products. For reference, Fig. 13 re-plots the top row of Fig. 12 with the sum of the ZAMS masses of the star-star mergers’ members along the abscissa, for the primordial-binary model. Here, of course, the deviations from the single-star ZAMS mass-remnant mass relations are less dramatic, especially for large combined ZAMS masses, and not all BHs beyond combined ZAMS are formed via PPSN.
Note that when a much more massive H-envelope, compared to what would be available from single stellar evolution, is present due to, say, a merger as considered here, a PPSN may shed only a part of the envelope resulting in a BH that is significantly more massive than what the collapse of the He-core alone (i.e., if the envelope is entirely expelled) would give. The PPSN/PSN recipe of B16 (see above and Sec. 2.2), as implemented in the current StarTrack and in the new BSE introduced in this work and which is applied also for star-star merger products, assumes complete stripping of the H-envelope, always resulting in a BH. This is in line with the hydrodynamic study of PPSN/PSN in single massive stars by Woosley (2017) where a small remaining H-envelope results in BHs up to as outcomes of PPSN in low-, massive stars. More massive BHs would potentially result from merger products where the PPSN is triggered in the presence of a substantially more massive H-rich envelope. This would, in turn, diminish the PSN mass gap (see Sec. 2.3 and references therein, Fig. 2). Also, due to the enhanced pressure from such an envelope, PPSN may be triggered “prematurely” when the He-core is of resulting in a longer PPSN phase with a larger number of mass loss episodes. At present, such stellar models are unavailable. Although beyond the scope of the present study, a more elaborate study of this regime will be presented elsewhere.
Of course, the extent of the exceeding of the single-star PPSN cutoff depends on the amount of “dynamical” mass loss during the star-star merger process. The physical process of merger among stars and the associated mass loss is poorly understood; especially so for massive stars beyond the main sequence. At present, the mass loss from star-star mergers in is at most moderate, ranging from % (energetic, equal-mass MS-MS mergers) to zero (non-MS or low-K.E. mergers; see above), the MS-MS merger mass loss being consistent with basic hydrodynamic (e.g., Lombardi et al., 2002) and “sticky-star” (e.g., Gaburov et al., 2008) treatments of the process. However, if the mass loss per merger is much higher, say, comparable to the mass of the secondary in all mergers, then the H-envelope in the merger products would be similarly massive as in single-stellar evolution. In such a case, the remnant mass-ZAMS mass relations and the PSN mass gaps would be much less affected by massive-stellar mergers from primordial binaries.
Fig. 14 demonstrates this for the case of F12-rapid remnant-formation including PPSN/PSN and low metallicity (). The panels in Fig. 14 show the remnant mass-ZAMS mass outcomes (black, filled squares) from different calculations of the primordial-binary model with increasing fraction, , of the secondary’s (the less massive merger participant’s) mass lost during a merger: from in energetic MS-MS mergers only ( default; left panel), through in all star-star mergers (middle panel), to in all star-star mergers (right panel) 999These enhanced s are implemented in NBODY7 for demonstration purposes only. It would, perhaps, be more realistic to have at least a mass condition for the value of which will be implemented in the near future. 101010The filled, black points won’t always map through the panels in Fig. 14 since varying alters the local gravitational potential and hence the dynamical interactions, so that a particular merger does not necessarily re-occur in the different N-body computations and the sequence of the mergers is also altered over the different runs.. Fig. 14 suggests that despite mergers of massive primordial binaries, (% loss of the total primary + secondary mass at the time of the merger; see de Mink et al. 2013 in this context) already restores the PPSN/PSN plateauing of BH mass as obtained from the latest core-collapse calculations (maximum BH mass of ; Woosley 2017) and (% loss of the total available mass budget) nearly restores the entire single-star remnant mass-ZAMS mass relation. In other words, if the single-star remnant mass-ZAMS mass relation is to be unaffected by massive primordial-binary mergers, then a substantial, up to %, mass loss would be necessary in such mergers.
4.1 Comparison with binary evolution
In Sec. 2.3, the comparison between the remnant outcomes of and were limited to single stellar evolution. In the present context, it would be worth comparing between the outcomes of binary evolution with and the updated . Such a comparison (as such, between any two binary-evolution programs) is a more complex affair than comparing single-stellar evolution due to the differences in the treatments of binary-interaction physics, namely, tidal interaction, mass transfer, CE evolution and as well in the treatments of star-star and star-remnant merger products. The difficulty of such a comparison exercise is augmented by the fact that these binary-interaction processes are themselves poorly understood (Ivanova et al., 2013).
At the moment, it is formidable to do a comparison while integrated within as done above for alone. Integration of within is far from being straightforward partly due to the difference in programming languages (C in versus Fortran in and ) and partly due to the rather rigid structure of ’s subroutines that connect to its N-body integration engine. A foreseeable solution would be to use a bridging framework such as AMUSE (Portegies Zwart & McMillan, 2018) but, currently, neither nor is adequately adapted for use in AMUSE.
Therefore, to preliminarily assess the differences in the outcomes of and in the context of remnant formation, which is the focus in this work, we evolve 12 selected massive binaries with updated and . The individual binary parameters and the corresponding remnant outcomes with and are summarized in Table 2. For both and evolutions, metallicity , CE efficiency parameter , F12-rapid+B16-PPSN/PSN remnant model, and star-star and star-BH merger mass loss fraction of of the secondary (less massive) member are applied. The mass transfer is set to be Eddington limited (Eddington limiting factor ).
The main difference between these and outcomes is that leads to more star-star mergers (especially during CE) than , leading to a single BH remnant. As evident from Table 2, all but one binaries yield just one remnant in StarTrack (including one no-remnant outcome due to PSN). In contrast, through evolution, 8 out of the 12 binaries lead to double remnants, 3 of which are BBHs and 2 of the BBHs merge within a short delay time.
The main cause for difference in binary-evolutionary outcomes between and , in these examples, is due to differences in the treatment of tidal energy loss among the interacting binary members. In particular, as described in Belczynski et al. (2008), employs 50 times stronger tidal energy loss than (Hurley et al., 2002) for convective stellar envelopes and the energy dissipation is weaker for radiative envelopes due to the use of tidal constants of Claret (2007). These differences cause the initial, post-SN, and post-CE-ejection detached (and eccentric) binaries to shrink (and circularize) at different rates, leading to subsequent binary interactions at different stellar-evolutionary ages. In particular, the more efficient tidal dissipation for convective envelopes would lead to faster orbit shrinking (during spin-orbit synchronization), leading to Roche lobe filling at a younger age when the member stars have smaller radii, thus leading to mergers in more of the systems. At the same time, for radiative envelopes, the reduced tidal efficiency leads to mergers by slowing down the orbital circularization, enabling Roche lobe overflow near the periastron at a younger age when the stars have swelled less.
Interestingly, if one considers only the most massive remnant produced from each binary evolution, they agree reasonably well among and and very closely for 7 out of the 12 binaries. This is demonstrated in Fig. 15; see also Table 2. The differences between and remnant masses are less visible in this way since, for most of the evolutions, the mass of the heavier among the two remnants or of the single remnant is determined mainly by the PPSN mass ceiling and star-star and star-BH merger mass loss which are taken to be identical for and , in this exercise.
The difference in treatments of tidal interaction and also in implementations and numerical approaches would lead to differences in outcomes of low-mass interacting binaries as well. Such a study is beyond the present scope; see Toonen et al. (2014) in this context.
5 Summary and outlook
The work presented above is summarized below where the next steps are also motivated.
- •
A reimplementation of the B10 wind model (Sec. 2.1) and the adoption of newer remnant-mass prescriptions such as F12-rapid and F12-delayed including the PPSN/PSN conditions of B16 (Sec. 2.2) in the widely-used population synthesis program reproduces the ZAMS mass-remnant mass relations of nearly perfectly over a wide range of ZAMS mass, remnant-formation scheme, and metallicity. (Sec. 2.3, Figs. 1, 2, 3). Such agreements, however, require smaller stellar-evolutionary time step parameters (the , , and parameters) than their default values while running (Sec. 2.3.1, Fig. 5), that slightly slows down the program. The new , therefore, now incorporates remnant-mass and stellar-wind prescriptions, in terms of detail, flexibility and the outcomes, at par with state-of-the-art population-synthesis programs such as and .
- •
These new developments in naturally get ported into the direct N-body code since the latter adopts as its stellar-evolution and remnant-formation driver. The adoption of the new routines (the functions and ) in yields a similar quality of reproduction of the ZAMS mass-remnant mass relations, for NSs and BHs formed from single stars within model star clusters (Sec. 3, Fig. 6).
- •
To deal with the retention of NSs and BHs in such model clusters at their birth, an explicit fallback modulation is implemented in the natal-kick generator of , which can be directly ported into the corresponding implementation in the standalone program (their routines; Sec. 3.1). In addition to this standard, fallback-controlled natal kick prescription, its variants describing convection-asymmetry-driven and collapse-asymmetry-driven natal kicks are prescribed and as well a neutrino-driven kick model (Sec. 3.2; Fig. 11).
The above implementations have now been tested to work very stably in an extensive set of long-term runs (Banerjee, S., in prep.). These new stellar-evolutionary ingredients and natal-kick models have now as well been ported to and tested to perform as stably (in prep.).
The various remnant and natal-kick prescriptions in the routines and can, at present, be selected based on switches built in these routines (for , also the stellar-evolutionary time step parameters), requiring recompilation of and . In the near future, arrangements will be made to provide these switches in the inputs, to make the options more accessible. Such a tested version of standalone BSE , with all the newly-implemented prescriptions described here, can be obtained freely. 111111https://tinyurl.com/v7uadf4
- •
For a given cluster (and given IMF), the amount of BHs retaining (i.e., remaining gravitationally bound to the cluster) due to low-enough natal kicks to become available for long-term dynamical processing depends on the remnant-formation mechanism (Fig. 7), the natal-kick mechanism (Fig. 9), and metallicity (Table 1). For all remnant-formation and natal-kick schemes (except for the collapse-asymmetry-driven case; see below), the retained BHs are generally top-heavier in mass spectrum compared to their birth mass spectrum (Figs. 8 & 10) due to the preferential escape of lower-mass BHs as a result of their higher natal kicks (or smaller fallback mass; Fig. 11). The F12-delayed remnant-formation model gives birth to the largest number of BHs (as opposed to the resulting retained BHs where B08 maximizes) but it also produces the largest number of escapees ultimately resulting in the smallest number of retainers.
- •
Table 1 gives the BH mass and number retention fractions in a moderately-massive stellar cluster as a function of remnant-formation scheme and metallicity, for the standard, fallback-controlled natal kick. The convection-asymmetry-driven kick mechanism would result in nearly the same retention fractions, the collapse-asymmetry-driven mechanism would lead to all retention fractions and the neutrino-driven mechanism would lead to all retention fractions (Fig. 9, 11). Hence, if BHs still retain in some present-day GCs as inferred when their observed photometric and internal-kinematic structures are compared with model computations (e.g., Askar et al., 2018; Kremer et al., 2018b, a), then the neutrino-driven mechanism can possibly be ruled out. If forthcoming observations suggest the presence of BHs in GCs and/or open clusters and/or young massive clusters, then this would provide a support to the collapse-asymmetry mechanism (Sec. 3.2.4, Fig. 10). For the case of the latter kick mechanism, a population of slow-escaping core-collapse NSs can also be expected in the field of view of Myr-aged young clusters (Sec. 3.2.4, Fig. 9), which may be identified in the future through radio observations. The BH retention fractions in Table 1 represent those in a moderate-escape-velocity system, essentially due to direct collapse or a substantial amount of material fallback alone: for massive GCs, their higher escape velocity would additionally aid retention, potentially resulting in somewhat higher retention fractions in all cases. These inferences encourage future, more elaborate studies of discriminatory signatures of natal-kick mechanisms in stellar clusters, especially in more massive systems (Leveque, A., et al., in preparation).
- •
If BH-progenitor stars exist in high primordial-binary fractions, as observed in young massive and open clusters, then, depending on the mass loss during star-star merger process, the BHs’ natal and cluster-retained mass spectra would become wider; at low metallicities, the BH mass can exceed well beyond even in the presence of PPSN, forming BHs bound to a stellar cluster well within the PSN mass gap (Sec. 4, Fig. 12). Since such BHs are still direct derivatives of stellar progenitors, they would potentially have very low spins (unless the star-star merger happens late in the evolution of both stars). A high primordial-binary fraction among BH progenitors would also lead to a moderately higher BH-retention fraction in both mass and number (Table 1) although the total number of BHs produced, under such circumstance, would be somewhat smaller due to the mergers (Fig. 12). The BH masses get altered due to pre-BH-formation stellar mergers and/or other stellar-hydrodynamical processes (e.g., CE envelope ejection, tidally-enhanced stellar wind). In this work, these modified BH masses are the outcomes of the simplistic treatment of these complex physical processes in the binary-evolution engine of and as well of the single-star-like treatment while forming remnant masses from stellar-merger products that would have altered stellar structure, which is why their values should be taken with caution.
The above effect, if real, would still preserve the (potentially-future) LIGO-detected BBH mergers comprising members as signatures of dynamically-triggered BBH mergers in clusters, except that such inferences need not, anymore, depend on low natal spins of BHs (see, e.g., Rodriguez et al., 2018). If stellar progenitors indeed produce low-spinning BHs irrespective of their own spin (Belczynski et al., 2017a), internal structure, and chemical composition, then such mass-gap BBH mergers could potentially exhibit low effective spin parameter. Whether a dynamically-active environment such as that inside a stellar cluster aid in inducing star-star mergers that lead to low-spin, PSN-mass-gap BHs, will be investigated in a separate study. If the single-star remnant mass-ZAMS mass relation is to be unaffected by massive primordial-binary mergers, then a substantial, up to %, mass loss would be necessary in such mergers (Sec. 4, Fig. 14). In this context, it would be worth bearing in mind that the PPSN cutoff mass is not strictly constrained (see B16). In particular, the recent measurement of a primary in the BBH merger event GW170729 (Abbott et al. 2019b, a; but see Fishbach et al. 2019) may call for a moderate revision of the PPSN BH mass limit in population-synthesis codes to increase it by as was already proposed in Belczynski et al. (2017a, see their Fig. 3). In fact, considering the full mass spectrum of BBHs detected to date by LIGO-Virgo (Abbott et al., 2019b), GW170729 is statistically consistent with being a first-generation merger (Chatziioannou et al., 2019; Kimball et al., 2020). Also, a more elaborate study of the influence of a dynamically-active environment on massive-binary evolution will be taken up in the near future.
Acknowledgements.
We thank the anonymous referee for their useful comments and suggestions which have improved the manuscript. SB acknowledges the support from the Deutsche Forschungsgemeinschaft (DFG; German Research Foundation) through the individual research grant “The dynamics of stellar-mass black holes in dense stellar systems and their role in gravitational-wave generation” (BA 4281/6-1; PI: S. Banerjee). SB acknowledges the pleasant hospitality of the Nicolaus Copernicus Astronomical Centre of the Polish Academy of Sciences (CAMK) where this work was initiated. SB acknowledges the generous support and efficient system maintenance of the computing teams at the Argelander-Institut für Astronomie (AIfA) and the CAMK where the computations have been performed. KB acknowledges support from the Polish National Science Center (NCN) grant Maestro (2018/30/A/ST9/00050). The work of CF was done at Los Alamos National Laboratory under project number 20190021DR. This work was supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – Project-ID 138713538 – SFB 881 (“The Milky Way System”, subproject A06) and by the Volkswagen Foundation under the Trilateral Partnerships grants No. 90411 and 97778. The special GPU accelerated supercomputer Laohu at NAOC has been used and we thank the Center of Information and Computing of NAOC for support. This work has been benefited from support by the International Space Science Institute, Bern, Switzerland, through its International Team program ref. no. 393 “The Evolution of Rich Stellar Populations and BH Binaries” (2017–18). The work of PB was also partially supported under the special program of the NAS of Ukraine “Support for the development of priority fields of scientific research” (CPCEL 6541230). RS and PB acknowledge support through the Silk Road Project at the National Astronomical Observatories (NAOC), Chinese Academy of Sciences, through the National Natural Science Foundation of China (NSFC, grant No. 11673032), and through the Strategic Priority Research Program (Pilot B) “Multi-wavelength gravitational wave universe” of the Chinese Academy of Sciences (No. XDB23040100). Visits of SB and KB have been supported by the Silk Road Project at NAOC. JH has been supported by the Kavli Visitor Program of Kavli Institute for Astronomy and Astrophysics, Peking University. LW acknowledges the support from Alexander von Humboldt Foundation and JSPS International Research Fellow (School of Science, The university of Tokyo). We thank Mirek Giersz and Diogo Belloni for pointing out an issue with white dwarf formation and white dwarf masses under our newly implemented schemes which has helped us to locate and fix an important bug in one of the routines in the updated code presented in this work. The “time step issue” in (Sec. 2.3.1) was hinted by Mirek Giersz during a lunch at CAMK.
Appendix A Summary of the changes implemented in the updated code
Below, we summarize the changes in the various subroutines of the updated described in this work with respect to the current publicly-available version of this program.
mlwind.f
Significant reorganization of the original mlwind function is done to implement the B10 wind recipe correctly; see Sec. 2.1 for the details. The original Hurley et al. (2000) wind and various variants of both B10 and Hurley et al. (2000) wind recipes can be activated by appropriate choice of the mdflag flag variable. mdflag selects the B10 wind which option is now recommended.
hrdiag.f
In the subroutine hrdiag, the “rapid” and “delayed” remnant-mass prescriptions of F12 are implemented utilizing the relevant analytical formulae; see Sec. 2.2 for the details. These formulae are inserted in the same parts of the code where the final baryonic remnant mass is evaluated from the older prescriptions. All the older remnant prescriptions are preserved and can be selected via choices of the input flag variable nsflag: currently nsflag (4) selects the F12-rapid (delayed) prescriptions.
The decision for undergoing PSN (giving a remnant of zero mass) or PPSN truncation according to the criteria in B16 is taken before the baryonic remnant mass evaluation and is implemented as an extension of the existing conditional statement that handles the cases of no remnant formation. PSN/PPSN can be enabled (disabled) with the choice of the input flag variable psflag.
Care is taken that these extensions are analogously coded in both the hydrogen rich-star and helium-star sections of hrdiag, as for the other remnant prescriptions.
Apart from these, two additional COMMON BLOCKs are introduced. One is
COMMON /FBACK/ FBFAC,FBTOT,MCO,ECS,
which enables the transport of the fallback fraction, fallback mass, CO core mass, and an indicator for the occurrence of an ECS to other subroutines where they would be useful such as kick (see below). The other is
COMMON /FLAGS2/ psflag,kmech,ecflag
to port relevant new input flags (ecflag for enabling ECS is now read from input just for symmetry and convenience).
Note that the ECS indicator, which is set to unity directly in the conditional statement in hrdiag that assigns the ECS-NS mass, is implemented to override the default approach of sensing ECS events (outside hrdiag) via NSs. This is an arrangement to address the fact that in the newly-implemented F12-rapid scheme, the core-collapse NSs’ mass can go down to so that a NS can as well be a core-collapse NS.
kick.f
The analytical formulae for modelling the various core-collapse SN natal kick mechanisms, as explored in this work (see Secs. 3.1 & 3.2), are all implemented in the subroutine kick. To evaluate the kick-velocity magnitude, the SN fallback fraction and the CO core mass are needed which are supplied by the subroutine hrdiag through the FBACK common block (see above). kick also uses the FLAGS2 common block to receive the value of the input flag variable kmech that toggles among the kick models.
Note that unlike mlwind and hrdiag, the kick subroutine is structurally quite different in , which also has to keep track of the identities of the members receiving the kick. The natal kick formulae have been analogously implemented in ’s version of kick. This subroutine is currently private and will be made public elsewhere.
star.f
The subroutine star differs slightly from its public version and is now identical with star.f that is supplied with the public version of .
sse.f
An extra READ statement is added in this main routine to read in the choices of the input flags psflag, kmech, and ecflag (see above).
bse.f
The same as in sse.f.
const_bse.h
In this header file, the common block
COMMON /FLAGS2/ psflag,kmech,ecflag
is added to carry the user-input values of psflag, kmech, and ecflag to other subroutines.
README_NEW
The descriptions of the new ingredients, including references, are elaborated in this companion manual. The functionalities of the new input flags (see above) and the amended formats of the runtime input files are also described here.
Stellar mass limit
Apart from the changes described above, the imposed stellar mass limit is, by default, commented out in star.f, hrdiag.f, mix.f, and evolv2.f. This facilitates the user to explore the very-massive-star regime where PPSN and PSN would occur and also the outcome of massive-stellar mergers. However, it should be borne in mind that extending beyond is a rather large extrapolation of the underlying stellar structure employed in .
The stellar-mass ceiling is by default suppressed in the public .
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Aarseth (2003) Aarseth, S. J. 2003, Gravitational N-Body Simulations, Cambridge University Press, Cambridge, UK, pp. 430. ISBN 0521432723, 430
- 2Aarseth (2012) Aarseth, S. J. 2012, MNRAS, 422, 841
- 3Abadie et al. (2010) Abadie, J., Abbott, B. P., Abbott, R., et al. 2010, Classical and Quantum Gravity, 27, 173001
- 4Abbott et al. (2016 a) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2016 a, Ap J, 818, L 22
- 5Abbott et al. (2016 b) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2016 b, Physical Review Letters, 116, 061102
- 6Abbott et al. (2019 a) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2019 a, Ap J, 882, L 24
- 7Abbott et al. (2019 b) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2019 b, Physical Review X, 9, 031040
- 8Abbott et al. (2017 a) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017 a, Physical Review Letters, 118, 221101