Populations of stellar mass Black holes from binary systems
Grzegorz Wiktorowicz Łukasz Wyrzykowski Martyna Chruslinska Jakub Klencki Krzysztof A. Rybicki Krzysztof Belczynski
astro-ph.HE


Populations of stellar mass Black holes]Populations of stellar mass Black holes from binary systems



   
     ^1 National Astronomical Observatories, Chinese Academy of Sciences, Beijing 100101, China

     ^2 School of Astronomy & Space Science, University of the Chinese Academy of Sciences, Beijing 100012, China

     ^3 Astronomical Observatory, Warsaw University, Al. Ujazdowskie 4, 00-478 Warsaw, Poland

     ^4 Institute of Mathematics, Astrophysics and Particle Physics, Radboud University Nijmegen, PO Box 9010, 6500 GL

     ^5 Nicolaus Copernicus Astronomical Center, Polish Academy of Sciences, Bartycka 18, 00-716 Warsaw, Poland 

            

  

	In large and complicated stellar systems like galaxies it is difficult to predict the number and characteristics of a black hole (BH) population. Such populations may be modelled as an aggregation of homogeneous (i.e. having uniform star formation history (SFH) and the same initial chemical composition) stellar populations. Using realistic evolutionary models we predict the abundances and properties of BHs formed from binaries in these environments. We show that the BH population will be dominated by single BHs (SBH) originating from binary disruptions and stellar mergers. Furthermore, we discuss how BH populations are influenced by such factors as initial parameters, metallicity, initial mass function (IMF), and natal kick (NK) models. As an example application of our results, we estimate that about 26 microlensing events to happen every year in the direction of the Galactic Bulge due to BHs in a survey like OGLE-IV. Our results may be used to perform in-depth studies related to realistic BH populations, e.g. observational predictions for space survey missions like Gaia, or Einstein Probe.	We prepared a publicly available database with the raw data from our simulations to be used for more in-depth studies.






§ INTRODUCTION


A BH is defined as a region in space from which nothing, even light, can escape <cit.>. BHs may be detected when interacting with other objects <cit.>, as gravitational wave sources <cit.>, when they interfere with radiation <cit.>, or in non-interacting binaries by observations of their companions <cit.>. It was also proposed that SBHs may be observed as X-ray novae <cit.>.

The so called stellar-mass BHs, i.e. stellar-origin BHs and BHs originating from mergers of stars and/or stellar-origin compact objects, with masses from ∼5 to possibly a few ×100 are final stages of massive stars evolution <cit.> and are distinguished from other subgroups such as super-massive BHs <cit.>, intermediate-mass BHs <cit.> and hypothetical primordial BHs <cit.>. Hereafter, we focus exclusively on stellar-mass BHs.

<cit.> provided a list of 59 BHs in transient X-ray binaries (XRB). 22 of them have dynamically confirmed mass estimates <cit.>. 5 BHs were detected in High-mass XRBs (defined as heaving a donor mass above 10). Recently, gravitational waves made it possible to discover first double BH merger events <cit.>, and a double NS merger GW170817, which may have formed a low-mass BH <cit.>. Interestingly, <cit.> observed a ∼25 star to disappear after a short brightening, which may be interpreted as a BH-formation event. Up to now, none SBHs were confirmed <cit.>.

Previous studies of BH populations are usually outdated and do not take into account recent progress in our understanding of massive star evolution <cit.>. The earliest estimations predicted ∼100 million BHs in the Milky Way (MW) galaxy <cit.> with as much as 45% of supernovae occurring in close binaries <cit.>. <cit.> provided an upper limit of 1.4×10^9 and <cit.> estimated 1.8×10^8 BHs in the galaxy taking into account  changes in star formation rate (SFR). <cit.> used the  population synthesis code in its earlier version to investigate BHs formed in star formation (SF) bursts. Their results indicate that most of BHs in the MW are actually single objects (not in binaries), whereas these which remain bound will have mostly main-sequence (MS) companions, provided that the stellar population is not very old. 

Recent studies do not provide the estimations for parameter distributions of the entire BH population. For example, <cit.> predicted ∼10^8 BHs in a MW-type galaxy, however, their approach does not take into account binary interactions directly and focuses on binary BHs (BH+BH) only. <cit.> performed a detailed cosmological simulation of the MW evolution (including e.g. changes in metallicity) and predicted, that ∼1 million of BH binaries have already merged in the galaxy, whereas ∼3 millions are still present. However, their study does not account for BHs with lower-mass companions, or disrupted systems, which may actually constitute a bulk of the BH population <cit.>.

Although the current population of BHs is predicted to be mainly single, their progenitors did not have to be single stars. Due to a NK <cit.>, a binary (especially with large orbital separation, or a low-mass companion) may be disrupted <cit.>. Recent observational studies <cit.> suggest that even ≳90% of massive stars (BH progenitors) are born in binaries (or multiple systems). The presence of a companion may significantly affect the evolution of a progenitor and, therefore, final properties of a BH. Some studies used a simplified approach. E.g., <cit.> encapsulated all the interactions in just two parameters. In this study we focused on a detailed consideration of binary interactions, however, we neglected any higher-order systems <cit.>.

A serious problem in obtaining estimations of BH populations is the lack of knowledge about the stellar environment. Taking the MW as an example, not even the total stellar mass is well constrained, let alone the SFH, or the distribution and motions of stars. What is more, a galaxy with complicated evolution and structure cannot be modelled with a uniform population of objects evolved from a single SF burst. Therefore, in our approach we take into account a range of evolutionary models with different initial parameters. Such results can be joined together according to SFH and chemical evolution of a stellar population in order to obtain more realistic distributions.  

Specifically, we propose a step-by-step approach in which we use results from modelling of homogeneous and scalable BH populations in order to build a complicated (non-homogeneous) stellar systems like galaxies. We provide the estimates for a simplified MW galaxy as an exemple (for detailed modeling of the MW galaxy see Olejak et al. in prep.). However, our results may be used to simulate virtually any large stellar population providing the dynamical interactions between stars may be neglected (it is not true e.g. for globular clusters).

The two most important factors that influence BH populations are the metallicity, and NKs. The metallicity of a star significantly affects the mass loss in stellar wind and, therefore, the final compact object's mass <cit.>. <cit.> showed that the maximal BH mass for solar metallicity environment is ∼15, whereas for lower metallicity environments may reach ∼80. However, their study did not consider the pair-instability supernova and pair-instability pulsation supernova <cit.> and was performed for single-stars only. Both, pair-instability supernovae and pair-instability pulsation supernova are believed to produce the second "mass gap" in the distribution of BH masses between ∼50135 <cit.>. However, in binary systems the masses of BHs can fill this range due to mass transfer (MT) phases and mergers <cit.>.

Uncertainties in the modelling of common envelope (CE) phase add additional significant error to the predictions of binary evolution <cit.>. The CE seems to be essential for the formation of XRBs where a strong reduction of initial separation is necessary to produce a Roche lobe overflow (RLOF), or significant accretion from wind. However, it was shown <cit.>, that different prescriptions for the CE phase, give similar results with respect to the population of XRBs. In contrast, the formation of close double compact objects (DCO), with time-to-merger smaller than the Hubble time, which also typically involves the CE phase, is highly influenced by the adopted CE model <cit.>. Nonetheless, as we show in this paper, only a small part of the total BH population resides in XRBs and close DCOs. Most stellar mergers of isolated binaries occur as a result of failed CE ejection, therefore, CE model may strongly influence the population of BHs originating form merger products. However, even though outcomes of mergers may constitute a significant fraction of a BH population, the poorly understood physics of stellar mergers and post-merger evolution make it impossible to reliably describe the population of these objects. In this paper, therefore, we do not include the analysis of CE models.

The fraction of binaries that remain bound after a supernova (SN) explosion as well as peculiar velocities of BHs are particularly sensitive to the assumptions about the NKs. The former impacts the predicted ratio of the number of SBHs to those found in binaries. The latter influences the spacial distribution of BHs in stellar systems (e.g. galaxies). Although the SN mechanism is not well understood, NKs are usually being connected to asymmetries arising in this process. Main candidates are: mass ejection <cit.>, gravitational waves <cit.>, and neutrinos <cit.>. <cit.> derived peculiar velocities of several BH binaries (BHB, i.e. binaries composed of a BH and a non-compact companion) obtaining values between 19÷144. However, it is usually difficult to assess NKs from present-day motions <cit.>. Recent analysis performed by <cit.> has shown that at least some of the BHs in the MW should have obtained high NKs (∼100 km/s), comparable to these of NSs. On the other hand, <cit.> showed that velocities higher than ∼80 km/s are not necessary, although cannot be ruled out. Additionally, <cit.> showed that XRBs with BH accretors have similar spacial distribution as theses with NS accretors, what suggests a similar NK magnitude at birth. See <cit.> and references therein for a recent discussion on the BH NKs.

In this paper, we performed the first simulation of the entire BH population for different evolutionary models. Our analysis is focused on SBHs from binary disruption events (Sec. <ref>), BHBs, DCOs (Sec. <ref>), and potential BHs originating from stellar mergers (Sec. <ref>). We pay a particular attention to the description of the results for the standard model (), which may depict the MW disk field population, and differences in relation to other models.

The main goals of this work are: 1) to describe the general characteristics of the BH population and its dependence on the most important model parameters (e.g., metallicity, NKs), and 2) to provide a reference point for forthcoming in-depth studies focusing on astrophysical problems involving BH populations like: XRBs, microlensing, gravitational wave sources, to name a few. Especially educating are possible comparison studies with the results of actual and future surveys focused on gravitational wave sources (e.g., aLIGO, Einstein Telescope, LISA), X-ray sources (e.g., XMM-Newton, Chandra, NuSTAR), or gravitational microlensing (OGLE, Gaia, LSST), which may result in obtaining better constrains for the evolutionary models. 

Particularly, in Sec. <ref> we estimate a number of microlensing events toward the Galactic bulge. The microlensing method seems to be a promising way to search for SBHs, which, as we show in this work, constitute a vast majority of all BHs in stellar populations.

In order to allow for easier and more flexible usage of the data introduced by this paper, we prepared a publicly available database where all the results are available for download: 
https://universeathome.pl/universe/bhdb.php
 






§ METHODOLOGY


We utilize a recent version of the  population synthesis code <cit.>. The code has been frequently used to study BHs in XRBs and DCOs <cit.>. Recent updates include, but are not limited to, new prescriptions for wind mass loss from massive stars <cit.>, pair-instability supernovae and pair-instability pulsation supernovae <cit.>.

