Black hole-neutron star mergers are unlikely multi-messenger sources

The promise by the LIGO/Virgo/Kagra (LVK) collaboration to detect black hole--neutron star (BH--NS) mergers via gravitational wave (GW) emission has recently been fulfilled with the detection of GW200105 and GW200115. Mergers of BH--NS binaries are particularly exciting for their multi-messenger potential, since the GW detection can be followed by an electromagnetic (EM) counterpart (kilonova, gamma-ray burst, afterglow) that can reveal important information on the equation of state (EOS) of NSs and the nature of the BH spin. We carry out a statistical study of the binary stars that evolve to form a BH--NS binary and compute the rate of merger events that can be followed by an EM counterpart. We find that $\gtrsim 50\%$ of the mergers can lead to an EM counterpart only in the case BHs are born highly spinning ($\chi_{\rm BH}\gtrsim 0.7$), while this fraction does not exceed about $30\%$ for stiff NS EOSs and a few percent for soft NS EOSs for low-spinning BHs ($\chi_{\rm BH}\lesssim 0.2$), suggesting that a high rate of EM counterparts of BH--NS would provide support for high natal BH spins. However, the possibilities that BHs are born with near-maximal spins and that NS internal structure is described by a stiff EOS are disfavored by current LVK constraints. Considering that these values only represent an upper limit to observe an EM counterpart due to current observational limitations, as in brightness sensitivity and sky localization, BH--NS mergers are unlikely multi-messenger sources.