cl
    
    Summary of models

     Model     difference in respect to standard model 
        standard (reference) model: 

         solar metallicity 

         distribution of initial periods

         P(log P)∼(log P)^-0.55

         distribution of initial eccentricities

         P(e)∼ e^-0.42

         BH/NS natal kicks are drawn from 

         Maxwellian distribution with σ=265 km s^-1

         BH natal kicks reduced due to fallback

         moderate slope for high-mass end of the IMF

         Γ=-2.3

        metallicity equal to 10%

        metallicity equal to 1%

        distribution of initial separations log(a)∼1

         distribution of initial eccentricities P(e)∼ e

         BH NKs are inversely proportional to the BH's mass 

        (Eq. <ref>)

        BH/NS natal kick proportional to 

         ratio of ejecta mass and remnant mass (Eq. <ref>)

        flat slope for high-mass end of the IMF (Γ=-1.9)

        steep slope for high-mass end of the IMF (Γ=-2.7)

    List of models. All the main parameters are provided only for  model. For other models only differences in respect to  model are given explicitly.



Here, we analyse 8 main models (Tab. <ref>) which differ in parameters that significantly affect the resulting BH population. Additionally, in the web database, we provide a grid of all 54 models which allow to investigate combinations of the main models. For each model we have simulated the evolution of 2×10^6 binaries from zero-age main-sequence (ZAMS) to disruption, merger, or reaching the age of 15. For mergers that are not DCOs we apply a simple formalism to estimate the endpoint of post-merger single-star evolution (see Sec. <ref>).

Initial primary masses () are drawn from a broken power-law distribution <cit.> with Γ=-1.3 for <0.5 and Γ=-2.2 for 0.5<<1.0. For initial masses >1.0 we chose Γ=-2.3 <cit.>, -2.7 <cit.>, or -1.9 <cit.>. Conclusions of <cit.> were recently revised by <cit.> who argued in favor of a much steeper Γ≈-2.11, or -2.15. Being aware of this, we leave the original value in order to emphasize the influence of IMF steepness on the population of BHs. The mass ratio q=/ (where  is the mass of the secondary – initially less massive star) is drawn from a uniform distribution between 0.08 / and 1. In  model, we assume that the initial distribution of periods is P(log P)∼(log P)^-0.55, and distribution of eccentricities is P(e)∼ e^-0.42 <cit.>. Even though, the results of <cit.> are for stars with initial masses between 15÷60, we extrapolate them to entire investigated range (0.08÷150). Additionally, in model , we test a distribution of initial separations that is flat in logarithm <cit.> with maximal value set as 10^5 and a thermal distribution of eccentricities <cit.>.

The effects of pair-instability supernovae and pair-instability pulsation supernovae <cit.> may significantly alter the final BH mass <cit.>. The instability leads to a significant mass loss for massive stars <cit.>, or disruption of the entire star <cit.>. Following <cit.>, we assume that for stars that form massive helium cores (M_ He=4565) all the envelope above the inner 45 is lost due to pulsations. Furthermore, we assume that stars which form heavier helium cores (M_ He=65135) are subject to the pair-instability supernovae and leave no remnant.

The CE phase is important for the evolution of many binaries and may lead to the formation of a much closer system, or a merger. In a general situation, we utilize the simple energy balance <cit.>. However, the donor type influences significantly the binary survival during this phase <cit.>. MS and Hertzsprung gap (HG) donors lack a clear core-envelope structure <cit.>, so the orbital energy is transferred to the entire star, instead of being transferred to the envelope only, what prevents envelop ejection and leads to a merger. It is not clear when the core-envelop boundary emerges (late HG, or red giant phase), thus two models for the treatment of HG donors in CE were introduced <cit.>. First one (model B), assumes that a CE phases with a HG donor always results in a merger. Second one (model A), producing significantly higher DCO merger rates <cit.> allows for survival in such situations. We adopt model B in the present study and analyse the influence of model A in Sec. <ref>.

We test models with 3 different metallicities. Most of the models have solar metallicity <cit.>. We also introduce two models with lower metallicity:  model (Z=0.002=10%), and   model (Z=0.0002=1%). Additionally, three different NK models are included in our simulations. Those models are described below.



 §.§ Natal kick models


<cit.> performed a study of 233 galactic pulsars proper motions and find out that for young (i.e. with characteristic age τ_ c=P/2Ṗ<3) pulsars the velocity distribution is well fitted with Maxwellian distribution with σ=265 km/s. Velocities of such a population may quite well resemble the NKs of NSs. However, it is not obvious, if the NKs of BHs follow the same distribution. Especially, if NKs are driven by the asymmetries in ejecta, the post-SN fall-back may significantly lower the BH's velocity. Proper motion measurements are available only for a few BHs residing in XRBs <cit.> and the connection between the current peculiar velocity and the NK is not straightforward <cit.>. On the other hand, the analysis of the positions of XRBs harbouring BHs may suggest that distributions of NKs of BHs and NSs are similar <cit.>. 

In the standard model, we lower the remnant's NK proportionally to the fall-back of material after a SN explosion. Precisely, the NK for a compact object is equal


    V_ kick=V_ kick, Maxwell (1-f_ FB),


where V_ kick, Maxwell has a Maxwellian distribution with σ=265 km/s and f_ FB indicates what fraction of ejected mass falls back onto a compact object <cit.>. This model produces both low-velocity and high-velocity BHs which is in agreement with observational estimates for XRBs <cit.>.

If NKs are driven by neutrino based mechanism, it is predicted that the imposed momentum will be the same for BHs and NSs <cit.>. Such NKs will be inversely proportional to the mass of a BH (V_ kick≈ M_ BH^-1). In  model we assume that NSs' NK distribution has a Maxwellian shape with σ=265 km/s, but BHs
NK distribution is lowered. Precisely, 


    [                      V_ kick,NS                               =               V_ kick, Maxwell,;                      V_ kick,BH                               = V_ kick, MaxwellM_max,NS/M_ BH, ]


where M_max,NS=2.5 is an adopted limit for a NS mass above which it collapses to a BH. It must be noted that, even for low-mass BHs (∼5–7), such a prescription will produce significantly lower NKs, than the NKs of NSs, what stands in contradiction with the results of <cit.>. Nevertheless, we included this model as a parameters study.

Finally, we incorporated a simple model proposed by <cit.>, who suggested that the NK distribution can be described as a linear function of a ratio between the mass of the ejecta (M_ ejecta) and the mass of the remnant (M_ remnant). We applied (model ) their fit to the observations provided by <cit.> according to which the NK in 3D is given by


    V_ kick=60 km/ s(M_ ejecta/M_ NS/BH) + 130 km/ s,

where M_ ejecta is mass ejected from the star during SN and M_ NS/BH is mass of a resulting NS, or a BH. Although <cit.> provided their fit only for NSs, we apply it also to BHs as a parameter study.




 §.§ Mergers of binary components


Stellar mergers may constitute a significant fraction of massive stars <cit.>, and therefore, progenitors of BHs. Head-on collisions may play a significant role in dense stellar systems like globular clusters <cit.>, however, stellar mergers in field populations are more probable due to orbital angular momentum loss (e.g. during CE phase). <cit.> predicted that ∼10% of 8<M<20 primaries merge with their companions before a SN and <cit.> showed that ∼8^+9_-4% of massive single stars may actually be products of binary mergers. Objects like Red Novae, or Blue Stragglers are thought to be observed during, or after the merging process <cit.>.
 
In case of non-DCO mergers (including mergers of a compact object with a non-compact object), a merger process is not well understood. Especially, the amount of mass ejected from the system, which is an important factor for post-merger evolution, is poorly constrained. A product of a merger of two MS stars may evolve similarly to a star which was always single, but if evolved stars are involved, the post-merger evolution is more complicated <cit.>. Especially, a merger outcome may evolve unlike any single star <cit.>.

We assume that non-DCO mergers occur in the following situations:

    
  * failed envelope ejection during a CE event <cit.>,
    
  * donor in a CE phase haven't developed a clear boundary between core and the envelope <cit.>, 
    
  * donor's radius exceeds two times its Roche lobe radius 
(R_ donor>2R_ donor, RL) during RLOF.

Due to the fact that the merger physics and post-merger evolution is poorly understood, we provide information on the binary parameters just before the merger to allow for different approaches to this conundrum.

As far as DCO mergers are concerned, many studies were devoted to an in-depth analysis of this phenomenon in stellar populations <cit.>. DCO mergers affect the BH populations, especially the BH mass distribution. As a part of our results, in the public database we provide information about DCO formation and estimate the time to merger using the formula of <cit.>. In order to give predictions and estimate the fraction of BHs which come from DCO mergers, we assumed that the mass loss during a DCO merger is negligible.

In this work, in order to estimate the population of SBH originating from binary mergers, we adopt a simplistic approach of Olejak et al. (in prep.). We include only main channels responsible for >95% of all the mergers. Specifically, for different evolutionary types of stars we use the following procedure:


    
  * MS+MS - Outcome is a MS star. We assume that half of the mass of the lighter component is lost in the process.
    
  * MS+HeS - Outcome is a helium star (HeS) star. We assume that half of the mass of the MS component is lost in the process.
    
  * HeS+HeS - Outcome is a HeS star. We assume that half of the mass of the lighter component is lost in the process.
    
  * NS+MS/HeS - We assume that half of the mass of the MS/HeS star is lost in the process and that it becomes a NS/BH if its final mass is lower/higher than M_ max,NS=2.5.
    
  * BH+MS/HeS - We assume that half of the mass of the MS/HeS star is lost in the process and that it becomes a BH.
    
  * BH+NS/BH - The outcome is a BH with a mass equal to the total mass of the binary before the merger.


We neglect other types of mergers as we found that they constitute only a small fraction (≲5%) of all mergers and predictions for their outcomes are even less certain. For post-merger MS and HeSs we calculate their further evolution assuming that they are on ZAMS, or zero-age helium main-sequence, respectively, in order to find BH predecessors and calculate compact object's mass. We note that presented approach is only quantitative, but helps to estimate the importance of stellar merger for BH populations.




 §.§ General simulation properties


In this work, we concentrate on a description of BH populations originating from massive binaries, which we define as binaries with a primary's (i.e., heavier star on ZAMS) initial mass >10. Our analysis does not include BHs in triples and higher-order systems <cit.>, as well as BHs in binaries formed due to stellar encounters in dense stellar systems <cit.>.