INTRODUCTION
Mergers of black hole-neutron star (BH-NS) binaries are particularly interesting since they can produce an electromagnetic (EM) counterpart associated with the gravitational (GW) signal at merger. The condition for an EM counterpart to occur is that during the merger phase, which follows the inspiral due to GW emission, the NS does not directly plunge into the BH, but rather is tidally disrupted. For these systems only, a post-merger phase is expected, during which NS matter debris can be ejected or accreted onto the BH producing a luminous event. BH-NS mergers that produce an EM counterpart can provide constraints on the BH spin and accretion process and unique information on the nuclear equation of state (EOS) (e.g., Pannarale et al. 2011;Tsang et al. 2012;Pannarale et al. 2015;Foucart et al. 2018;Ascenzi et al. 2019;Hinderer et al. 2019;Zappa et al. 2019;Fragione & Loeb 2021;Tiwari et al. 2021).
Hundreds of merging BH-NS binaries are expected to be detected in the next few years. The second Gravitational Wave Transient Catalog by the LIGO/Virgo/Kagra (LVK) Collaboration reports tens of BH-BH mergers, two NS-NS mergers, and a candidate BH-NS merger (GW190426; Abbott et al. 2020). Recently, the first two confirmed BH-NS mergers, GW200105 and GW20115, from the third observational run have been publicly released (Abbott et al. 2021a).
The astrophysical origin of BH-NS mergers remains highly uncertain. The most promising scenario is represented by BH-NS mergers produced as a result of the evolution of field binaries. Within the uncertainties of stellar evolution models, this channel predicts merger rates consistent with the empirical LVK estimate (e.g., de Mink & Mandel 2016;Giacobbo & Mapelli 2018;Kruckow et al. 2018;Broekgaarden et al. 2021;Shao & Li 2021). BH-NS binaries do not efficiently form in a dense star cluster, resulting in a merger rate several orders of magnitude smaller than LVK estimated rates (e.g., Clausen et al. 2013;Bae et al. 2014;Belczynski et al. 2018;Fragione & Banerjee 2020; Arca Sedda 2020; Ye et al. 2020), which can increase under favorable initial conditions (Rastello et al. 2020). BH-NS mergers from massive triples as a result of the Lidov-Kozai mechanism have been proposed as an alternative channel, but merger rates consistent with the LVK results are obtained only if very low natal kicks for BHs and NSs are assumed (e.g., Fragione & Loeb 2019a,b). Finally, BH-NS binaries can be assembled and merge in AGN disks at high rates, but there are still major uncertainties in the models (e.g., Yang et al. 2020;Tagawa et al. 2021).
The detection of a short gamma-ray burst (GRB) (GW170817A; Abbott et al. 2017a), of an ultraviolet-opticalinfrared transient (AT2017gfo; Abbott et al. 2017c), and of an off-axis jet (e.g., Margutti et al. 2017;Troja et al. 2017) associated with GW170817 (Abbott et al. 2017b) has opened the long-awaited era of multi-messenger astronomy for double NS mergers. The question whether BH-NS mergers are typically expected to have EM counterparts remains still to be answered. Many detailed follow-up observations for the candidate and confirmed LVK BH-NS mergers have shown no associated EM counterpart (Hosseinzadeh et al. 2019;Goldstein et al. 2019;Coughlin et al. 2020;Thakur et al. 2020;Anand et al. 2021;Alexander et al. 2021;Kilpatrick et al. 2021). On one hand, the lack of detections could simply be due to the fact that present searches are not sensitive enough (BH-NS mergers can be observed via GW emission at larger distances with respect to binary NSs) and sky localization is not optimal. On the other hand, there is the possibility that EM counterparts are intrinsically missing because the NS typically plunges into the BH, leaving behind no debris to accrete.
In this Letter, we use population synthesis of BH-NS systems from field binaries to study the rate of expected EM counterparts. We adopt different models for natal kicks, efficiencies for common-envelope ejection, BH birth spins, and NS EOSs. We show that a significant fraction of BH-NS mergers is expected to be followed by an EM counterpart only if BHs are born highly-spinning and NSs have stiff EOSs, both of which are currently disfavored by LVK constraints.
This Letter is organized as follows. In Section 2, we describe our models and methods to build a cosmologicallymotivated population of BH-NS mergers that can be followed by an EM counterpart. In Section 3, we discuss the expected distributions and rates of EM counterparts to BH-NS mergers. Finally, in Section 4, we summarize our findings and draw our conclusions.
2. METHOD 2.1. Binary population synthesis We sample the initial mass m 1 of the primary from a Kroupa (2001) in the mass range [20 M -150 M ], appropriate for BH progenitors. To determine the initial secondary mass (m 2 ), we adopt a flat mass ratio distribution (Sana et al. 2012;Duchêne & Kraus 2013;Sana 2017). Finally, we extract orbital periods (in days) from f (log 10 P) ∝ (log 10 P) −0.55 in the range [0.15-5.5] (Sana et al. 2012) and assume a thermal distribution for the eccentricity. The binaries we sample are evolved using the stellar evolution code BSE (Hurley et al. 2000(Hurley et al. , 2002. We use the latest version of BSE from Banerjee et al. (2020a), which includes the most up-to-date prescriptions for stellar winds and remnant formation, and produces remnant populations consistent with those from STARTRACK (Belczynski et al. 2008).
We model the distribution of natal velocity kicks (due to recoil from an asymmetric supernova explosion) imparted to compact objects at birth as a Maxwellian distribution with velocity dispersion σ. In our main model, we assume σ = 265 km s −1 , consistent with the distribution inferred from proper motions of pulsars by Hobbs et al. (2005). However, the value of σ is uncertain. For example, Arzoumanian et al. (2002) found a bimodal distribution with characteristic velocities of 90 km s −1 and 500 km s −1 based on the velocities of isolated radio pulsars, while Beniamini & Piran (2016) found evidence for a low-kick population ( 30 km s −1 ) and a high-kick population ( 400 km s −1 ) based on observed binary NSs. Therefore, we also run additional models with σ = 10 km s −1 , 40 km s −1 , 150 km s −1 . For BHs, we sample natal kicks from the same distribution as for NSs, scaling them down with increasing mass fallback fraction (Repetto et al. 2012;Janka 2013). We also self-consistently keep track of the spin-orbit misalignment (I BH−NS ) produced as a result of natal kicks, by extracting them from BSE and computing the tilt of the binary orbit (whenever the orbit remains bound, as outlined in Sect. 2 in Fragione et al. (2021)). We do not model other possible sources of spin-orbit misalignment, as gas torques due to accretion during common-envelope events.
We consider three different models for BH spins. In the first model, we assume that the BH spin is χ BH = 0.9, which is consistent with the prescriptions of the GENEVA stellar evolution code over the mass range relevant for BH-NS mergers (Eggenberger et al. 2008;Ekström et al. 2012). In the second model, we assume χ BH = 0.1, consistent with the results of the MESA stellar evolution code (Paxton et al. 2011(Paxton et al. , 2015. Indeed, MESA includes a treatment for the magnetic field that makes the outwards angular momentum transport from the core much more efficient by forming a Tayler-Spruit magnetic Observational constraints on the mass-radius relationship compared with NS equilibrium sequences. The dotted-gray line represents constraints (95% confidence region) for low-mass and high-mass NSs from LVK using data from the NS-NS merger GW170817, while the black-dashed line and the pink-dot-dashed line represent constraints (95% confidence region) from NICER and NICER+XMM, respectively, using PSR J0740+6620. dynamo, resulting in a small BH spin. Recent detailed studies of massive stars have suggested that the Tayler-Spruit magnetic dynamo can essentially extract all of the angular momentum, birthing BHs with extremely low spins (Fuller & Ma 2019). Therefore, we also consider a model where the initial spin of BHs is assumed to vanish. For full details see Banerjee et al. (2020b) and Banerjee (2021).

Electromagnetic counterpart
To compute whether a BH-NS merger produces an EM signature, we compute the remnant baryon mass (M rem ) outside the BH after merger. We assume that if M rem > 0 a disk is formed and there is EM emission, otherwise the NS plunges directly into the BH leaving no signatures other than the GW inspiral (e.g., Foucart 2012;Foucart et al. 2013). To estimate the remnant disk mass (in units of the initial mass of the NS), we usê where (α, β, γ, δ) = (0.406, 0.139, 0.255, 1.761) are fitting parameters to numerical relativity simulations (Foucart 2012;Foucart et al. 2013Foucart et al. , 2018, C NS = GM NS /(R NS c 2 ) is the NS compaction, which depends on the EOS, η = Q/(1 + Q) 2 , with Q = M BH /M NS , is the symmetric mass ratio, and (in units is the innermost stable circular orbit (ISCO) radius (Bardeen et al. 1972), with In the previous equations, χ BH, is the aligned component of the BH spin (χ BH, = χ BH cos I BH−NS ) with respect to the orbital angular momentum. For the NS EOS, we consider six different equilibrium sequences, namely APR4 (Akmal et al. 1998), ENG (Engvik et al. 1996), MPA1 (Müther et al. 1987), MS0 (Müller & Serot 1996), SLy4 (Chabanat et al. 1998), and WFF1 (Wiringa et al. 1988). We show these NS equilibium sequences in Figure 1 along with current observational constraints from LVK (Abbott et al. 2018) and NICER(+ XMM) (Miller et al. 2021). Note that, while NICER constraints allow NSs with large radii (stiff EOSs), this portion of the parameter space is excluded by LVK constraints at 95% confidence, favoring softer EOSs.

Distribution of BH-NS mergers with electromagnetic
counterpart To place our BH-NS population in a cosmological context, we assign to each binary a formation time t form by sampling the formation redshift z form from the cosmic star formation history of Madau & Dickinson (2014) Then, we convolve the delay time t delay of the BH-NS mergers with the distribution of formation times and discard the binaries that merge later than the present day. Each BH-NS that is not discarded is then assigned a weight that accounts for the cosmic distribution of metallicity, which we assume is described by a log-normal distribution Π(z, Z), with mean given by (Madau & Fragos 2017) log Z/Z = 0.153 − 0.074z 1.34 (8) and a standard deviation of 0.5 dex (Dvorkin et al. 2015). This weighting procedure provides us with the underlying astrophysical distribution of sources at a given redshift interval per comoving volume.
2.4. Merger rates of BH-NS mergers with electromagnetic counterpart For a specific EOS, we compute the rates of BH-NS mergers with electromagnetic counterpart as where f b = 0.5 is the fraction of stars in binaries (e.g., Raghavan et al. 2010;Duchêne & Kraus 2013;Sana 2017), f IMF = 0.115 is a correction factor that accounts for our truncation of the primary mass distribution ≥ 20 M (assuming a Kroupa (2001) initial mass function), t lb is the look-back time at redshift z 1 . In Eq. 9, Φ is the merger efficiency at a given metallicity of binaries that produce an EM counterpart where M tot (Z) is the total simulated mass at metallicity Z and N merger (z, Z) the total number of BH-NS mergers at redshift z originating from progenitors at metallicity Z, and  is the fraction of merging systems that have an EM counterpart (M rem > 0), assuming a given EOS. To compute BH-NS merger rates (considering both systems with and without EM counterparts), Γ BH−NS (z), we simply impose S (EOS) EM (z, Z) = 1.
3. RESULTS For a given EOS, whether a NS plunges directly onto a BH or is tidally disrupted producing an EM counterpart depends crucially on the BH spin and its orientation with respect to the orbital angular momentum. In Figure 2, the region below each line covers the portion of the mass (M NS -M BH ) parameter space that allows the tidal disruption of a NS by a BH for each combination of NS EOSs and BH spin models considered in this work. We also report the 90% credible interval for the LVK candidate event GW190426, the two confirmed LVK BH-NS mergers GW201005 and GW20015, and the re-analysis by Mandel & Smith (2021) of GW200115 (GW200115 M21). Note that, for the BH masses of interest here, χ BH = 0.9 and χ BH = 0.1 correspond to the GENEVA model and the MESA model, respectively, while χ BH = 0.0 to the model of Fuller & Ma (2019) (see Figure 3 in Baner-  jee 2021). In the top panel, we show the case I BH−−NS = 0 • , that is the BH spin is aligned to the orbital angular momentum. Under this assumption, a BH of ∼ 10 M could lead to a tidal disruption of a NS of a typical mass of 1.3 M for even the softest EOSs here considered. Moreover, we find that the stiffer the NS EOS is and the higher the BH spin is, the larger is the parameter space that leads to the tidal disruption of the NS, eventually producing an EM counterpart. We find that LVK events could have had (within mass uncertainty limits) an EM counterpart in the case of highly-spinning BHs. However, current constraints on the effective spin of these events suggest a very low value for the BH spin, rendering the likelihood of an EM counterpart small. The previous trend is highly affected by the orientation of the BH spin with respect to the orbital angular momentum. As an illustrative case, we plot the critical lines in the case I BH−−NS = 140 • in the bottom panel, which is representative of the median misalignment of the BH spin in GW200115 (Abbott et al. 2021a), even though recent analysis by Mandel & Smith (2021) cast doubt on this result. In this case, the region of the mass parameter space that allows the tidal disruption of a NS shrinks significantly and does not critically depend on the magnitude of the BH spin. In this case, APR4 and WFF1, the two softest EOS here considered, always lead to a direct plunge, while the other EOSs allow the tidal disruption of a NS only for BH masses 6 M . Only MS0, the stiffest EOS we consider, would predict an EM counterpart for GW190426 and GW200115 (within mass uncertainty limits), regardless of the BH spin. Note that, however, the possibility that all BHs are born with near-maximal spins is becoming disfavored (Abbott et al. 2021b) and very stiff EOSs are excluded at 95% confidence by LVK constraints on the NS EOS (Abbott et al. 2018). Moreover, their kilonova brightness would be too faint to be detected for present follow-up search campaigns, justifying the lack of detections so far (e.g., Zhu et al. 2021).
To understand whether BH-NS mergers are typically expected to be multi-messenger sources, we illustrate in Figure 3 the contours (68%, 95%, 99.7%) of component masses of BH-NS mergers that produce an EM counterpart in our cosmologically-motivated model for different EOSs and BH spin models. In this case, we fix σ = 260 km s −1 and α CM = 1 in our population synthesis models; other combinations of σ and α CM lead to qualitatively similar results. Note that we account for the cosmological distribution of inclination angles, I BH−NS , (see Fragione et al. 2021), which is critical to determine whether a NS directly plunges onto or is tidally disrupted by the BH. We also report in Figure 3 the median mass for the LVK candidate event GW190426, the two confirmed LVK BH-NS mergers GW201005 and GW20015, and the re-analysis by Mandel & Smith (2021) of GW200115. Using similar considerations used in the previous plot, we find that a considerable fraction of our cosmologically-motivated population of BH-NS mergers can eventually produce an EM counterpart only in the case BHs are born highly spinning. Under this assumption, current LVK systems would lie in the 95%-density region, with the exception of GW201005 in the case we consider WFF1, our softest EOS. If BHs are born slowly spinning, only MS0, our stiffest EOS, allows the tidal disruption of NSs, and the eventual EM counterpart, in a nonnegligible fraction of the mass parameter space.
Finally, we compute the fractional rate in the local Universe of BH-NS systems that could lead to an EM counterpart with respect to the total rate of BH-NS mergers, for different values of χ BH . We show our results in Figure 4, for different values of σ and assuming α CM = 1. We report our results for different values of α CM in Appendix A. We find that 50% of the mergers can lead to an EM counterpart only in the case BHs are born highly spinning (χ BH 0.7). In the case BH are born with small spins (χ BH 0.2), the fraction of systems that could have an EM counterpart does not exceed about 30% for MS0, our stiffest EOS. These trends appear to be qualitatively the same for different values of α CM and σ, since high natal kicks imply larger spin-orbit misalignment (see Fig. 2 in Fragione et al. 2021), but, at the same time, they may disrupt the binary system, preventing the BH-NS binary to eventually merge, and affect the shape of the delay-time distribution (time between binary formation and merger). As the EOS of NSs gets better constrained, our results predict that a high rate of detected EM counterparts to BH-NS mergers would support the case that BHs are preferentially born with high spins. The possibility that all BHs are born with near-maximal spins is becoming disfavored (Abbott et al. 2021b) and very stiff EOSs are excluded at 95% confidence by LVK constraints on the NS EOS (Abbott et al. 2018). Under these assumptions, we find that the fraction of BH-NS mergers that can possibly have an EM counterpart does not exceed about 10%. In any case, the future measured EM counterpart fraction could be translated sensitively into constraints on the BH spin using the results of our Figure 4. However, note that these values only represent an upper limit of observable EM counterparts, since, for example, they might be too faint to observe their kilonova and sky localization could not be optimal for a possible follow-up observation (e.g., Zhu et al. 2021).
We finally note that in our models we have taken into account only the spin-orbit misalignment produced as a result of natal kicks and have not modeled other possible sources of spin-orbit misalignments, as gas torques due to accretion during common-envelope events. The detailed prescriptions of core-collapse physics may also play a role (Román-Garza et al. 2021). Also, we note that there might be other types of EM counterparts associated with BH-NS mergers, as shockpowered radio precursors (Sridhar et al. 2021) and resonant shattering flares (Neill et al. 2021). 4. CONCLUSIONS BH-NS mergers are very interesting since they could provide crucial information on their origin, the nature of the BH spin, and the NS internal structure. However, this relies on the fact that BH-NS mergers are followed by an EM counterpart, which happens whenever the NS does not directly plunge onto the BH.
We have carried out a broad statistical study of the field binary stars that form merging BH-NS binaries and have evaluated the fraction that can eventually be associated with an EM counterpart. We have considered different NS EOSs and BH spin models, and we have taken into account the uncertainties on the natal-kick magnitudes and efficiency of commonenvelope ejection for compact objects. We have found that 50% of BH-NS mergers can lead to an EM counterpart only in the case BHs are born highly spinning (χ BH 0.7), otherwise (χ BH 0.2) this fraction does not exceed a few percent for soft EOSs and 30% for stiff EOSs. However, current LVK constraints tend to disfavor the scenario where all BHs are born with near-maximal spins (Abbott et al. 2021b) and exclude very stiff EOSs at 95% confidence (Abbott et al. 2018). Considering that these values only represent an upper limit due to current limitations of EM follow-up observations to GW detections, we conclude that BH-NS mergers are un-likely multi-messenger sources.