The total simulated stellar mass is 4.8×10^8 for the  model and 1.1×10^9 (2.4×10^8) for the () model. Throughout this paper, we assume 50% binary fraction on ZAMS for low-mass stars (M_ ZAMS<10) and 100% binary fraction for heavier stars (M_ ZAMS≥10). The simulated stellar mass corresponds to about ∼0.8% (1.9% for  model, or 0.4% for  model) of the stellar mass of the MW equivalent galaxy <cit.>.

For the purpose of presentation, all results are scaled to the stellar mass of a MWEG <cit.> and the SFH is chosen to be constant through the last 10 (SFR=6 M_⊙ yr^-1). Such a model, although simple, allows us to draw general conclusions. Raw results may be scaled by using any realistic total stellar mass and any SFH may be applied. In general, the changes in SFH are not as important as the total mass of formed stars, as there is only a small delay time between the formation of a binary (ZAMS) and formation of a BH (typically less than a  few ×10), which are, typically, single. We note, that for a part of the BH population the SFH may, actually, be important, e.g. for DCOs, which have steep time-to-merger distribution <cit.>. Additionally, if the spacial distribution of BHs is concerned, the change of SFH may play a role, because it affects the time which BHs have to spread throughout the space. What is more, BHs formed earlier have more time to be ejected from the stellar systems due to their NKs, or as a result of dynamical interactions. Deeper analysis of changes in spacial distributions, SFH, and chemical evolution are left for separate studies (e.g. Olejak et al. in prep.).









§ RESULTS
 

lcrlrlrlrl
    
    Number of BHs in the Milky Way Equivalent Galaxy
     Model     N_ BH,tot     2cN_ dSBH     2cN_ BHB     2cN_ BH,DCO     2cN_ mSBH
        1.18     5.47     (49.7%)     1.56     (1.4%)     5.16     (4.7%)     4.87     (44.2%) 

        1.18     3.47     (30.6%)     2.06     (1.8%)     1.27     (10.8%)     6.37     (56.8%) 

        1.28     3.77     (30.0%)     1.56     (1.2%)     1.97     (15.4%)     6.67     (53.4%) 

        1.18     6.67     (58.8%)     3.46     (3.0%)     6.96     (6.1%)     3.67     (32.1%) 

        1.18     6.47     (58.1%)     5.94     (0.1%)     1.14     (0.0%)     4.67     (41.8%) 

        1.18     6.57     (58.5%)     5.44     (0.0%)     2.73     (0.0%)     4.67     (41.4%) 

        2.78     1.48     (52.6%)     3.26     (1.2%)     1.37     (4.9%)     1.18     (41.3%) 

        3.47     1.67     (46.8%)     5.65     (1.6%)     1.66     (4.7%)     1.67     (46.8%) 

    Number of BHs for different tested models (see Tab. <ref>). Values are presented for a simple MW model with a total simulated stellar mass of 6×10^10 and constant star formation during the last 10. Column headers stand for: N_ BH,tot - total estimated number of BHs; N_ dSBH - number of single BHs from disrupted binaries; N_ BHB number of binaries harbouring a BH and a non-compact companion; N_ BH,DCO - number of BHs in DCOs (BH+BH is counted as two BHs; see Tab. <ref>); M_ mSBH - a rough estimate of the number of stellar mergers (including DCO) which will be massive enough to form a BH (see Sec. <ref>). 
    


The main results of the simulations are summarized in Tab. <ref>. Most of the BHs are predicted to exist as SBHs either as a result of a binary disruption (dSBH), or a stellar merger (mSBH; with a BH involved, or producing a star massive enough to form a BH). The population of BHs in binaries (both in DCOs and in BHBs) is about an order of magnitude smaller and in the case of models with increased BH NKs ( and ) even two orders of magnitude smaller. Only a very small fraction of the BHBs are the interacting ones (see Sec. <ref>).



The initial separation influences significantly the fate of a massive binary. In general picture, if the separation is very large (a≳ 3,000) the binary will mostly disrupt and produce SBHs (Sec. <ref>). A low separation (≲ 30) will lead frequently to a non-DCO merger (Sec. <ref>), which may be massive enough to form a SBH. Objects with medium initial separations have much more uncertain fate and may merge, become disrupted, or form a bound binary with at least one BH inside (Sec. <ref>). In the case of a binary harbouring a BH, the system after some time may still merge (e.g. as a double compact object merger, Sec. <ref>), or be disrupted during the formation of a second compact object (e.g.,  formation route, Sec. <ref>). 

ll
    
    Typical formation routes of BHs
     Route     evolutionary route 
    2cSingle BHs from disrupted binaries

        MT1(1/2/4/5-1/2/4/5) SN1 Disruption SN2 

        SN1 SN2 Disruption


    2cBH XRBs

        CE1(4/5-1;7/8-1) SN1 MT2(14-1/2/3/4/5/6) 

        CE1(4/5-1;7/8-1) SN1 MT2(13-1/2/3) AICBH1

         MT2(14-1/2/3)


    2cDCOs

        SN1 SN2

        MT1(1/2/4/5-1/2/4) SN1

         CE2(13/14-4/5;13/14-7/8) SN2


    2cMergers

        CE1(2/3/4-1;7-1) MT1(7-1) Merger

    Schematic representations of the evolutionary routes for typical BH formation channels. Wide BHBs (Sec. <ref>) are not included in the table as their evolution does not include any interactions and their evolutionary route is 'SN1' by definition. Details of the post-merger evolution are not included in our study. Symbols in evolutionary routes represent: SN1/2 - supernova of the primary/secondary; MT1/2 - mass transfer (primary/secondary is a donor); CE1/2 - common envelope (primary/secondary is a donor; numbers in parenthesis represent typical primary evolutionary type (left) and secondary (right); first two numbers represent initial types (prior to CE), whereas, last two numbers represent the final types (after CE)); AICBH1 - accretion induced collapse of a NS primary to a BH (assumed to occur after a NS obtain the mass of M≥ M_max,NS=2.5). Evolutionary types (numbers inside parenthesis) represent: 1-main sequence, 2-Hertzsprung gap, 3-red giant, 4-core helium burning, 5-early asymptotic giant branch, 6-thermal pulsing asymptotic giant branch, 7-helium main sequence, 8-evolved helium star, 13-neutron star, 14-black hole.
    









 §.§ Single BHs from disrupted binaries


If a system hosting a BH progenitor becomes disrupted, what may occur before, after, or during the BH formation, a compact object becomes an dSBH (we use name dSBH to represent a SBH from a disrupted binary, as SBH can also form due to mergers; see Sec. <ref>). Binaries may become disrupted due to a NK which one component obtains after a SN explosion, due to a Blaauw kick <cit.>, or interaction with a third star. Although in our simulations we included only isolated binaries, we note that interactions, even in sparse stellar systems like the Galactic disk, may significantly alter the binary evolution <cit.>.


lccccl
    
    Single BHs from disrupted binaries
     Model     Route     Number     a_ ZAMS[]  
            3.6×10^7 (67%)     28÷400 

             1.7×10^7 (32%)     >4400 

            2.0×10^7 (58%)     20÷800 

             1.2×10^7 (33%)     >7700 

            1.5×10^7 (40%)     17÷180 

             2.0×10^7 (55%)     >1600 

            4.0×10^7 (61%)     28÷1000 

             2.4×10^7 (36%)     >6200 

            4.0×10^7 (62%)     28÷400 

             2.2×10^7 (35%)     >5200 

            4.1×10^7 (63%)     28÷400 

             2.2×10^7 (34%)     >5200 

    
    
            9.4×10^7 (66%)     28÷430 

             4.7×10^7 (33%)     >5400 

            1.1×10^7 (68%)     28÷400 

             5.1×10^6 (31%)     >5000 

    Number of single BHs originated from disrupted binaries. Only main evolutionary routes are included). The initial separation (a_ZAMS) is the main differing factor between the main routes. See Sec. <ref> for details.



We found out that the tested initial parameters distributions have little effect on the properties of the resulting population of dSBHs. The steepness of the IMF changes number of resulting dSBH by Δ N_ dSBH≈70%/160% due to the decreased/increased relative number of BH progenitors for / models. The effect of metallicity and NK models is not as strong (Δ N_ dSBH≈40% and 20%, respectively).




  §.§.§ Typical evolutionary routes


We found that two evolutionary routes (see Tab. <ref>) dominate most of the dSBH production and are common for all models. Number of dSBHs produced through these routes and for different models is presented in Tab. <ref>. Both channels involve a binary initially composed of two massive stars that will undergo SN explosions, or collapse directly into BHs. The former case may lead to a binary disruption (due to NK and/or mass loss). The main difference between the routes comes from the initial separation. In the case of , the initial distance between stars is lower than ∼1000, whereas for  it is usually larger then ∼5000. Consequently, the former route leads through an interaction phase (typically the MT), whereas a binary evolving through the latter one experiences no interactions.



For a typical binary evolving through  route (Fig. <ref>, upper plot), the primary is about 26 on ZAMS, whereas secondary is slightly lighter (∼21). The separation is modest (∼300). The primary evolves faster and after ∼6.7, while expanding as a HG star, fills its Roche lobe and a MT begins onto the secondary, which is still on its MS. Although the MT is relatively short (∼4), the primary loses its hydrogen envelope and becomes a ∼7 helium star. Half of the expelled envelope is accreted by the secondary which grows to ∼29 still being on its MS. After about 1 primary goes through SNIb/c and becomes a ∼1.5 NS. A NK disrupts the system. Due to its significant mass, the velocity which the secondary obtains is very low (∼5). The massive secondary needs only ∼2 to become a BH with a mass of ∼8. The NK in a case of such a BH is negligible, therefore, its velocity doesn't change significantly. It is noteworthy that most of the BHs formed through  route originate from secondary stars, i.e. less massive ones on ZAMS.

The  channel is similar to the one proposed for the formation of single massive stars (potential BH progenitors) by <cit.>. Using different population synthesis code, they found that 86^+11_-9% (depending on the model parameters) of binaries evolving through this channel will become disrupted during the first SN.

In the other typical route (; Fig. <ref>, bottom plot) the primary is about 25 on ZAMS, whereas the secondary is ∼10. The initial separation is as large as ∼6000. As a consequence, their Roche lobes are huge and no interaction is possible throughout their evolution. The heavier star evolves faster and in ∼8 forms a ∼7.5 BH with a small NK. The separation grows to ∼10,000 what results from a loss of mass in stellar wind from both stars. After additional ∼17 secondary explodes as a SN and forms a NS. A strong NK disrupts the system. Due to a large mass ratio between compact objects (q≈7) the post-disruption velocity of the BH is relatively small (∼5).




  §.§.§ Mass distribution of dSBHs




The dSBHs' mass distribution is presented in Fig. <ref>. The shape is similar to the distribution of BH masses in single star evolution <cit.>. The main peak (∼7–8 for solar metallicity ( model; Z=), or ∼15–30 for lower metallicities (models  and ; Z=10% and Z=1%, respectively)) relates to BHs formed from ∼20–35 progenitors through failed SN explosion <cit.>. The second peak (∼15 for , or ∼40 for  model) comes from BHs originating from most massive stars (≳100) on ZAMS, which lose a large part of their mass in stellar wind. In the case of low metallicity environments (model ) the peak results from the pair-pulsation SN <cit.>, which prevents the formation of heavier BHs in single star evolution. The third peak (∼22 for , or the tail extending to ∼45 for  and ) is formed through binary interactions and is not present in distributions for single star evolution <cit.>.

The similarity of the mass distributions of BHs for single stars and binaries is a direct consequence of the fact that massive binaries, which later on became disrupted, usually interact only little (mass of BH progenitor is not changed significantly; ), or not at all (), so in most cases the binary components evolve as in isolation. 

A notable difference is the presence of the high-mass peak and tail in the mass distributions. It is most prominent for the models with solar metallicity where the tail extends up to ∼22, whereas the heaviest BHs forming at that metallicity from single stars are only ∼15 <cit.>. These massive SBHs constitute only a tiny fraction of the population (<1%). Typically, the initial total mass of a system in which such over-massive BHs are formed is above 200 and a mass ratio is ∼1. The primary, being slightly heavier, is the first to expand and fill its Roche lobe while being on a HG. At the moment of interaction both stars have masses of about 60 due to strong mass losses in stellar wind. When the mass transfer ceases, the primary becomes a ∼25 core helium burning (CHeB) star. The secondary being now much heavier (∼80), appears rejuvenated <cit.>, evolves much faster and after ∼400 forms a BH first with a mass of ∼21. The primary needs additional 100 to become a 7.3 BH. During the second SN, NK easily disrupts a wide (∼15,000) binary.




  §.§.§ Significance of metallicity




The relation between the number of dSBH and metallcity was found to be non-monotonic. The largest number of dSBHs in respect to metallicity is produced in the  model (5.4×10^7), whereas the smallest is in the  model (3.4×10^7). For  model the number of dSBHs is slightly higher (3.7×10^7) than in the  model.  

The number of dSBHs produced through the  route increases monotonically with metallicity. The higher is the metallicity, the stronger is the line-driven stellar wind. Only a fraction of this mass may be accreted by a companion, thus bulk of it leaves the system what results in orbital expansion. Consequently, the stronger is the stellar wind, the wider the binary becomes. Simultaneously, the stronger mass loss also makes the companion lighter during the NS formation. Wider orbits and lighter companions make the system easier to disrupt, what as a result, produces more SBHs. 

A different mechanism results in the non-monotonic relation between the metallicity and the number of BHs produced through the  route. Metallicity affects the evolutionary expansion of the stars which happens to be the strongest in  environment and the smallest in the  environment (Fig. <ref>). Only binaries with separations high enough to avoid Roche lobe filling may evolve without interactions. The number of such systems will be smaller if the nuclear expansion is on average higher as in the case of  model.
 




  §.§.§ Initial BH velocities and high velocity dSBHs


When a SBH forms (after disruption of the binary, or BH formation, whatever happens later) its velocity is not only a result of the motion in the gravitational potential (e.g. of a galaxy), but also preceding evolutionary processes, like NK, Blaauw kick, and binary disruption. The velocity which is the result of the latter processes calculated at the moment of the binary disruption, or the formation of the BH (whatever happens later) we call here an initial BH's velocity (). Such a definition in which we do not involve the motion in the gravitation potential allows for application of our results to any gravitational potential (e.g. different galaxies, or different models of the MW galaxy). In this work, we provide the magnitudes of 3D velocities (i.e. lengths of the velocity vectors) and assume that the distribution of NKs is isotropic. We note that in a realistic situation, the velocity of a dSBH will be modified by the gravitational potential of a stellar system and interactions with other stars. Although these effects may also affect the binary evolution, the short dSBH formation timescale (≲50) allow us to assume that the effect is usually negligible.

The highest values of  are obtained for BHs with the highest NKs. We found that majority of dSBHs with the highest initial velocities (>300 km/s) within the highest 10% were formed through  route and involve lowest BH masses and the shortest pre-disruption periods. We note, that <cit.> found that the majority of the secondaries (∼90%) after the disruption will have low velocities (≲20 km/s) what agrees with our predictions for standard NK model.

 distributions for all tested NK models are presented in Fig. <ref>. The maximal  of a dSBH may reach ∼500 for  model and nearly ∼1000 for  model. The higher velocities for lower metallicities are the result of a formation of closer binaries through  route, what together with higher masses of BHs, give larger orbital speeds <cit.>.



A distribution of  for standard NK prescription <cit.> has two distinctive parts (Fig. <ref>). The low-velocity peak is connected with BHs formed with no NK due to a heavy fallback, which, as we assume, decreases the NK, or in a direct collapse associated with no NK. The wide high-velocity component comprises BHs which are much lighter, thus formed with no significant fallback. Models  and  lead to a comparable number of dSBHs (6.4×10^7 and 6.5×10^7) as the  model (5.4×10^7), however, the shape of  distribution is significantly different.  model lacks a low-velocity peak as the NK is independent of the fallback. Therefore, the average NK is much higher than in the  model and typically equals ∼125 km/s for a typical BH mass of 7.5. For  model there is typically full fallback associated with the formation of BHs with masses of ∼7–8 and no NK is applied. Even more dissimilar is the  distribution for the  model, where NKs are proportional to the ratio of ejected mass to compact object mass (v_ NK∼ M_ ej/M_ BH; eq. <ref>) and typically, for 7.5 BHs, are ∼130 km/s.







 §.§ Bound systems


The number of BHs which are bound in binaries is between 5.74 – 2.17 (depending on the model, see Tab. <ref>) per MWEG. The fraction of BHs which reside in binaries for most of the tested models is ≲17%. For models  and , the fraction of BHs in binaries is ≲0.1%, therefore, the model of NKs highly influences the survival of BH binaries. 

The maximum mass for a BH is ∼22 (or ∼45 for lower metallicity models). The binary with such a massive BH exists for ∼100 before it disrupts after the second SN. Afterwards, these BHs become the most heavy dSBHs (Sec. <ref>).

Most of the binaries harbouring BHs are wide, i.e. formed without any interactions between stars. BHs reside  predominantly in DCOs (Sec. <ref>). The rest is mostly accompanied by WDs (such binaries we treat as non-DCOs, i.e. BHBs; Sec. <ref>). A small fraction  of BHBs could be observed during their MT phase, although this fraction is strongly dependent on the adopted SFH (Sec. <ref>).




  §.§.§ Double compact objects
 

l|ccrcrcrcrcrcr
    
    Double compact objects containing black holes
     Model     N_BH,DCO     2cBH+BH     2c     2c     2ca[]<10^3     2ct_merge<10
        5.16     4.96     (95.9%)     2.46     (90.1%)     1.64     (0.6%)     2.54     (0.9%)      1.74     (0.6%) 

        1.27     1.27     (96.6%)     2.66     (42.6%)     4.05     (3.4%)     1.76     (27.3%)      8.15     (13.0%) 

        1.97     1.97     (97.5%)     5.76     (57.8%)     1.86     (18.1%)     2.76     (27.0%)      1.96     (19.0%) 

        6.96     6.66     (95.6%)     3.46     (93.1%)     2.44     (0.7%)     3.04     (0.8%)      2.34     (0.6%) 

        1.14     5.13     (45.5%)     8.82     (10.3%)     1.63     (19.1%)     6.83     (79.4%)      3.23     (36.8%) 

        2.73     2.52     (9.5%)     2c–     2c–     2.53     (100%)      1.32     (5.0%) 

        1.37     1.37     (96.6%)     6.26     (90.6%)     3.84     (0.6%)     6.04     (0.9%)      3.84     (0.6%) 

        1.66     1.56     (95.3%)     7.75     (90.7%)     3.73     (0.4%)     5.73     (0.7%)      3.73     (0.4%)
    Results for double compact objects harbouring BHs (i.e. BH+BH and BH+NS systems). N_BH,DCO - number of BHs in double compact objects; BH+BH - number of BHs in double BH systems (the rest reside in BH+NS systems);  and  - number of BHs formed through specific evolutionary routes; a[]<10^3 - number of BHs residing in close DCOs, i.e. systems with separations lower than 10^3; t_merge<10 - number of BHs in systems with time to merger after second supernova smaller than 10.
    


In this study, by DCOs we understand only these harbouring at least one BH (i.e. BH+BH and BH+NS), totally neglecting NS+NS systems. BH+WD systems are included in non-DCO binaries (BHBs; see Sec. <ref>). The fact that the majority of binaries harbouring BHs are actually DCOs is a direct consequence of the adopted uniform distribution of mass ratios on ZAMS - if the primary is a BH progenitor (≳20) the most likely companions are NS and BH progenitors (≳8). NSs generally receive higher NKs than BHs and many potential BH+NS binaries get disrupted during the NS formation <cit.>, therefore, BH+BH systems typically dominate (∼95–98%; see Tab. <ref>) over BH+NS binaries in most tested models. The exception are models with higher average NKs ( and ) where most of DCOs are compact BH+NS systems. Many groups studied the formation of DCOs harbouring BHs over the years <cit.>. However, the main focus of those studies was lied on the merger rates and properties of gravitational wave sources and, hence, they performed no deeper analysis of wider DCOs. Here, we present the analysis of the entire DCO population as a part of the total BH population. Therefore, our results may be applied not only to the study of merging DCOs, but also to other phenomena, e.g. microlensing by DCOs. 

In Tab. <ref> we present detailed results for the main models. The number of BHs in DCOs is typically between 1.56 – 1.97 (note that BH+BH counts as two BHs), although for  and  models, the number drops to 5.13 and 2.52, respectively. The low number of DCOs in the models with higher average NKs is connected to the fact that most of binaries with massive components which are potential progenitors of DCOs are wide (route ; see below), thus easily disrupted by even moderate NKs.



There are two main evolutionary routes leading to the formation of DCOs harbouring BHs (Tab. <ref>). The main feature that distinguish these routes is the initial separation. It is connected with the lack (), or presence () of the pre-DCO formation interactions (MT and CE), which results in the formation of wide (a≳10^4), or close (a≲10^3) DCOs, respectively (see Fig. <ref>). Among close DCOs (mainly route ) are potential DCO merger progenitors <cit.> and their number and relative fraction strongly depends on metallicity and NKs as previously noted by <cit.>.



Route  is the main formation route of DCOs if standard NK model is concerned. In this formation channel none of strong interactions are present (Fig. <ref>). In a typical situations the stars on ZAMS are not very massive (∼25 and ∼22) and the separation is large (∼12,300). After about 7.7 the primary explodes and forms a ∼7.5 BH. In additional ∼1, the secondary forms a ∼7.5 BH in a SN explosion. The final separation is ∼37,500.

The other formation route () is responsible for much smaller fraction of DCOs (≲21%), but systems which form through this route have much shorter separations and, thus, may merger after less than 10 due to the emission of gravitational waves. In a typical case, the binary is massive on ZAMS (≈57 and ≈38) and the separation is small (≈600). The primary evolves fast and after ∼4.1 fills its Roche lobe while expanding on the HG. The MT commences and results in the loss of nearly entire hydrogen envelope by the primary, which shortly after becomes a ∼16 HeS. Half of the primaries envelope is accreted by the secondary which grows to ∼49 while still on its MS. After additional ∼600, the primary forms a ∼8 BH in a SN explosion. The secondary starts to expand after leaving the MS and while being a CHeB star fills its Roche lobe and a CE commences. As a result the separation shrinks from ∼4000 to about 5.2 and the secondary loses its entire envelope and becomes a ∼15 HeS. A second SN occurs ∼600 later and ∼7.7 BH forms. The NK is not very strong, thus the separation remains small (∼7.3). Such a binary is estimated to merge after ∼420. This channel was previously thoroughly analysed in the context of DCO mergers by e.g. <cit.>.

BH+NS formation through  and  routes is qualitatively similar to the formation of BH+BH systems, although the ZAMS secondary masses are significantly lower (∼3× for  and ∼1.5× for ). Additionally, the progenitors are more easily disrupted during the second SN and have a smaller chance of CE survival due to more extreme mass ratio. This results in a significantly lower fraction of BH+NS systems among DCOs than BH+BH systems, in spite of a higher abundance of potential progenitors in initial (ZAMS) populations[The ratio of potential BH+NS progenitors to potential BH+BH progenitors (i.e. binaries with components' masses high enough to form a NS, or BH in single star evolution) may be calculated as

    f=∫_M_ NS,ZAMS,max^150P_ BH+NS(M_1)  IMF(M_1) dM_1/∫_M_ NS,ZAMS,max^150P_ BH+BH(M_1)  IMF(M_1) dM_1,

where P_ BH+NS/BH+BH(M_1) is the probability that a companion of a BH progenitor with mass M_1 is a NS/BH progenitor. Assuming flat mass ratio distribution for companions (P(q)=1/(M_1-M_ NS,ZAMS,min)),

    P_ BH+NS(M_1)=M_ NS,ZAMS,max-M_ NS,ZAMS,min/M_1-M_ NS,ZAMS,min
 and 
    P_ BH+BH(M_1)=M_1-M_ NS,ZAMS,max/M_1-M_ NS,ZAMS,min,

where M_ NS,ZAMS,min/max=8/22 are the minimal/maximal ZAMS masses for a NS progenitor (assuming no interactions with the companion) for solar metallicity (). IMF(M_1)∝ M_1^Γ is the initial mass function. 150 is the upper limit for the ZAMS mass adopted in our calculations. Assuming Γ=-2.3 and solar metallicity ( model) f≈1.1 meaning more BH+NS than BH+BH progenitors in initial populations. This value is slightly affected by metallicity, as Z=1% ( model; where M_NS,ZAMS,max≈19) gives f≈1, and significantly affect by the steepness of the IMF, as Γ=-1.9 ( model) gives f≈0.9, i.e. more BH+BH progenitors.]. There are typically (except  and  models) ≳20 times more BHs in BH+BH binaries than in BH+NS. We observe no significant difference between the shapes of BH+BH and BH+NS distributions of orbital separations at the moment of of DCO formation (Fig. <ref>).

We have found that the number of BHs in DCOs is 3–4 times higher in lower metallicity environments than in the  model. The number of wide DCOs () is only slightly affected by metallicity (by a factor of ≲2.5), whereas the number of close DCOs () is significantly affected (by more than two orders of magnitude). The latter is an effect of higher BH masses (typically ∼15) in low metallicity environments in comparison to  environments (typically ∼7–8). Heavier BHs in the standard NK model usually obtain lower NKs due to fallback, or no NK when formed in a direct collapse. Therefore, DCOs in lower metallicity environments, where BHs form with higher masses, are less frequently disrupted <cit.>.

Significantly different situation concerns the relation for close systems (), which in lower metallicity environments are about two orders of magnitude more numerous than in the  model. This results from an interplay between mass loss in stellar wind and stellar expansion (Fig. <ref>). The former is proportional to metallicity and leads to orbital widening due to mass loss from the system. On the other hand, solar metallicity stars do not expand as significantly as these in  model. This means that massive binaries in high-metallicity environment (e.g.  model) more frequently have separations too large to go through a CE phase which is an essential step in the  route. Consequently, lower-metalicity models have more close DCOs, than high metallicity ones. A similar result was recently obtained by <cit.> with the use of SEVN population synthesis code.



Distributions of component masses are shown in Fig. <ref>. For  formation channel (wide DCOs) distributions for primaries and secondaries are joined together, except BH+NS systems for which only BH masses are shown. The results resemble these for BHs from disrupted binaries (Fig. <ref>). For  channel (close DCOs) the total mass of a binary is provided. Such information is more important from the point of view of microlensing surveys, because a close binary will most probably act as a single lens. For lower metallicites the total mass may reach ∼70–80 what matches the current range of observed mass in BH+BH mergers <cit.>.

Within our models up to ∼15% of BHs are found in DCOs (Tab. <ref>) which are predominantly wide <cit.>. For the standard NK model, up to 19% of them have merger times (t_ merge) smaller than 10, what mainly depends on metallicity <cit.>. In the  and  models, although the fraction of merging systems may be much higher (∼37% for  model), the number of these systems is significantly smaller. It is a consequence of higher average NKs (∼125–130 km/s in contrast to ∼10 km/s in  model), what frequently disrupts wide binaries. 



Typical center-of-mass initial velocities (v_ COM,0; i.e. binary velocities just after the second SN not including velocity resulting from the motion in a gravitational potential) are below ∼20 km/s (Fig. <ref>). Only less than 15% of DCOs (mainly compact; route ) possess velocities exceeding 20 km/s. Therefore, vast majority of the DCOs remains in the vicinity of their birth places. Highest velocities (v_ COM,0≳200 km/s) are obtain by BH+NS systems, whereas systems with lower velocities are dominated by BH+BH binaries, what results from typically higher NKs obtained by NSs, than BHs. Typical velocities of merging systems (t_ merge<10) are below ∼50 km/s for  model, so most of the mergers are predicted to occur in the vicinity of the formation places <cit.>. However, some systems may obtain also high velocities (≳200 km/s). These are mainly BH+NS systems, which we predict to merge away from their birth environments. For  and  models, DCOs obtain typically intermediate velocities (∼50–130 km/s; mainly route ; typically, t_ merge≲10), although some systems (mainly BH+NS) may have velocities as high as ∼400 km/s.




  §.§.§ Wide non-DCO binaries


l|ccc|ccc
    
    Wide and mass-transferring binaries
     Model     3cWide BHB     3cMTBHB 
     N     a[]     Typ. comp.      N     a[]     Typ. comp.
        1.56      ≳10,000     CO WD (86%)     8.51     ≲115      MS (75%) 

        7.55      ≳8,300      CO WD (79%)     4.23     ≲110     MS (85%) 

        1.26      ≳5,800      CO WD (85%)     4.42     ≲84      MS (91%) 

        3.36      ≳8,300      CO WD (74%)     5.72     ≲108      MS (73%) 

        3.13      ≳4,800      CO WD (89%)     1.72     ≲90      MS (86%) 

        1.0-4     ≳14,000     CHeB (100%)      1.82     ≲75      MS (82%) 

        3.16      ≳10,000     CO WD (86%)     1.32     ≲113      MS (73%) 

        5.55      ≳10,000     CO WD (86%)     3.21     ≲116      MS (70%) 

    Results for wide and mass transferring binaries with one of the components being a BH. N - number of binaries (also number of BHs); a - typical separation; Typ. comp. - typical evolutionary type of the companion: MS - main sequence; CHeB - core helium burning; CO WD - carbon-oxygen white dwarf.
    


Throughout the paper we use the term wide BHB referring to systems harbouring a BH with a non-compact companion and experiencing no strong interaction during their entire evolution (the separation between the stars is large and the stars never fill their Roche lobes). In this section, we focus on wide BHBs, excluding wide DCOs, which were discussed in previous section. 

Tab. <ref> summarizes our results for the main models. The typical number of wide BHBs is between 5.55 and 3.36. Only for models with higher average NKs ( and ) the predicted number of such binaries is much lower (3.13 and 1.0-4, respectively). The probability for a binary to obtain a low NK (≲20 km/s) in  and  models is much lower than in  model (c.f. Fig. <ref>) and typically NKs are ≳100 km/s. Therefore, in these models wide binaries are more easily disrupted during the formation of a BH. In contrast, standard NK model <cit.> predicts mostly low, or negligible NKs for BHs which allows for more frequent survival of wide binaries.

The number of wide BHBs is higher in  model than in  and  models mainly because of the more significant orbital expansion resulting from mass loss in stellar wind, which is stronger for higher metallicities. Additionally, stellar expansion is higher in  model, thus futher increasing the chance of interactions (Fig. <ref>). In  model the stellar expansion is smaller than in the  model, therefore, results in higher number of wide BHBs.

BH companions in wide systems are predominantly CO WDs (74–89%). WDs, being the final evolutionary stage of low-mass (≲8) stars' evolution, are simultaneously the longest evolutionary stage of the evolution of stars with initial masses 2≲ M≲8 during the Hubble time. For comparison, low mass stars (M≲1) spend most of this time on the MS, whereas, heavier stars (M≳8) quickly (≲50) end their evolution and form a NS, or a BH. The typical separations of wide BHBs are a≳5000.



In an exemplary evolution of a typical wide BHB, the progenitor system on ZAMS is composed of a ∼25 primary and a 6.0 companion (see Fig. <ref>). 
The initial separation is large (∼16,000), which precludes any interactions. 
The separation further increases due to wind mass-loss from the stars and just before the SN reaches ∼34,000. The primary, being heavier, evolves faster and in ∼7.8 becomes a BH. The NK and the Blaauw kick are crucial at that point. Only a small total kick will allow for a survival of the binary. The secondary is still on its MS for another ∼60.  At the age of about 76 a CO WD forms. The final separation may reach 56,000.



Mass distribution of BHs in wide BHBs (Fig. <ref>) resembles this for dSBHs as most dSBHs are not interacting significantly prior to disruption, thus frequently originate from wide BHBs. Two distributions differ in the presence (dSBHs), or the lack (wide BHBs) of the high mass peak (∼22) for  model. Also the high-mass tail (above ≳15 for  model and ≳40 for  and  models), which originates from wind accretion by a BH, or its progenitor, in the case of non-interacting binaries is less pronounced, than in distributions for interacting binaries. The maximal BH masses are ∼16 for  model and ∼46 in lower metallicities. The distribution of companion masses is dominated by WDs (carbon-oxygen and oxygen-neon), which form a characteristic peak between ∼1 and the Chandrasekhar mass (M_ Ch=1.44). The long tail on this distribution is mainly composed of MS stars, which may have masses up to ∼80 for solar metallicity, or ≳100 for lower metallicity models. The massive companions (M≳8) are going to explode in SN explosion which may disrupt the BHBs (either through a NK, or a Blaauw kick) populating dSBHs, or form a wide DCOs (discussed earlier in Sec. <ref>).

Although in rare cases the initial center-of-mass velocity (v_ COM,0) of wide BHBs may reach ∼80 km/s, typically  (≳95% of cases) they are negligible (≲10 km/s). Due to the fact that orbital velocities are typically also small (≲ 20 km/s), the binaries rarely leave their birth places if evolved in isolation. We note, that in realistic situation a wide BHB is frequently influenced by fly-by encounters <cit.>.




  §.§.§ Mass Transferring binaries with BH accretors


Mass transferring BHBs (MTBHB) have different evolutionary routes and properties than the typical BHBs, which are dominated by non-interacting systems (wide BHBs; Sec. <ref>). Although we note that BHs were detected also in high-mass XRBs where MT rates high enough to fuel an accretion disk may occur also through companion's stellar wind <cit.>, in this study only RLOF systems were included.

For any adopted SFH the fraction of BHs accreting mass from their companions through RLOF is small. However, this small fraction is more a result of the brevity of MT phase in comparison to the evolutionary time-scales, than the rarity of BHBs that are close enough to commence RLOF. Here, for the sake of presentation, we assumed a constant star formation for the last 10. In such a case, only <1% (see Tab. <ref>) of all BHBs contain a donor which is filling its Roche lobe and transferring mass onto a BH. MTBHBs can be perceptible as an X-ray binaries <cit.>. We found out that in our models 2 – 28% of BHs have gone through a MT phase which lasts typically ≲40. A BH accretes on average 0.7 – 3.1 in that time. 

Our simple model predicts 32 – 4,200 MTBHBs per MWEG. The number is the highest in lower metallicities (4.42 and 4.23 systems for  and  models, respectively) where it is easier to start RLOF MT <cit.>. Also in  model we see ∼8 times higher predicted number of MTBHBs, than in the  model. It is a result of assumed in this model thermal distribution of eccentricities (P(e)∼ e), what significantly increases average initial eccentricity, thus lowering the initial periastron distance, what makes RLOF easier.

There are two main evolutionary routes (Tab. <ref>) typical for all models (comprising 79÷99% of all MTBHBs). Both channels differ mainly in the mass of an accreting BH. In the  route, typical BH masses are above ∼5, whereas in  route, they are around ∼3. This distinctness stems from the different BH formation mechanisms. In , a BH forms directly after a SN, whereas in , a NS forms first and after a period of mass accretion reaches a critical mass (M_ max,NS) and collapses to a BH. The critical mass is in general not well constrained and depends on the applied equation of state and rotation <cit.>. In our simulations a value of M_ max,NS=2.5 was adopted. 



An exemplary system which becomes a MTBHB through  route begins its evolution having ZAMS masses of 42 and 5.4 (Fig. <ref>, upper plot). The separation is moderate, a_ ZAMS≈5000. In ∼4.9, the primary evolves off the MS and commences a CE phase. Afterwards, the separation is reduced to ∼13. After additional 0.5, the primary becomes a BH receiving a small NK. The companion needs additional ∼60 to expand due to nuclear evolution and fill its Roche lobe. The MT prolongs for 16 during which the companion loses ∼90% of its mass and evolves off the MS. The BH accretes about 1.5.

As far as the  route is concerned (Fig. <ref>, lower plot), the initial masses on ZAMS are much smaller and in a typical case equal ∼23 and ∼3.3. The initial separation is large (a≈9400). Similarly to  the primary evolves off the MS and commences a CE phase while being a ∼10 CHeB star. The outcome is a compact binary (a≈40) composed of a ∼7.9 HeS and ∼3.5 MS star. The primary is not heavy enough to form a BH. Instead, in 100 forms a heavy (∼1.8) NS with a strong NK. The orbit becomes highly elongated (a≈200, e≈0.96). The secondary, which is now the more massive star, expands due to nuclear evolution and after ∼200 fills its Roche lobe during periastron passage and a MT begins. In such a situation, we assume that the orbit is immediately circularized. In ∼500, the NS accretes >0.7, the critical mass is reached, and it collapses to a low-mass (2.5) BH. Afterwards, the MT restarts when the secondary is evolving off the MS. Finally, the BH grows to a mass of ∼3.



Fig. <ref> shows distribution of masses for BHs and companions. The leftmost peak (both for BHs and companions distributions) is composed of BHs formed through  route, whereas, heavier BHs and companions represent  route. The maximal BH masses, naturally, exceed these for wide BHBs, or expected from single star evolution <cit.>. The heaviest companions are usually accompanied by the heaviest BHs what allows to avoid extreme mass ratios, which in our calculations are assumed to lead to dynamical instability during MT.

Models  and  predict larger populations of MTBHBs than  model. In lower metallicity environments it is easier to produce heavier BHs, thus companions Roche lobes are relatively smaller. Also, the nuclear expansion is slower what allows for longer phases of stable MT.

Typical initial center-of-mass velocities (v_ COM,0) are smaller than 20 km/s. It is a result of the fact that both heavy BHs (M≳40) and heavy NSs (M>1.8) present in MTBHBs formed through channels  and , respectively, and have low, or negligible NKs. The exception are models with higher average NKs ( and ) which give much higher v_ COM,0 (50–100 km/s). We note, that velocities around 100 km/s are also attainable in models with standard NK distribution (e.g. ) in systems with lighter BHs (M_ BH≈6).




 §.§ Mergers




Mergers are a frequent outcome in population synthesis calculations of isolated binaries due to a failed CE ejection. Our results show that most typically (more than 95% of cases) mergers occur when the components are on MS, HG, or are HeS (Fig. <ref>). Although our results are not representative for the entire binary star population, as we include only binaries with initial primary's mass in the range 10150, we include all potential progenitors of SBHs from merging binaries for which +≥20. We note that even if some lower-mass binaries may produce a BH predecessor after a merger, the poorly understood merger physics does not allow to include such cases in the following analysis. Throughout the paper we use a designation mSBH to distinguish SBHs formed from merger products from these originating from disrupted binaries.



In order to include merger products in the BH population we have implemented a simplified formalism of Olejak et al. (in prep.; see Sec. <ref>). After merger, the products were evolved as single stars till the formation of a compact object, or reaching the age of 15. This way, we were able to roughly estimate the population of mSBHs. According to our simulations, the minimal ZAMS mass that produce a BH in single-star evolution is between 19–22 depending on the metallicity. Similarly, the minimal zero-age helium main-sequence mass that produce a BH is ∼9 <cit.>. In a typical case (MS+MS and HG+MS mergers) there are no interactions prior to merger.

The predicted numbers of BHs for all tested models are provided in Tab. <ref>. The model dependence is very low with ∼3.66.67 BHs originating from mergers, except  and  models where the contrast is significant (1.18 and 1.67, respectively) due to different (higher, or lower, respectively) number of massive stars on ZAMS.



The mass distribution (Fig. <ref>) shows a major peak at ∼7.5 (∼15 for  and  models). The distribution noticeably differs from distributions of BH masses obtained through other formation channels. Particularly interesting is the presence of BHs inside the Mass Gap (∼2.5 – 5) were no compact objects have been detected through observations <cit.>. These BHs mainly originate from mergers of NSs with HeS thus exceeding the maximal mass for a NS (M_max,NS). In a typical case, the components masses on ZAMS are ∼10 and ∼9.6. The primary, being slightly heavier, evolves faster and in ∼24 fills its Roche lobe while expanding on a HG. The mass transfer leads to a mass reversal with primary becoming a 2.2 HeS, whereas secondary grows to 14 still remaining on a MS. In about 4 primary explodes as an electron capture supernova and forms a NS with a mass of 1.26 and negligible NK. The separation is 770 at the moment. Secondary needs additional 1 to fill its Roche lobe during a CHeB phase. The mass ratio is extreme with the donor being ten times heavier then the accretor, what leads to a CE and, in a consequence, merger of both stars. Finally, according to our prescription a 2.9 BH forms. Such low-mass free-floating BHs can be detectable by microlensing methods <cit.>. We note that the formation of low-mass BHs through this channel heavily depends on our assumptions concerning mergers (see Sec. <ref>) and, particularly, the fraction of the HeS which is actually accreted onto the NS.

Another striking feature of the distribution is the presence of high mass tails exceeded up to ∼50, ∼80, or ∼130 for , , and  models, respectively (Fig. <ref>). Typically, in  model, progenitors of such massive BHs are formed from a merger of a ∼25 HeS with a 46 MS star which leaves a ∼43 HeS as an outcome. Such a massive HeS collapses directly to a BH with ∼40 mass. Evolution for massive BHs in other metallicities are analogous. These massive stars, although rare, may produce significantly stronger microlensing events than those with typical masses (∼7–8). W note that <cit.> obtained even heavier single BHs as a result of stellar mergers using different prescription for post-merger stellar parameters. In contrast to our results, their heavies BHs (up to ∼145 for metallicity Z=0.0001) originate from mergers of two MS stars.

BHs originating from mergers will have typically very low initial velocities (). According to our results more than 75% of the post-merger BHs will have no additional velocity to that resulting from the movement in the galactic potential (≈0). Only BHs originating from double NS mergers may obtain high  (above 100), however, they constitute only <0.1% of all post-merger BHs. We note, that asymmetries in the merger process and mass ejections may be a source of significant NKs, but are not involved in this study.









§ APPLICATION: MICROLENSING OF BULGE STARS BY ISOLATED BHS


Gravitational microlensing effect is an amplification of the light of a distant object by another object (luminous, or not) aligned in the line-of-sight toward the source <cit.>. Unfortunately, lens mass measurement, which is necessary to conclude on its stellar, or BH nature, is typically prone to parameter degeneracies <cit.> and BH lenses can be found only in particular circumstances <cit.>. On the other hand, population studies allow us to estimate the expected number of microlensing events with BH lenses, what we do in the following analysis. Other similar attempts include studies by <cit.> and <cit.>.

Due to the high star number density and relative proximity, the Galactic bulge is a frequent direction of observations aimed at detecting microlensing events <cit.>, therefore, here we also concentrate on this particular direction.

We make a series of simplifying assumptions. First of all, we assume that the source is always located in the bulge, which is represented by a flat disk of size equal θ_ bulge=31 square degrees, possessing N_⋆, bulge=1.5×10^8 stars. The tangential velocity distribution (radial velocity can be ignored in microlensing) is on average zero and has a dispersion σ_ bulge,z=σ_ bulge,y=80 km/s in both directions: parallel (y-axis) and perpendicular (z-axis) to the galactic plane <cit.>.

Secondly, we assume that a BH lens is always located in the galactic disk, which has a total mass of 5×10^10 <cit.> and its stellar number density (we assume that BHs have the same distribution as stars in general) is described by 


    ρ_ star∝ e^|z|/0.3 kpc-R/2.6 kpc,

where z is the distance from the galactic plane and R is the distance from the galactic centre in the galactic plane, following <cit.> and <cit.>. 
We assume that components of the tangential velocity of stars in the disk are on average v_ disk,mean,z=0 in the direction perpendicular to the galactic plane and v_ disk,mean,y=200 km/s in the galactic plane. The tangential velocity dispersion is assumed to be σ_ disk,z=40 km/s perpendicular to the galactic plane and σ_ disk,y=55 km/s parallel to galactic plane.

In contrary to other studies, here we use the BHs from disrupted binaries and stellar mergers only, as BHs in binaries represent small fraction (<20%) of the total population (Tab. <ref>). Additionally, microlensing by binary lens is more complex (caustic crossing, high amplification, etc.) and will be investigated in a separate study. Also, due to the assumed 100% binary fraction for stars heavier than 10 on ZAMS, we do not expect any BHs formed through single-star evolution.

Stellar mergers may potentially constitute a significant part of the BH population (see Sec. <ref>), so it is necessary to include them in the calculation of the microlensing rate. In such a case, not only the number of BHs is important, but also their masses and proper motions. We calculate the final BH mass following the formulas of Olejak et al. (in prep.; see Sec. <ref> for the post merger mass and evolutionary phase). As far as the proper motions are concerned, we assume that the  of the post merger star is equal to the pre-merger center-of-mass velocity of the binary, i.e. merger occurs without a kick.

We use the Monte Carlo method to sample the spacial and velocity distributions of both lenses and sources. In our analysis we include only lenses which were localized between the observer, assumed to be located at x_⊙=8.3 kpc from the Galactic centre, and the Galactic bulge, ignoring less likely lenses and sources from the far Galactic disk.


An angular Einstein radius of a lens depends on its mass (M_L) and distances between the observer and the lens (D_L), or the source (D_S):

    θ_E = √(κ M_Lπ_rel),     π_rel = 1/D_L - 1/D_S,

where κ=4G/c^2≈ 8.1  mas/. 
If the relative proper motion between the lens and the source is μ⃗_rel = μ⃗_L-μ⃗_S, the source crosses the Einstein radius in time called Einstein Radius crossing time, or simply, event's time-scale, t_E = θ_E/|μ_rel|.


The relative proper motion is computed using the velocities of the lens and the source as:

    μ_ rel,y/z=v_ rel,y/z/x_⊙-x,

and

    v_ rel,y/z=v_ BH,y/z-v_ Earth,y/z
    
        -(v_ source,y/z-v_ Earth,y/z)(1-x/x_⊙),


where v_ rel,y/z is a component of the relative velocity along y-/z-axis. v_ BH,y/z, v_ source,y/z, and v_ Earth,y/z, are the velocities of the BH, the source, and the observer, along y-/z-axis, respectively. The observers velocity was assumed to be v_ Earth,y/z=230/15.5 km/s. x=x_⊙-D_ L is the distance of the lens from the Galactic centre. 

Finally, we calculated an estimated number of microlensing events (E(N_ ML)) during t=1 yr as (assuming θ_ E is small), 

    E(N_ ML)=∑_i (2θ_ E,iv_ rel,i/x_-x_ i t + πθ_ E,i^2)/Ω_ bulgeN_⋆, bulge,

where, θ_ E,i is the Einstein radius generated by i-th BH, v_ rel,i=√(v_ rel,x,i^2+v_ rel,y,i^2) is the relative tangential velocity of the i-th BH, and Ω_ bulge is the considered area of the bulge in steradians. The summation goes over all the included BHs, i.e localized between the observer and the bulge.

If we assume that the stellar mass of the Galactic disk is ∼5×10^10 <cit.> and we take the number of stars in the bulge as N_⋆, bulge= 1.5× 10^8 from the number of monitored stars by the OGLE-III survey covering 31 sq.deg of the bulge <cit.>, we predict that there might be as many as ∼14 yr^-1 microlensing events due to BHs.
In the so-far analysed OGLE-IV data <cit.>, which covers 160 sq.deg and monitors 400 million sources, the number of expected events due to BHs grows to ∼26 yr^-1. 
Here we include only events with characteristic time t_E in range 1–1000 days which includes most of the events (>95%). Although BHs from disrupted binaries have higher velocities and, therefore, higher chance for microlensing, BHs from disrupted binaries (dSBH; Sec. <ref>) dominate the rate (∼55% of events) due to higher masses and greater abundance.



For a given source and lens distances and their velocities we computed θ_E and μ_rel and hence the distribution of time-scales (t_E, Fig. <ref>). The typical values of t_E are on average ∼97 days for BH originating from disrupted binaries (dSBH; Sec. <ref>) and ∼108 days for these originating from stellar mergers (mSBH; Sec. <ref>. 

Events with such time-scales are detectable in current microlensing surveys with efficiency of more than 75% <cit.>. This means in the OGLE-III data covering 8 years of monitoring of 150 million of stars there should be about 84 events due to BHs. Since, as shown above, these tend to have longer time-scales, the annual Earth motion (the parallax effect) is likely to cause the distortion to the standard microlensing light curve, therefore, it is possible that some fraction of those 84 events are not in the <cit.> sample of standard events. Indeed, <cit.> have identified 59 OGLE-III events with significant parallax effect, among which 13 are likely due to NSs and BHs. The improved analysis of <cit.> has increased the number of NS and BH lenses to 18 in that sample. The remaining undetected BH lenses might still await discovery in the OGLE-III data with moderate and small annual parallax signatures. 




§ DISCUSSION




 §.§ Detectability of single BHs


Although the predicted number of SBHs in MWEG (both dSBH and mSBH) is very significant (3.27 – 2.58), these objects still evade detection. Already <cit.> estimated that BHs in the Galaxy should be counted in millions. <cit.> estimated that there are even 1.4×10^9 BHs in the Galaxy. <cit.> calculated the BH population to be 1.8×10^8. They included the changes in SFR, but excluded binary stars from their analysis. It was proposed that SBHs may be detected when interacting with interstellar matter <cit.>, however, the estimated X-ray emission from an accreting SBH is very small <cit.>. Such faint sources may be potentially detected only in our Galaxy. <cit.> suggested that the accretion onto a SBH from interstellar medium may not be necessarily spherical and result in the formation of an accretion disk. In such disks the same instabilities may arise as in X-ray Novae resulting in a transient behavior with outburst luminosities reaching ∼10^38.

SBHs can be also detected as a gravitational microlenses <cit.>. As shown in this paper, SBHs may have very high velocities what increases a chance for microlensing event. On the other hand, a significant part of DCOs ( channel) have separations small enough to act as a single, very massive lens, with the mass as the sum of the masses of individual components. In combination with small proper motions of such binaries, they may generate long-lasting (t_E∼1) events, however, we note that undisrupted compact DCOs make less than 15% of all BHs in the Milky Way.

Although none SBHs were detected so far, a few observational methods were proposed. The predicted number of observations in future X-ray missions is strongly dependant on the BH velocities and accretion models. Future missions may help to improve the constrains on this physics. Nevertheless, detectability of SBHs in X-rays is strongly spoiled by background AGNs, hard coronal emitters and CVs <cit.>. Another possibility, to which more attention is put in this work, is the detection of SBHs with microlensing <cit.>. First candidates were proposed in <cit.> and <cit.> and were followed-up with X-ray observations, however, no signal from ISM accretion was detected <cit.>. Recently, several candidates were proposed <cit.>, but the lens's mass estimates are heavily degenerated due to the lack of measurement of relative velocities. Consequently, WDs or NSs cannot be fully ruled out.  For the currently on-going events, there is an opportunity for degeneracy breaking thanks to astrometric measurements from optical interferometry <cit.> or Gaia<cit.>, however, this will only be possible for the brightest events. In near future with the LSST it should be possible to observe thousands of microlensing events due to BHs and, therefore, to confront the stellar evolution predictions with the observed parameters of Galactic BHs. 




 §.§ Detectability of BHs in binaries


All dynamically confirmed galactic BHs were detected in XRBs. We note that a BH candidate was also detected  in a non-interacting binary <cit.> and many more are expected to be seen by Gaia <cit.>. In XRBs a compact object accretes mass from its companion which results in a formation of an accretion disk and production of highly energetic radiation. If the companion is observable, radial velocities can be measured and the mass function can be estimated. Unfortunately, the majority of XRBs containing BHs are transient systems, thus are visible only during outbursts and only a few recurring systems were observed. What is more, typical transient systems are characterized by low mass donors (M_don≲2), therefore, undetectable from extra-galactic distances and outshined by the disk during an outburst. In 19 XRBs, BHs were dynamically confirmed <cit.>, however, there are many more candidates <cit.>.  

Gravitational waves emission from double BH mergers may be perceived as a detection of BHs. Up-to-now, 10 binary BH merger events were detected by LIGO/Virgo <cit.>. Using different pipeline <cit.> found 6 additional sources in publicly available data from the second observing run. Notably, in spite of one double NS merger detection <cit.>, no BH+NS mergers were detected. The masses of pre-merger BHs range from ∼7.7 to ∼50.6 with the total binary mass as high as ∼85.1 <cit.>. Many tens of detections are expected from the next observational runs <cit.> and several candidates have already been found <cit.>.

BHs in binaries may be detected also through microlensing. They may lens differently depending on the compactness of the binary. If the separation is significantly smaller than the total mass Einstein radius, the system will act as a single lens with a mass equal to the sum of the components' masses. On the other hand, if the system is very wide (many Einstein radii), a BH will act as a single lens unaffected by a remote companion. In-between these two options, we have systems that act as a binary lens with all the effects connected with caustics crossing. Up-to-now, no microlensing events with BH lenses were confirmed. Although there are several candidates <cit.>, parameter degeneracies make it impossible to infer mass, thus distinguish between a BH and other low-luminosity objects like WDs.




 §.§ Comparison with previous estimates of microlensing rates by BHs


Our prediction of ∼14 microlensing events due to BHs in a year stays in contrast to previous estimates. <cit.> analytically calculated that only 1% of microlensing events in the direction of the Bulge are due to BH lenses, what gives about 30 events in OGLE-III (3-4 per year). There are several reasons for this contradiction. First of all, <cit.> concentrates on BH lenses in the bulge, whereas our lenses are in the disk. π_ rel (see Eq. <ref>) is typically smaller for bulge lenses, thus also the Einstein radius, which is proportional to the probability of microlensing. Secondly, <cit.> assumed 40 as a minimal ZAMS mass producing a BH remnant, whereas in our calculations we take ∼22 (for solar metallicity). Consequently, they obtained much less BHs from the same initial population. Moreover, BHs in <cit.> follow the velocity distribution of stars, whereas in our simulations BHs obtain additionally a kick (including NK, Blaauw kick, and the kick from binary disruption). Thirdly, the binary systems were ignored in the work of <cit.>, while here we show that the binary systems could potentially contribute significantly to the population of lensing BHs.

A more recent work, <cit.> used the same population synthesis  code to estimate the microlensing rate events due to BHs. However, they focus on astrometric microlensing events, i.e. events in which not only photometric observation is possible, but also astrometric signal can be detected with Gaia space mission. They estimated there should be a few such events per year. About an order of magnitude of difference between our results stems from the fact that only ∼5% of bulge stars (5×10^6 in their work) are bright enough to potentially generate a detectable astrometric microlensing event. In contrast, we included all bulge stars as potential sources (1.5×10^8), thus 30 times more, which results in significantly higher estimated rate. Moreover, <cit.> used averaged, or typical values to calculate the rates, whereas in this paper we populate the Galaxy with stars, evolve them to the formation of a BH and calculate the rates directly.

We note that in our approach to estimate the microlensing rates we made a number of simplification. Especially, we used only one model from 8 analysed in this papers. Additionally, we chose a simple model of BH distribution in the Galaxy. For example, model of <cit.>, which was used in <cit.>, gives microlensig rate by BH of ∼9 yr^-1, so nearly two times smaller than the rate obtained for the distribution described by Eq. <ref>. We also note that the rate may change by a factor of ∼5 depending on the chosen evolutionary parameters, e.g. IMF, natal kicks, metallicity. 



 §.§ Hertzspurng Gap donors in common envelope evolution


lcrlrlrlrl
    
    Number of BHs in the Milky Way Equivalent Galaxy (model A)

     Model     N_ BH,tot     2cN_ dSBH     2cN_ BHB     2cN_ BH,DCO     2cN_ mSBH
        1.08     5.47     (54.2%)     1.56     (1.5%)     6.16     (6.1%)     3.87     (38.2%) 

        1.28     3.57     (29.9%)     2.06     (1.7%)     2.27     (18.8%)     5.87     (49.6%) 

        1.28     3.77     (30.7%)     1.56     (1.2%)     2.17     (17.4%)     6.17     (50.6%) 

        1.08     6.67     (65.0%)     3.46     (3.3%)     8.16     (8.0%)     2.47     (23.6%) 

        1.08     6.47     (63.8%)     6.04     (0.1%)     2.15     (0.2%)     3.67     (35.9%) 

        1.08     6.57     (64.3%)     5.54     (0.1%)     8.74     (0.1%)     3.67     (35.6%) 

        2.58     1.48     (56.0%)     3.26     (1.3%)     1.67     (6.4%)     9.17     (36.4%) 

        3.17     1.67     (50.9%)     5.65     (1.8%)     1.96     (6.0%)     1.37     (41.3%) 

    The same as Tab. <ref>, but for a models in which HG donors were allowed to survive the CE phase (see Sec. <ref> for details). 



It is still not known if binaries undergoing CE events with HG donors will merger due to the lack of clear core-envelope boundary, or they can survive, because the boundary develops during the HG phase. In order to check the importance of our assumption (model B - CE with HG donor always leads to a merger), we have tested also an opposite case (model A) where binaries are allowed to survive the CE phase even if donor is on its HG. The results are presented in Tab. <ref>.

The difference in the total number of BHs is small (≲10%), as most of the BHs form without a CE phase during the binary evolution. More pronounced is the higher number of DCOs (particularly, the close ones), which may exceed two orders of magnitude in comparison to model B for higher NK models ( and ). A part of close DCOs progenitors evolve through a CE phase and in model A donors are frequently HG stars.

The number of mergers leading to the formation of a BH is also affected. In model B, CE events with HG donors are treated as mergers, what is not always a case in model A. In a consequence, there are significantly fewer (up to 33%) BHs originating from stellar mergers in model A than in model B.

As the total number of BHs is affected only slightly by the treatment of HG donors in CE events, we predict a negligible impact on the microlensing observations. However, it may be significant if events caused by close systems are considered. 




 §.§ Binary fraction


In this work, we assumed that all stars heavier than 10 are born in binaries. Such a high binary fraction (100%) for massive stars is supported by observations, which suggest their binary fraction to be higher than 90% <cit.>. For lower mass stars we assumed equal number of binaries and single stars on ZAMS (binary fraction equal 50%).

Nonetheless, we agree that the binary fraction may not change in such a drastic way when the mass of the primary increases. More probably, the binary fraction increases continuously. Therefore, a fraction of massive stars, which are heavy enough to produce BHs (M_ ZAMS≳20) may actually form without companions. If included in our simulations, these stars could increase the relative fraction of single BHs in the entire BH population. However, the total number of BHs is expected to be lower due to the fact that flat mass ration distribution in binaries gives higher average mass of two stars in a binary, thus more BH progenitors, than two single stars following the IMF relation. We note, that mergers of binary systems of two massive stars (both heavy enough to form a BH in single stars evolution), may decrease the total number of BHs, but mergers of NS progenitors (M_ zams≈8–20) may become heavy enough to form a BH in posterior evolution, thus increasing the total number of BHs.

Many of the massive binaries may be part of triples and higher-order systems <cit.>. Even ∼50% of massive (OB type) stars may exist in triples <cit.>, but the inner binary in a hierarchical triple system may evolve effectively isolated from the third star, thus estimating the significance of higher-order systems on BH populations is complicated. Of particular interest are mergers induced by Kozai-Lidov mechanism in stellar triples <cit.>. However, higher-order systems are not understood yet well enough to be included in population synthesis modeling.







§ SUMMARY

   
In this paper we analysed the general properties of synthetic BH populations in different stellar environments represented by models which differ in the metallicity, initial parameter distributions, and natal kick prescriptions. We note that the results are applicable for further studies like predictions for present and future survey missions, or in-depth analysis of specific BH populations (e.g. the Milky Way galaxy). We particularly focused on BHs originating from disrupted binaries and stellar mergers, which were frequently omitted in previous studies. We find that those BHs constitute a majority of the total BH population.

Particularly, we show that single BHs dominate the total population of BH (≳83% of all BHs), even though massive stars form predominantly in binaries. Both binary disruptions and stellar mergers are important and the predicted number of BHs is only slightly affected (up to a factor of ∼3) by a chosen model. Although BHs in binaries constitute only a small part (≲17%) of the population, their number is heavily dependent (about two orders of magnitude) on the adopted model parameters (especially natal kick prescription and star formation history), therefore, new observations may significantly help to constrain these parameters and better understand evolutionary processes <cit.>.

Using our results, we calculated the expected rate of microlensing events with BH lenses in the direction of the Galactic bulge. We expect that as many as ∼14 such events per year with average crossing times around 100 days in the OGLE-III footprint and about 26 in the OGLE-IV data. Only some of these events may be observable from Earth due to technical limitation (e.g., low luminosity, extinction). In this estimate we have neglected BHs which remained in binaries, because, low fraction (≲17% of all BHs) and slow velocities (typically, ≲20 km/s) allow to suppose that their influence is small.

A grid of 54 models (including 8 main models analysed in the paper) is available in a free-access database[https://universeathome.pl/universe/bhdb.php]. In our future work we plan to utilize these datafiles in order to study microlensing by single BHs in the Galaxy and provide detailed predictions for surveys like Gaia, or Einstein Telescope.




We are thankful to the anonymous referee who helped to improve the paper and thousands of volunteers, who took part in the Universe@Home project[https://universeathome.pl/] and provided their computers for this research. We acknowledge the Polish NCN grant HARMONIA 2015/18/M/ST9/00544 to ŁW, OPTICON H2020 EC grant no. 730890, MNiSW grant DIR/WK/2018/12 and SONATA BIS 2 DEC-2012/07/E/ST9/01360. G.W. is partly supported by the President’s International Fellowship Initiative (PIFI) of the Chinese Academy of Sciences under grant no.2018PM0017 and by the Strategic Priority Research Program of the Chinese Academy of Science Multi-waveband Gravitational Wave Universe (Grant No. XDB23040000). KB and GW acknowledge NCN grant MAESTRO 2018/30/A/ST9/00050. MC and JK acknowledge support from the Netherlands Organisation for Scientific Research (NWO). We acknowledge the Polish NCN grant PRELUDIUM 2017/25/N/ST9/01253 to KR. This work is partly supported by the National Key Program for Science and Technology Research and Development (Grant No. 2016YFA0400704), and the National Natural Science Foundation of China (Grant No. 11690024 and 11873056).


apj




















