# Search for a Black Hole Binary in Gaia DR3 Astrometric Binary Stars with
Spectroscopic Data

Ataru Tanikawa Department of Earth Science and Astronomy, College of Arts and
Sciences, The University of Tokyo, 3-8-1 Komaba, Meguro-ku, Tokyo 153-8902,
Japan  Kohei Hattori National Astronomical Observatory of Japan, 2-21-1 Osawa,
Mitaka, Tokyo 181-8588, Japan  Institute of Statistical Mathematics, 10-3
Midoricho, Tachikawa, Tokyo 190-8562, Japan  Norita Kawanaka Center for
Gravitational Physics and Quantum Information, Yukawa Institute for
Theoretical Physics, Kyoto University, Kitashirakawa Oiwake-cho, Sakyo-ku,
Kyoto 606-8502, Japan  Tomoya Kinugawa Institute for Cosmic Ray Research,
University of Tokyo, Kashiwa, Chiba 255-8582, Japan  Minori Shikauchi
Department of Physics, the University of Tokyo, 7-3-1 Hongo, Bunkyo, Tokyo
113-0033, Japan  Research Center for the Early Universe (RESCEU), School of
Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033,
Japan  Department of Physics and Astronomy, the University of British
Columbia, 6224 Agricultural Road, Vancouver, BC, V6T 1Z1, Canada  Daichi Tsuna
Research Center for the Early Universe (RESCEU), School of Science, The
University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan

###### Abstract

We report the discovery of a candidate binary system consisting of a black
hole (BH) and a red giant branch star from the Gaia DR3. This binary system is
discovered from 64108 binary solutions for which both astrometric and
spectroscopic data are available. For this system, the astrometric and
spectroscopic solutions are consistent with each other, making this system a
confident candidate of a BH binary. The primary (visible) star in this system,
Gaia DR3 5870569352746779008, is a red giant branch whose mass is quite
uncertain. Fortunately, albeit the uncertainty of the primary’s mass, we can
estimate the mass of the secondary (dark) object in this system to be $>5.68$
$M_{\odot}$ with a probability of $99$ %, based on the orbital parameters. The
mass of the secondary object is much larger than the maximum neutron star mass
($\sim 2.0$ $M_{\odot}$), which indicates that the secondary object is likely
a BH. We argue that, if this dark object is not a BH, this system must be a
more exotic system, in which the primary red giant branch star orbits around a
quadruple star system (or a higher-order multiple star system) whose total
mass is more than $5.68$ $M_{\odot}$. If this is a genuine BH binary, this has
the longest period ($1352.22\pm 45.81$ days) among discovered so far. As our
conclusion entirely relies on the Gaia DR3 data, independent confirmation with
follow-up observations (e.g. long-term spectra) is desired.

Stellar mass black holes – Astrometric binary stars – Spectroscopic binary
stars

††software: Matplotlib (Hunter, 2007); NumPy (van der Walt et al., 2011);
Astropy (Astropy Collaboration et al., 2013); SciPy (Virtanen et al., 2020)

##  1 Introduction

Stellar-mass black holes (BHs) are the final state of massive stars with
several $10$ $M_{\odot}$ (e.g. Woosley et al., 2002). BHs are not just dark,
especially when they are members of close binary stars. Thus, they have been
discovered as X-ray binaries (e.g. Casares et al., 2017) and gravitational
wave transients (Abbott et al., 2019, 2021; The LIGO Scientific Collaboration
et al., 2021). Nevertheless, since such BH populations are rare, just a
handful of BHs are known. So far, $\sim 100$ BHs have been detected as X-ray
binaries in the Milky Way, (Corral-Santana et al., 2016), while there should
be $\sim 10^{8}$ BHs in the Milky Way (e.g. Shapiro & Teukolsky, 1983; van den
Heuvel et al., 1992). This is because BHs are bright in X-rays only when they
have close companion stars: binary periods of less than about $1$ day.

Great efforts have been made to discover a variety of BHs in binary stars
(hereafter BH binaries). Many spectroscopic observations have reported BH
binaries with periods of $1$–$100$ days (Thompson et al., 2019; Liu et al.,
2019; Rivinius et al., 2020; Jayasinghe et al., 2021, 2022b; Lennon et al.,
2022; Saracino et al., 2022). However, many concerns have been raised for
these reports (Abdul-Masih et al., 2020; El-Badry & Quataert, 2020; van den
Heuvel & Tauris, 2020; El-Badry & Quataert, 2021; Eldridge et al., 2020;
Irrgang et al., 2020; Tanikawa et al., 2020; Safarzadeh et al., 2020;
Bodensteiner et al., 2020; Shenar et al., 2020; El-Badry & Burdge, 2022; El-
Badry et al., 2022a, b). Several BH binaries (Giesers et al., 2018; Shenar et
al., 2022) still survive, despite such harsh environment for BH binary
searchers.

Gaia have monitored more than $10^{9}$ stars and their astrometric and
spectroscopic motions during 34 months (Gaia Collaboration et al., 2016,
2018a, 2021, 2022a), and have published $\sim 3\times 10^{5}$ astrometric and
spectroscopic binary stars in total in Gaia Data Release 3 (GDR3; Gaia
Collaboration et al., 2022b; Holl et al., 2022a, b; Halbwachs et al., 2022).
Before GDR3, many studies have predicted that Gaia discovers a large amount of
compact objects in binary stars, such as white dwarfs (WDs), neutron stars
(NSs), and BHs, from Gaia’s astrometric data (Mashian & Loeb, 2017; Breivik et
al., 2017; Yalinewich et al., 2018; Yamaguchi et al., 2018; Kinugawa &
Yamaguchi, 2018; Shahaf et al., 2019; Shao & Li, 2019; Andrews et al., 2019,
2021; Shikauchi et al., 2020, 2022; Chawla et al., 2022; Janssens et al.,
2022). Starting with Gaia Collaboration et al. (2022b), many research groups
have searched for WD, NS, and BH binaries in spectroscopic binaries (Gomel et
al., 2022; Jayasinghe et al., 2022a; Fu et al., 2022) and astrometric binaries
(Andrews et al., 2022; Chakrabarti et al., 2022; El-Badry et al., 2023a;
Shahaf et al., 2023) just after GDR3.

GDR3 has presented several $10^{4}$ binary stars with both of astrometric and
spectroscopic data. However, previous studies have focused on either of
astrometric or spectroscopic data. In this paper, we first search for BH
binaries from binary stars where both data are available, taking into account
both of astrometric and spectroscopic data. In other words, we first make a
comparison between astrometric and spectroscopic mass functions (see Eqs. (1)
and (3), respectively) to search for BH binaries.

We eventually find a promising BH binary candidate whose source ID is GDR3
5870569352746779008. After we posted this work to arXiv, El-Badry et al.
(2023a) have independently pointed out that the BH binary candidate is
promising, and El-Badry et al. (2023b) have confirmed it as a genuine BH
binary by follow-up observations. This shows that our search is helpful and
efficient to narrow down BH binary candidates. Although we recognize that El-
Badry et al. (2023b) call it “Gaia BH2”, we call it “BH binary candidate” in
this paper. This is not because we disagree with their confirmation, but
because we regarded it as a BH binary candidate when we posted this work to
arXiv (September 2022).

The structure of this paper is as follows. In section 2, we describe how to
select a sample of binary stars from GDR3, and how to list up BH binary
candidates. Finally, we find one BH binary candidate. In section 3, we analyze
the BH binary candidate in detail. In section 4, we discuss the BH binary
candidate, comparing it with BH binary candidates listed by previous studies.
In section 5, we summarize this paper.

##  2 Sample selection

###  2.1 Search for BH binaries with $m_{2}>3M_{\odot}$

We select GDR3 binary stars with astrometric and spectroscopic data (Gaia
Collaboration et al., 2022b). There are three types of such binary stars. The
orbital solutions of the first type are obtained from astrometric and
spectroscopic data. They have a nss_solution_type name of “AstroSpectroSB1” in
the non-single star tables of GDR3 (nss_two_body_orbit). We call them
AstroSpectroSB1 binary stars. The second type has an orbital solution derived
only from astrometric data, and additionally has the total amplitude in the
radial velocity time series called “rv_amplitude_robust”. Such binary stars
have a nss_solution_type name of “Orbital”, and satisfy the following two
conditions. First, they are bright stars; they have Gaia RVS magnitude less
than and equal to 12. Second, their radial velocities are computed more than
twice. For the third type, binary stars have two nss_solution_type names of
“Orbital” and “SB1” independently. Such binary stars also have a
non_single_star value of 3. Hereafter, the second and third types are
collectively called Orbital binary stars simply. We can extract such a sample
of binary stars from GDR3 with following ADQL query:

    select nss.*, gs.*
    from gaiadr3.nss_two_body_orbit as nss,
    gaiadr3.gaia_source as gs
    where nss.source_id = gs.source_id
    and (nss.nss_solution_type = ’AstroSpectroSB1’
         or (nss.nss_solution_type = ’Orbital’
             and gs.rv_amplitude_robust IS NOT NULL)
         or gs.non_single_star = 3)

The line numbers 5, 6-7, and 8 in the above ADQL query try to pick up the
first, second and third types, respectively. However, the line number 8 picks
up binary stars not only the third type binaries but also many other binaries,
for example, binary stars with nss_solution_type names of “acceleration7” and
“SB1”. We exclude them later. Finally, the total number of binary stars is
64108 consisting of 33467 “AstroSpectroSB1” and 30641 “Orbital” binary stars,
where the numbers of the second and third types are 30629 and 12,
respectively.

We search for BH binary candidates from the above sample, using astrometric
and spectroscopic mass functions ($f_{\rm m,astro}$ and $f_{\rm m,spectro}$,
respectively). We express these mass functions as follows:

$\displaystyle f_{\rm m,astro}$
$\displaystyle=(m_{1}+m_{2})\left|\frac{m_{2}}{m_{1}+m_{2}}-\frac{F_{2}/F_{1}}{1+F_{2}/F_{1}}\right|^{3}$
(1) $\displaystyle=1\left(\frac{a_{1}}{{\rm
mas}}\right)^{3}\left(\frac{\varpi}{{\rm mas}}\right)^{-3}\left(\frac{P}{{\rm
yr}}\right)^{-2}\;[M_{\odot}],$ (2)

and

$\displaystyle f_{\rm m,spectro}$
$\displaystyle=(m_{1}+m_{2})\left(\frac{m_{2}}{m_{1}+m_{2}}\right)^{3}$ (3)
$\displaystyle=3.7931\times 10^{-5}\left(\frac{K_{1}}{{\rm
km~{}s^{-1}}}\right)^{3}\left(\frac{P}{{\rm yr}}\right)$
$\displaystyle\times(1-e^{2})^{3/2}\sin^{-3}i\;[M_{\odot}],$ (4)

where $m_{1}$ and $m_{2}$ are the primary and secondary stars of a binary
star, $F_{2}/F_{1}$ is the flux ratio of the secondary star to the primary
star, $a_{1}$ is the angular semi-major axis of the primary star, $K_{1}$ is
the semi-amplitude of the radial velocity of the primary star, and $\varpi$,
$P$, $e$, and $i$ are the parallax, period, eccentricity, and inclination
angle of the binary star, respectively. We define a primary star as a star
observed by astrometry and spectroscopy, and a secondary star as a fainter
star than the primary star. The secondary star is an unseen star if
$F_{2}/F_{1}=0$. We can get $a_{1}$, $\varpi$, $P$, $e$, and $i$ from
astrometry, and $K_{1}$ from spectroscopy. We have to remark that $f_{\rm
m,spectro}$ is similar to but different from the spectroscopic mass function
ordinarily defined (hereafter $\hat{f}_{\rm m,spectro}$), since we obtain
$f_{\rm m,spectro}$, dividing $\hat{f}_{\rm m,spectro}$ by $\sin^{3}i$. We can
know the inclination angle, $i$, thanks to astrometric observation, and thus
mainly refer to $f_{\rm m,spectro}$, not $\hat{f}_{\rm m,spectro}$.

Practically, we calculate $f_{\rm m,spectro}$ of AstroSpectroSB1 binary stars
as

$\displaystyle f_{\rm m,spectro}=\left[\left(\frac{C_{1}}{{\rm
au}}\right)^{2}+\left(\frac{H_{1}}{{\rm
au}}\right)^{2}\right]^{3/2}\left(\frac{P}{{\rm
yr}}\right)^{-2}\sin^{-3}i\;[M_{\odot}],$ (5)

where $C_{1}$ and $H_{1}$ are Thiele-Innes elements (Binnendijk, 1960; Heintz,
1978), derived by spectroscopic observation. On the other hand, we calculate
$f_{\rm m,spectro}$ of Orbital binary stars, substituting half
rv_amplitude_robust into $K_{1}$.

We regard binary stars as BH binary candidates if they satisfy the following
two conditions:

$\displaystyle 0.5\leq f_{\rm m,spectro}/f_{\rm m,astro}\leq 2,$ (6)
$\displaystyle f_{\rm m,astro}\geq 3M_{\odot}.$ (7)

We adopt the first condition expressed by Eq. (6) for the following reason.
When a binary star is a BH binary, the secondary star is an unseen star;
$F_{2}/F_{1}=0$. Substituting $F_{2}/F_{1}=0$ into Eq. (1), we find $f_{\rm
m,astro}=f_{\rm m,spectro}$. Thus, BH binaries should satisfy $f_{\rm
m,astro}\simeq f_{\rm m,spectro}$. By the second condition of Eq. (7), we can
select binary star candidates with $m_{2}\geq 3M_{\odot}$ irrespective of
$m_{1}$. Such binary stars are likely to be BH binaries, since the maximum
mass of neutron stars is expected to be $\sim 2M_{\odot}$ (Kalogera & Baym,
1996).

Figure 1: Scatter plots of $f_{\rm m,astro}$ and $f_{\rm m,spectro}/f_{\rm
m,astro}$ for AstroSpectroSB1 (top) and Orbital (bottom) binary stars. The
color scale represents the square root of the relative density of binary
stars. Shaded regions satisfy the two conditions of BH binary candidates
expressed as Eqs. (6) and (7). The BH binary candidate found in this work
(GDR3 5870569352746779008) is emphasized as a star in the top panel.

Figure 1 shows $f_{\rm m,astro}$ and $f_{\rm m,spectro}/f_{\rm m,astro}$ of
all the samples. The shaded region in this figure corresponds to the two
conditions imposed in this study (Eqs. (6) and (7)). Only one binary star
satisfies these two conditions. Its basic parameters are summarized in Table
1. We analyze this BH binary candidate in later sections.

In general, we have $f_{\rm m,spectro}\geq f_{\rm m,astro}$ for any binary
stars, which can be easily confirmed from their definitions in Eq. (1) and Eq.
(3). However, Figure 1 shows that the distribution of $f_{\rm
m,spectro}/f_{\rm m,astro}$ spreads under 1. There can be two reasons. First,
$f_{\rm m,spectro}$ is underestimated for the second type of binary stars. For
these binary stars, we adopt rv_amplitude_robust for $K_{1}$ in Eq. (4).
However, the observed radial velocities may not fall at the right phase to
fully sample the orbit’s maximum and minimum radial velocities. Second, some
of binary stars contain large errors of either $f_{\rm m,spectro}$ or $f_{\rm
m,astro}$, while they have $f_{\rm m,spectro}\geq f_{\rm m,astro}$ in reality.
In fact, such binary stars may hide BH binaries. However, in this paper, we
conservatively select binary stars with $f_{\rm m,spectro}\simeq f_{\rm
m,astro}$ as BH binary candidates. This is because the small discrepancy
between $f_{\rm m,spectro}$ and $f_{\rm m,astro}$ is anticipated for a binary
system in which the secondary star is much fainter than the primary star
($F_{2}/F_{1}\simeq 0$).

It is a bit strange that the $\log(f_{\rm m,spectro}/f_{\rm m,astro})$ values
are centered on zero for both the AstroSpectroSB1 and Orbital binary stars.
Typically, binary stars should have luminous secondary stars (e.g. Sana et
al., 2012), and thus should have $F_{2}/F_{1}$ close to 1, and large $f_{\rm
m,spectro}/f_{\rm m,astro}$ (or small $f_{\rm m,astro}$). The reason for this
discrepancy might be that Gaia preferentially select binary stars with faint
secondary stars.

Figure 2: Distributions of the $\log f_{\rm m,astro}$ and $\log f_{\rm
m,spectro}$ dispersions for AstroSpectroSB1 binary stars. The red and blue
dashed lines indicate the $\log f_{\rm m,astro}$ and $\log f_{\rm m,spectro}$
dispersions, respectively, for the BH binary candidate (GDR3 source ID
5870569352746779008).

Figure 2 shows the distributions of the $\log f_{\rm m,astro}$ and $\log
f_{\rm m,spectro}$ dispersions for AstroSpectroSB1 binary stars. In order to
obtain these dispersions, we generate $10^{3}$ Monte Carlo random draws of the
covariance matrix of the BH binary candidate in the GDR3 nss_two_body_orbit
table. Note that the number of Monte Carlo random draws is sufficient, since
the distributions for $10^{3}$ random draws are similar to those for $10^{2}$
random draws. We cannot calculate the $\log f_{\rm m,spectro}$ dispersions of
Orbital binary stars, since the covariance matrix does not include the data of
rv_amplitude_robust. The $\log f_{\rm m,spectro}$ dispersions of Orbital
binary stars should be larger than that of AstroSpectroSB1 binary stars. This
is also the reason why $\log(f_{\rm m,spectro}/f_{\rm m,astro})$ values spread
more widely in Orbital binary stars than in AstroSpectroSB1 binary stars. The
peak of the $\log f_{\rm m,spectro}$ dispersion is at $\sim 0.2$. On the other
hand, $\gtrsim 100$ binary stars have the $\log f_{\rm m,spectro}$ dispersions
of $\gtrsim 1$. This should affect the presence of binary stars with $\log
f_{\rm m,spectro}/f_{\rm m,astro}<0$.

If a binary star system is actually a triple star system, $f_{\rm m,spectro}$
will be significantly overestimated, and consequently $f_{\rm
m,spectro}>f_{\rm m,astro}$ for the following reason. Astrometric observation
is sensitive to the outer binary’s motion: the relative motion between the
inner binary and third star. On the other hand, spectroscopic observation is
sensitive to the inner binary’s motion, since its motion velocity is much
larger than the outer binary’s motion velocity. Thus, we will calculate
$f_{\rm m,spectro}$ in Eq. (4), using $P$, $e$, and $i$ of the outer binary
and $K_{1}$ of the inner binary. This $f_{\rm m,spectro}$ will be larger than
the actual $f_{\rm m,spectro}$ of both the inner and outer binaries. For
obtaining the inner (outer) binary’s $f_{\rm m,spectro}$, the adopted $P$
($K_{1}$) is larger than the actual inner (outer) binary’s. This may happen in
the Orbital binary stars more frequently.

###  2.2 Some comments on rejected binaries

Before analyzing the BH binary candidate in detail, we review our search. In
particular, we focus on binary stars which look like BH binaries at a glance,
but which our search rejects. GDR3 provides the binary_masses table including
the masses of primary and secondary stars estimated from the PARSEC isochrone
models111http://stev.oapd.inaf.it/cgi-bin/cmd (Bressan et al., 2012). We can
obtain such binary stars with following ADQL query:

    select nss.*, gs.*, bm.*
    from gaiadr3.nss_two_body_orbit as nss,
    gaiadr3.gaia_source as gs,
    gaiadr3.binary_masses as bm
    where gs.source_id = nss.source_id
    and bm.source_id = nss.source_id
    and (nss.nss_solution_type = ’AstroSpectroSB1’
    Ψ or (nss.nss_solution_type = ’Orbital’
             and gs.rv_amplitude_robust IS NOT NULL)
    Ψ or (gs.non_single_star = 3))

We just add the binary_masses table to the ADQL query in section 2.1. Note
that our samples are the ones obtained with the ADQL query in section 2.1
unless otherwise stated. Not all our samples are listed in the binary_masses
table, because the mass estimation is only applied to primary stars in the
main sequence (MS) on the color-magnitude diagram. In the binary_masses table,
there are 6 AstroSpectroSB1 and 3 Orbital binary stars containing secondary
stars with $>3M_{\odot}$. In spite of their secondary masses, none of them are
regarded as BH binary candidates by our search.

As for the 6 AstroSpectroSB1 binary stars, they are rejected, because all of
them have too large $f_{\rm m,spectro}/f_{\rm m,astro}$ ($>10$). This means
that, although these binary stars have main-sequence primary stars with
$1$–$2$ $M_{\odot}$, they have secondary stars with $>3$ $M_{\odot}$ and
smaller (but non-zero) luminosity than the primary stars. It is difficult to
interpret these binary stars as BH binaries. Thus, we remove them from our
list of BH binary candidates. The 3 Orbital binary stars are ruled out, since
they have too small $f_{\rm m,spectro}/f_{\rm m,astro}$ ($<0.01$).
Incomprehensibly, their $F_{2}/F_{1}$ values are negative. Astrometric or
spectroscopic results might not be appropriate. In fact, all of them have
large goodness-of-fit values ($>5$), where the goodness-of-fit is expected to
obey the normal distribution if astrometric parameters are correctly derived.
When Andrews et al. (2022) search for NS and BH binaries, they rule out binary
stars with goodness-of-fit values more than 5 from NS and BH binary
candidates.

The second condition expressed by Eq. (7) may be too strict to complete a
search for BH binaries from our sample. This condition means that the
secondary mass is more than $3$ $M_{\odot}$ for any primary masses. We convert
this condition to $m_{2}>3$ $M_{\odot}$, where $m_{2}$ is drawn from the lower
limit of $m_{2}$ (m2_lower) in the GDR3 binary_masses table. By this
conversion, we can relax our search for BH binaries, since the secondary mass
can be more than $3$ $M_{\odot}$ even for $f_{\rm m,astro}<3$ $M_{\odot}$ if
the primary mass is larger than a certain value. However, we find no other BH
binary candidate. Although the two conditions expressed by Eqs. (6) and (7)
are slightly strict, we confirm that there is only one BH binary candidate
(GDR3 source ID 5870569352746779008) in GDR3 astrometric binary stars with
spectroscopic data.

##  3 Analysis of a BH candidate

We summarize the basic parameters of the BH binary candidate in Table 1. For
the right ascension, declination, BP-RP color, reddening of BP-RP color,
[M/H], and surface gravity ($\log g$), we adopt the mean values in the GDR3
gaia_source table. The galactic longitude and latitude are derived from the
right ascension and declination. We obtain the mean value of the extinction in
G band ($A_{\rm G}$) from the EXPLORE G-Tomo scientific data application
(Lallement et al., 2022; Vergely et al., 2022)222https://explore-platform.eu/,
while the value in the parentheses is the mean value of the GDR3 gaia_source
table. Hereafter, we adopt the former value for the extinction. We obtain the
goodness-of-fit value from the GDR3 nss_two_body_orbit table. In order to
calculate the mean values and one standard deviation intervals of the
distance, period ($P$), physical semi-major axis ($a_{1}/\varpi$),
eccentricity ($e$), inclination ($i$), radial velocity semi-amplitude
($K_{1}$), astrometric mass function ($f_{\rm m,astro}$), and spectroscopic
mass function ($f_{\rm m,spectro}$), we generate $10^{4}$ Monte Carlo random
draws of the covariance matrix of the BH binary candidate in the GDR3
nss_two_body_orbit table333The number of Monte Carlo random draws is
sufficiently large, since the results are similar if we adopt $10^{3}$ for the
number of the random draws.. In this method, we also obtain $f_{\rm
m,astro}>5.68$ $M_{\odot}$ and $f_{\rm m,spectro}>6.57$ $M_{\odot}$ at a
probability of $99$ %. Note that the distance is calculated from the parallax
in the GDR3 nss_two_body_orbit table, not in the GDR3 gaia_source table.
According to Gaia Collaboration et al. (2022b), the parallax in the former
table is more accurate than in the latter table. We get the absolution
magnitude in G band ($M_{\rm G}$) from the mean of apparent magnitude in the
GDR3 gaia_source, and the mean of the distance derived above.

Table 1: Basic parameters of the BH binary candidate.

Quantities | Values
---|---
(1) Source ID | 5870569352746779008
(2) Orbital solution | AstroSpectroSB1
(3) Right ascension | $207.5697^{\circ}$
(4) Declination | $-59.2390^{\circ}$
(5) Galactic longitude | $310.4031^{\circ}$
(6) Galactic latitude | $2.7765^{\circ}$
(7) Absolute magnitude in G band ($M_{\rm G}$) |  $1.95$ mag
(8) Extinction in G band ($A_{\rm G}$) |  $0.5628$ mag ($0.70$ mag)
(9) BP-RP color |  $1.49$ mag
(10) Reddening of BP-RP color |  $0.37$ mag
(11) Surface gravity ($\log g$, logg_gspphot) |  $3.25$ [cgs]
(12) $[{\rm M/H}]$ (mh_gspphot) |  $0.0066$ dex
(13) Goodness-of-Fit | $3.07$
(14) Distance |  $1164.41\pm 25.16$ pc
(15) Period ($P$) |  $1352.22\pm 45.81$ day
(16) Physical semi-major axis ($a_{1}/\varpi$) |  $4.5194\pm 0.1305$ au
(17) Eccentricity ($e$) | $0.5323\pm 0.0153$
(18) Inclination ($i$) | $35.15\pm 0.99^{\circ}$
(19) Radial velocity semi-amplitude ($K_{1}$) |  $27.0\pm 1.0$ km s-1
(20) Astrometric mass function ($f_{\rm m,astro}$) |  $6.75\pm 0.51$ $M_{\odot}$
(21) Lower bound in $f_{\rm m,astro}$ (99%) | $f_{\rm m,astro}>5.68M_{\odot}$
(22) Spectroscopic mass function ($f_{\rm m,spectro}$) |  $8.85\pm 1.13$ $M_{\odot}$
(23) Lower bound in $f_{\rm m,spectro}$ (99%) | $f_{\rm m,spectro}>6.57M_{\odot}$
(24) Probability of $f_{\rm m,spectro}>f_{\rm m,astro}$ |  $95$ %
(25) Probability of $f_{\rm m,spectro}<f_{\rm m,astro}$ |  $5$ %

Note. — From row 3 (right ascension) to 13 (goodness-of-fit) except for row 8
(Extinction in G band), we show the mean value in GDR3. For row 8 (Extinction
in G band), the value is obtained from the EXPLORE G-Tomo scientific data
application, and the value in the parentheses is obtained from the mean value
in GDR3. From row 14 (distance) to 20 ($f_{\rm m,astro}$) as well as in row 22
($f_{\rm m,spectro}$), we show the mean value and one standard deviation
interval. In rows 21 and 23, we show the 99% confidence level of $f_{\rm
m,astro}$ and $f_{\rm m,spectro}$, respectively. In rows 24 and 25, we show
the probabilities of $f_{\rm m,spectro}>f_{\rm m,astro}$ and $f_{\rm
m,spectro}<f_{\rm m,astro}$, respectively (see section 3 for more detail).

The goodness-of-fit value, $3.07$, is relatively low, since Andrews et al.
(2022) consider that NS and BH binary candidates should have the goodness-of-
fit value less than $5$. Note that the goodness-of-fit value for reliable
sources should be normally distributed with a mean of zero. Thus, we are not
going to arguing that the BH binary candidate is reliable only from the
goodness-of-fit value. Nevertheless, we have to remark that the goodness-of-
fit value is typical of the AstroSpectroSB1 binary stars as described later
(see Figure 7). Although the goodness-of-fit value largely deviates from zero,
it would not directly mean that this BH binary candidate is unreliable. We
find that the ratios of means to standard deviation intervals are high for
$f_{\rm m,astro}$ and $f_{\rm m,spectro}$ ($13.2$ and $7.83$, respectively).
They should be relatively well-measured. Additionally, the $\log f_{\rm
m,astro}$ and $\log f_{\rm m,spectro}$ dispersions are small, compared to
those of other AstroSpectroSB1 binary stars, as seen in Figure 2. This should
be another evidence that the parameters of this BH candidate are well-
measured. Moreover, at a probability of $99$ %, $f_{\rm m,astro}>5.68$
$M_{\odot}$ and $f_{\rm m,spectro}>6.75$ $M_{\odot}$. These values are
unlikely to fall below $3$ $M_{\odot}$. A concern is that $f_{\rm m,spectro}$
is systematically larger than $f_{\rm m,astro}$, which we discuss in section
4.

Figure 3: Color magnitude diagram of stars in the GDR3 gaia_source table.
These stars are filtered in the same way as those in figure 6c of Gaia
Collaboration et al. (2018b). The color scale represents the square root of
the relative density of stars. The star and diamond points indicate the BH
binary candidate (GDR3 source ID 5870569352746779008), where the star and
diamond points are not corrected and corrected by its extinction and
reddening, respectively. We define the regions of MS and RGB stars below and
above the dashed line, respectively. The line is expressed as Eq. (10).

Figure 3 shows the color-magnitude diagram of the primary star of the BH
binary candidate, and GDR3 stars whose G-band absolute magnitudes and BP-RP
colors are well-measured. MS and red giant branch (RGB) regions are defined as
regions below and above the dashed line. The dashed line is expressed as

$\displaystyle M_{\rm G}=\left\\{\begin{array}[]{ll}3.14\left({\rm BP}-{\rm
RP}\right)-0.43&({\rm BP}-{\rm RP}<1.41)\\\ 4&({\rm
otherwise})\end{array}\right..$ (10)

The first case of Eq. (10) is the same as in Andrews (2022). We induce the
second case to avoid regarding low-mass MS stars as RGB stars. As seen in
Figure 3, the primary star of the BH binary candidate is likely to be a RGB
star. This is consistent with its small surface gravity ($\log g=3.25$). The
primary star indicated by the star point (not corrected by its extinction and
reddening) is redder than RGB stars on the color-magnitude diagram. It suffers
from interstellar reddening, since it is located in the Galactic disk
($b=2.7765^{\circ}$). In fact, the primary star indicated by the diamond point
(corrected by its extinction and reddening) is on the RGB.

Generally, BH binary candidates are thought dubious when their primary stars
are RGB stars. This is because such primary stars can easily outshine
companion stars even if the companion stars are more massive than the primary
stars. Moreover, it is difficult to estimate the masses of RGB stars in binary
systems. Such RGB stars can be in so-called Algol-type systems (El-Badry et
al., 2022b). They can be luminous but low-mass (say $\sim 0.1$ $M_{\odot}$) if
they experience mass transfer. These types of problems frequently happen in BH
binary candidates with only spectroscopic data, or usual spectroscopic mass
function, $\hat{f}_{\rm m,spectro}$ (not $f_{\rm m,spectro}$), and
$\hat{f}_{\rm m,spectro}$ is $\sim 1$ $M_{\odot}$. In order to conclude that
their secondary stars are $>3$ $M_{\odot}$ compact objects, we need to
estimate the primary stars’ masses and inclination angles of the binary stars.
As an illustration, let us consider a spectroscopic binary characterized by
$\hat{f}_{\rm m,spectro}=1$ $M_{\odot}$ and inclination angle $i=60^{\circ}$
(from which we obtain $f_{\rm m,spectro}=1.54$ $M_{\odot}$). If the the
primary star’s mass is $1.2$ $M_{\odot}$, the secondary star’s mass is $3$
$M_{\odot}$. In this case, the $3$ $M_{\odot}$ secondary star is highly likely
a BH. In contrast, if the the primary star’s mass is $0.2$ $M_{\odot}$, the
secondary star’s mass is $1.9$ $M_{\odot}$. In this case, we cannot exclude
the possibility that the $1.9$ $M_{\odot}$ secondary star is a main-sequence
star outshined by the primary RGB star. From this simple illustrative example,
we can see that the mass estimation of primary stars critically affects
whether their secondary stars are BH or not.

Fortunately, these types of problems do not happen in our BH binary candidate.
We know the inclination angle $i$ of the binary star from the astrometric
data, and get $f_{\rm m,spectro}$ in a model-independent way. Moreover, this
BH binary candidate has $f_{\rm m,astro}>5.68$ $M_{\odot}$ and $f_{\rm
m,spectro}>6.75$ $M_{\odot}$ at a probability of 99 %. The secondary mass is
more than $5$ $M_{\odot}$, even if this BH binary candidate is an Algol-type
system, or the primary RGB mass is close to zero. The primary RGB star cannot
outshine the $>5$ $M_{\odot}$ secondary star even if the secondary star is in
the main-sequence phase, or the faintest among $5$ $M_{\odot}$ stars in any
phases except a BH. This point is described in detail below. Thus, the
secondary star is likely to be a BH.

We examine the possibility that the secondary star of the BH binary candidate
may be a single object except a BH, or multiple star systems. When a stellar
mass is fixed, a MS star is the faintest except stellar remnants like WD, NS,
and BH. If a MS star with the same mass as the secondary star is more luminous
than the primary star, the possibility that the secondary star is a single
object except a BH can be ruled out. When the total mass of a multiple star
system is fixed, a multiple star system with equal-mass MS stars is the least
luminous. This is because MS stars become luminous more steeply with their
masses increasing. If an $n$-tuple star system with equal-mass MS stars has
the same mass as the secondary star, and larger luminosity than the primary
star, the possibility that the secondary star is any $n$-tuple star systems
can be rejected. Thus, we compare the luminosity of the primary star with the
luminosities of a single MS star or multiple MS star systems with equal
masses.

Figure 4: G-band absolute magnitude of multiple star systems with equal-mass
MS stars whose ages are $10^{7}$, $10^{8}$, $10^{9}$, and $10^{10}$ yrs. The
total mass of the multiple star systems is $5.68$ $M_{\odot}$, the lower bound
mass of the secondary star of the BH binary candidate at a probability of 99
%. The component mass and the number of stars are shown in the lower and upper
$x$-axis, respectively. We show only MS stars defined in Eq. (10). That is the
reason why the curves of $10^{9}$ and $10^{10}$ yrs cut off in the middle. We
obtain G-band absolute magnitude and BP-RP color at each mass and age, using
the PARSEC code (Bressan et al., 2012). The metallicity is set to the solar
metallicity, the same as the primary star of the BH binary candidate. It seems
that there is no publicly available spectroscopic survey data that provides
reliable metallicity for the primary star of our BH binary candidate. The
dashed line indicates the G-band absolute magnitude of the primary star, which
is corrected by the G-band extinction. The dotted line indicates the G-band
absolute magnitude of a star half as luminous as the primary star.

Figure 4 shows the G-band absolute magnitude of multiple star systems with
equal-mass MS stars. The total mass of the multiple star systems is $5.68$
$M_{\odot}$, the lower bound mass of the secondary star of the BH binary
candidate at a probability of 99 %. We can rule out single, binary and triple
stars with the total mass of $5.68$ $M_{\odot}$. They would outshine the
primary star if they were the secondary star. A quadruple star system with
each stellar mass of $1.4$–$1.5$ $M_{\odot}$ is as luminous as the primary
star. However, such a quadruple star system should be detected by Gaia itself.
A quintuple star system with each stellar mass of $1.1$ $M_{\odot}$ has
luminosity twice less than the primary star, and might not be observed by
Gaia. Except for multiple star systems with MS stars, the secondary star can
be a triple NS star system or a quadruple WD star system, where the maximum
mass of NS and WD are about $2.0$ and $1.4$ $M_{\odot}$, respectively. Such
systems may be more valuable than a single BH, since they have never been
discovered to our knowledge. In any case, the secondary star should be
quadruple or higher-order star systems except for a single BH. Moreover, the
size of the system should be more compact than the pericenter distance of the
primary star, $\sim 2.4$ au. It is unclear that such multiple systems are
stable under the perturbation of the primary star.

Figure 5: Bottom left: Scatter plots of $f_{\rm m,astro}$ and $f_{\rm
m,spectro}/f_{\rm m,astro}$ for RGB stars in AstroSpectroSB1. The color scale
represents the square root of the relative density of binary stars. Contours
indicate $\sigma$ levels of $1$, $2$, $\cdots$, and $7$ from the inner to the
outer. The shaded region is considered to calculate the p-values in Table 2.
The p-values are calculated by a kernel-density estimate with the kernel
bandwidth of Scott’s rule (Scott, 1992). The star point indicates the BH
candidate (GDR3 source ID 5870569352746779008). It is not included in the
samples with which the p-values are calculated. Top and left: $f_{\rm
m,astro}$ and $f_{\rm m,spectro}/f_{\rm m,astro}$ distributions, respectively.
The histograms indicate the sample distribution, and the curves indicate the
projected distributions derived by the kernel-density estimate.

In order to assess whether the BH binary candidate is coincidentally located
on the $f_{\rm m,astro}$–$f_{\rm m,astro}/f_{\rm m,spectro}$ plane, we
calculate the p-values of a $f_{\rm m,astro}$–$f_{\rm m,astro}/f_{\rm
m,spectro}$ region around the BH binary candidate. We adopt a kernel-density
estimate with a kernel bandwidth of Scott’s rule (Scott, 1992). The bandwidth
is $N_{\rm sample}^{-1/6}$, where $N_{\rm sample}$ is the number of samples.
At first, we select RGB primary stars from AstroSpectroSB1 as samples for the
kernel-density estimate. The number of samples is $9047$. Note that the BH
binary candidate is excluded from the samples. Figure 5 shows the kernel-
density contours of $1$, $2$, $\cdots$, and $7\sigma$ levels from the inner to
the outer. We calculate the p-value in the shaded region. The p-value is
$9.6\times 10^{-12}$, and the $\sigma$ level is $6.1$. The position of $f_{\rm
m,astro}$ and $f_{\rm m,astro}/f_{\rm m,spectro}$ of the BH binary candidate
is unlikely to be coincident.

Table 2: P-values.  Sample | Number | p-value |  $\sigma$ level | Remark
---|---|---|---|---
All | $64095$ | $2.4\times 10^{-12}$ | $7.0$ |
All in AstroSpectroSB1 | $33466$ | $9.1\times 10^{-12}$ | $6.8$ |
Low-error in AstroSpectroSB1 | $28188$ | $1.0\times 10^{-11}$ | $6.8$ | Exclude samples with top 10 % large errors
Low-error in AstroSpectroSB1 | $17614$ | $1.1\times 10^{-11}$ | $6.8$ | Exclude samples with errors more than 0.2 in log-scale
RGBs in AstroSpectroSB1 | $9047$ | $9.6\times 10^{-12}$ | $6.1$ | The same samples used in Figure 5
Low-error RGBs in AstroSpectroSB1 | $8626$ | $7.5\times 10^{-10}$ | $6.2$ | Exclude samples with top 10 % large errors
Low-error RGBs in AstroSpectroSB1 | $5395$ | $1.2\times 10^{-9}$ | $6.1$ | Exclude samples with errors more than 0.2 in log-scale

We select samples for the kernel-density estimate in different ways in order
to investigate whether the p-values depend on the choice of samples. We
summarize the choices of samples and their results in Table 2. The first
column indicates the choice of samples. Note that the BH binary candidate is
not included in any choices. For “All”, we choose all the samples selected in
section 2. For “All in AstroSpectroSB1”, we choose all the samples in
AstroSpectroSB1. For “RGBs in AstroSpectroSB1”, we extract only the RGB
primary stars in the samples of “All in AstroSpectroSB1”. This is the samples
shown in Figure 5. We also make samples, excluding samples with large errors
of $f_{\rm m,astro}$ and $f_{\rm m,spectro}$ from “All in AstroSpectroSB1” and
“RGBs in AstroSpectroSB1”. We calculate the errors in the same way as the one
standard deviation of the BH binary candidate in Table 1, where we generate
$10^{3}$ Monte Carlo random draws for each sample for calculation cost
savings. We adopt two cases to exclude samples. In the first case, we exclude
10 % samples with largest errors in either of $f_{\rm m,astro}$ and $f_{\rm
m,spectro}$. In the second case, we exclude samples with errors larger than
$0.2$ in log-scale for either of $f_{\rm m,astro}$ and $f_{\rm m,spectro}$.
Note that $0.2$ is similar to the bandwidth of the kernel-density estimate. In
any cases, the p-values are small, and the $\sigma$ levels are high. The
position of $f_{\rm m,astro}$ and $f_{\rm m,astro}/f_{\rm m,spectro}$ of the
BH binary candidate is unlikely to be coincident, independently of the choices
of samples for the kernel-density estimate.

We search for the BH binary candidate in several database. The GDR3
variability table (Eyer et al., 2022) and the All-Sky Automated Survey for
SuperNovae (ASAS-SN; Kochanek et al., 2017) do not include the BH binary
candidate as a variable star. Its light curve is available on the ASAS-SN
Photometry Database (Shappee et al., 2014; Jayasinghe et al.,
2019)777https://asas-sn.osu.edu/photometry. The BH binary candidate is
observed in V and g band over $\sim 3000$ days. We do not find any periodic
feature. The Transiting Exoplanet Survey Satellite (TESS; Ricker et al., 2015)
has performed twice high-cadence observations during about 30 days according
to data downloaded from TESScut444https://mast.stsci.edu/tesscut/ for the BH
binary candidate. The duration is too short to detect its periodic variability
due to its binary orbit if any, since it has a period of about 1000 days.
Wide-field Infrared Survey Explorer (WISE; Wright et al., 2010) also observes
the BH binary candidate over $\sim 4000$ days according to ALLWISE Multiepoch
Photometry Table555https://irsa.ipac.caltech.edu/cgi-bin/Gator/nph-
scan?submit=Select&projshort=WISE and NEOWISE-R Single Exposure (L1b) Source
Table666https://irsa.ipac.caltech.edu/cgi-bin/Gator/nph-dd. We do not
recognize any periodic variability. At the time of September 2022, the BH
binary candidate is not listed in the following data base:
SIMBAD888http://simbad.cds.unistra.fr/simbad/, the ninth catalog of
spectroscopic binary orbits (SB9; Pourbaix et al., 2004), RAdial Velocity
Experiments (RAVE; Kunder et al., 2017), the Galactic Archaeology with HERMES
(GALAH; Buder et al., 2021), the Large sky Area Multi-Object fiber
Spectroscopic Surveys (LAMOST; Cui et al., 2012), and the Apache Point
Observatory Galactic Evolution Experiment (APOGEE; Majewski et al., 2017).
High-energy telescopes, such as the Fermi gamma-ray space telescope (Atwood et
al., 2009), the Swift Burst Alert Telescope (Swift BAT: Barthelmy et al.,
2005) XMM-Newton (Strüder et al., 2001), the Chandra observatory (Weisskopf et
al., 2000), and the Galaxy Evolution Explorer (GALEX; Martin et al., 2005), do
not observe it as far as we see Aladin
lite999https://aladin.u-strasbg.fr/AladinLite/. ESO
archive101010http://archive.eso.org/scienceportal/home does not list it. In
summary, we do not find any positive nor negative evidence for the BH binary
candidate.

##  4 Discussion

First, we compare the BH binary candidate with other BH binary candidates by
previous studies, and assess whether our BH binary candidate is similar to
others rejected before. As described in section 3, BH binary candidates tend
to be rejected when their primary stars are RGB stars. It is difficult to
estimate the masses of RGB stars, and such binary stars can be Algol-type
systems in which primary stars are low-mass (say $\sim 0.1$ $M_{\odot}$).
Since such BH binary candidates have $\hat{f}_{\rm m,spectro}\sim 1$
$M_{\odot}$, the mass estimate of RGB stars severely affect the secondary
mass. However, our BH binary candidate has $f_{\rm m,astro}>5.68$ $M_{\odot}$
and $f_{\rm m,spectro}>6.75$ $M_{\odot}$ at a probability of 99 %. In this
case, the secondary mass is more than $\sim 5$ $M_{\odot}$ even if the primary
mass is nearly zero. Note that the secondary mass increases monotonically with
the primary mass increasing when $f_{\rm m,astro}$ or $f_{\rm m,spectro}$ is
fixed. Thus, the secondary star is likely to be a BH, even if the BH binary
candidate is an Algol-type system.

Gaia Collaboration et al. (2022b) listed up BH binary candidates with $\sim 2$
$M_{\odot}$ MS stars and $\sim 3$ $M_{\odot}$ BHs. However, El-Badry & Rix
(2022) pointed out possibility that they are Algol-type systems consisting of
$\sim 0.2$ $M_{\odot}$ stripped stars and $\sim 2$ $M_{\odot}$ MS stars. The
reason for this discrepancy is as follows. Gaia Collaboration et al. (2022b)
thought that $\sim 2$ $M_{\odot}$ MS stars dominate the luminosity
(photometry) and radial-velocity motion (spectroscopy) of the binary stars. On
the other hand, El-Badry & Rix (2022) claimed that $\sim 2$ $M_{\odot}$ MS
stars dominate the luminosity, while $\sim 0.2$ $M_{\odot}$ stripped stars
dominate the radial-velocity motion. This interpretation better explains their
spectral energy distribution and spectroscopic mass function ($\hat{f}_{\rm
m,spectro}\sim 1.5$ $M_{\odot}$) more naturally. We do not expect that similar
things happen in our BH binary candidate for the following reason. If a hidden
star dominates the radial-velocity motion, we replace $m_{2}$ with $m_{1}$ in
Eq. (3). Since $f_{\rm m,astro}\sim f_{\rm m,spectro}$, we obtain
$m_{1}=4f_{\rm m,astro}(1+F_{2}/F_{1})^{2}$ $M_{\odot}$ and $m_{2}=4f_{\rm
m,astro}(1+F_{2}/F_{1})^{2}(1+2F_{2}/F_{1})$ $M_{\odot}$. Thus, the RGB
primary mass should be at least $4f_{\rm m,astro}$ ($\sim 23$) $M_{\odot}$.
However, its luminosity requires its mass much less than $23$ $M_{\odot}$.
Thus, a hidden star does not dominate the radial-velocity motion of our BH
binary candidate in contrary to the BH binary candidates in table 10 of Gaia
Collaboration et al. (2022b).

Gaia Collaboration et al. (2022b) also show another table of BH binary
candidates (their table 9) in which BH binary candidates belong to SB1, and
have high $\hat{f}_{\rm m,spectro}$ ($>3$ $M_{\odot}$). Hereafter, we call
them “Gaia’s table 9 candidates”. Although these candidates have secondary
stars with more than $3$ $M_{\odot}$ for any primary masses, Gaia
Collaboration et al. (2022b) cannot rule out that the secondary stars consist
of multiple star systems, similarly to our description in section 3. We remark
that our BH candidate will be better-constrained than all of Gaia’s table 9
candidates. Our BH binary candidate has larger mass function and smaller
luminosity than Gaia’s table 9 candidates except for GDR3 source IDs
4661290764764683776 and 5863544023161862144. GDR3 source ID
4661290764764683776 has high $\hat{f}_{\rm m,spectro}$ ($=13.67$ $M_{\odot}$),
however its primary star has high luminosity, $-6.707$ mag in G band. Since
the primary star can be more luminous than a $\sim 13$ $M_{\odot}$ MS star, it
is difficult to confirm that the secondary star is a BH. GDR3 source ID
5863544023161862144 shows eclipses, and consequently its secondary should not
be a BH (Gaia Collaboration et al., 2022b). In summary, we can easiest rule
out the possibility that the secondary star of our BH candidate consists of a
multiple star system.

Figure 6: Comparison of orbital parameters (Period, eccentricity, radial
velocity semi-amplitude, and $\hat{f}_{\rm m,spectro}$) between GDR3 and SB9.
Triangle, circle and star points indicate GDR3 SB1, AstroSpectroSB1 with a
period of $\leq 10^{3}$ day, and AstroSpectroSB1 with a period of $>10^{3}$
day, respectively. The dashed lines show the mean values of the orbital
parameters of the BH binary candidate (GDR3 source ID 5870569352746779008).

Pourbaix et al. (2022) and Jayasinghe et al. (2022a) compared orbital
parameters in GDR3 with those in SB9 (Pourbaix et al., 2004), in particular
for spectroscopic binary stars with either one component being parameterized
(SB1). They found that Gaia’s and SB9’s periods are inconsistent for periods
of more than $10^{3}$ days in SB9. Since they did not investigate
AstroSpectroSB1, we investigate both of SB1 and AstroSpectroSB1. We find 304
SB1 and 109 AstroSpectroSB1 in common between GDR3 and SB9. Our BH binary
candidate is not included in SB9 as described in the previous section. In
Figure 6, we make comparison between orbital parameters of binary stars in
GDR3 and SB9. Note that the $x$-axes in Figure 6 adopt GDR3 values, while
Pourbaix et al. (2022) and Jayasinghe et al. (2022a) adopt SB9 values for the
$x$-axes in their figure 7.41 and figure 6, respectively. Similarly to
Pourbaix et al. (2022) and Jayasinghe et al. (2022a), we find that periods in
GDR3 are largely different from those in SB9 for SB1 with periods of more than
$10^{3}$ days. However, for AstroSpectroSB1, their periods do not deviate up
to periods of a few $10^{3}$ days. The other parameters in GDR3 are also in
good agreement with those in SB9 for AstroSpectroSB1, in particular around the
mean values of the orbital parameters of the BH binary candidate. This does
not directly show that the spectroscopic data of the BH binary candidate is
reliable, since most of binary stars in SB9 are brighter than our BH binary
candidate. Nevertheless, this means that GDR3 values of AstroSpectroSB1 binary
stars may be reliable even if the binary stars have periods of a few $10^{3}$
days.

Bashi et al. (2022) compared GDR3 SB1 with the database of LAMOST (Cui et al.,
2012) and GALAH (Buder et al., 2021), and found that GDR3 SB1 with periods of
less than $10^{1.5}$ days may be refuted. Although our BH binary candidate
belongs to AstroSpectroSB1 (not to SB1), it has a period of $\gtrsim 10^{3}$
days, much larger than $10^{1.5}$ days. Our BH binary candidate may not be
refuted by the criteria of Bashi et al. (2022).

Andrews et al. (2022) and Shahaf et al. (2023) independently presented lists
of NS and BH binary candidates in GDR3. Their lists do not include our BH
binary candidate. This is because they focus on binary stars with primary MS
stars. The masses of MS stars can be estimated less model-dependently than
those of RGB stars. The masses and natures of secondary stars can be derived
robustly. Thus, they avoided binary stars with primary RGB stars. On the other
hand, although the primary star of our BH binary candidate is a RGB star, we
can call it a “BH binary candidate”, because its $f_{\rm m,astro}$ and $f_{\rm
m,spectro}$ are high; $\gtrsim 5.68$ and $6.75$ $M_{\odot}$, respectively, at
a probability of 99 %. Its secondary mass is more than $5$ $M_{\odot}$,
regardless of the primary RGB mass.

Figure 7: The ratio of $f_{\rm m,spectro}$ to $f_{\rm m,astro}$ as a function
of goodness-of-fit values for AstroSpectroSB1 (top) and Orbital (bottom)
binary stars. The color scale represents the square root of the relative
density of binary stars. The star point indicates the BH binary candidate
(GDR3 source ID 5870569352746779008).

Conversely, we examine the lists of Andrews et al. (2022) and Shahaf et al.
(2023) from our conditions. We focus on binary stars with $m_{2}>2$
$M_{\odot}$ in their lists. Note that the maximum NS mass can be $\sim 2$
$M_{\odot}$. Our sample selected in section 2 does not include BH binary
candidate in the Andrews’s list. The candidates do not have spectroscopic
data. This may be partly because astrometric binary stars with spectroscopic
data (i.e. our sample) has systematically large goodness-of-fit values. Figure
7 shows that goodness-of-fit values in AstroSpectroSB1 and Orbital are
centered at $\sim 3$ and $\sim 5$, respectively. Actually, this can be seen in
the middle panel of figure 4 in Andrews et al. (2022). Their figure 4 includes
all the Orbital binary stars with and without spectroscopic data, and
indicates the second peak around the goodness-of-fit value of $\sim 5$. The
second peak should consist of Orbital binary stars with spectroscopic data. We
do not know the reason for the systematic upward shift. We have to remark that
bright binary stars (G band magnitude of $<13$), i.e. those with spectroscopic
data, have the systematically higher goodness-of-fit values, while faint
binary stars (G band magnitude of $>13$) typically have lower goodness-of-fit
values111111The exact value of $13$ is obtained by El-Badry et al. (2023b)..
In any case, our sample does not include the list of Andrews et al. (2022),
because they avoid including binary stars with goodness-of-fit values of more
than 5 in their list.

Our sample includes Shanaf’s three BH binary candidates (GDR3 source IDs:
3263804373319076480, 3509370326763016704, and 6281177228434199296). However,
we do not list up them as BH binary candidates. This is because their $f_{\rm
m,spectro}/f_{\rm m,astro}$ are small ($0.25$, $0.0053$, $0.0017$,
respectively) for our first condition as seen in Eq. (6). We do not intend to
reject the three BH candidates completely. The three BH candidates may suffer
from large errors of spectroscopic data, and consequently have small $f_{\rm
m,spectro}/f_{\rm m,astro}$. We suspect this possibility, because two of the
three BH candidates are not included in AstroSpectroSB1 binary stars despite
the fact that they have spectroscopic data. Our sample selected in section 2
does not include the other 5 BH binary candidates because of the absence of
spectroscopic data.

A few days after we posted this work on arXiv, El-Badry et al. (2023a)
reported one promising BH binary candidate different from our BH binary
candidate. They made follow-up spectroscopic observation, and showed that
Gaia’s astrometric data is consistent with their spectroscopic data. Since
their BH binary candidate has a shorter period ($185.6$ days) than our binary
candidate ($1352.22$ days), they can finish their follow-up observation in a
short period of time. They also mentioned our BH binary candidate, and did not
conclude whether our BH binary candidate is genuine because of the absence of
follow-up spectroscopic observations. Their argument is in good agreement with
ours. Note that we analyze our BH binary candidate in detail.

Several BH binary candidates can be rejected for exceptional reasons. Although
Gaia Collaboration et al. (2022b) found that GDR3 source ID
2006840790676091776 has high $\hat{f}_{\rm m,spectro}$, they did not include
it in their list of BH binary candidates. This is because it is close to a
bright star, whose apparent magnitude is 3.86 mag in G band. There is no such
bright stars close to our BH binary candidate. Nearby stars have apparent
magnitude of at least 13 mag in G band. The reason for this rejection can not
be applied to our BH binary candidate. Andrews et al. (2022) removed GDR3
source ID 4373465352415301632121212This object is later confirmed as a BH
binary (also known as Gaia BH1) by El-Badry et al. (2023a)., since its period
($\sim 186$ days) is roughly 3 times Gaia’s scanning period ($63$ days). Our
BH binary candidate has a period of $1352$ days, not integer multiple of
Gaia’s scanning period.

Figure 8: Bottom left: 2-dimensional probability distribution of $f_{\rm
m,astro}$ and $f_{\rm m,spectro}$ for the BH binary candidate (GDR3 source ID
5870569352746779008). The cross point means the typical values of $f_{\rm
m,astro}$ and $f_{\rm m,spectro}$. Contours indicate $\sigma$ levels of $1$,
$2$, and $3$ from the inner to the outer. The dashed line shows $f_{\rm
m,astro}=f_{\rm m,spectro}$. Top and left: $f_{\rm m,astro}$ and $f_{\rm
m,spectro}$ distributions, respectively.

Hereafter, we make several concerns. First of all, we mostly rely on GDR3
astrometric and spectroscopic data, which are already largely processed. We do
not assess correctness of the data of our BH binary candidate. Aside from
this, we find that the BH binary candidate has $f_{\rm m,spectro}>f_{\rm
m,astro}$ and $f_{\rm m,spectro}<f_{\rm m,astro}$ at probabilities of $95$ and
$5$ %, respectively (see Table 1). Although $f_{\rm m,spectro}=f_{\rm
m,astro}$ is possible, $f_{\rm m,spectro}$ is always larger than $f_{\rm
m,astro}$ in the $1\sigma$-level region as seen in Figure 8. For comparison,
we calculate the probabilities of $f_{\rm m,spectro}>f_{\rm m,astro}$ and
$f_{\rm m,spectro}<f_{\rm m,astro}$ for GDR3 source ID 5136025521527939072
which is in AstroSpectroSB1, and suggested as a NS binary candidate by Gaia
Collaboration et al. (2022b)131313The NS binary candidate is more likely to be
a WD binary according to El-Badry et al. (2023a).. They are $72$ and $28$ %.
Both of $f_{\rm m,spectro}>f_{\rm m,astro}$ and $f_{\rm m,spectro}<f_{\rm
m,astro}$ are in the $1\sigma$-level region, in contrast to our BH binary
candidate. The $f_{\rm m,spectro}$ and $f_{\rm m,astro}$ of our BH binary
candidate are not as similar as those of GDR3 source ID 5136025521527939072.
Nevertheless, we may regard $f_{\rm m,spectro}=f_{\rm m,astro}$, since our BH
binary candidate may contain some systematic error in either of spectroscopic
or astrometric data.

Another concern is that the primary star of the BH binary candidate is a RGB
star. Theoretical studies (e.g. Shikauchi et al., 2020, 2022) expected that a
BH binary with a $\gtrsim 10$ $M_{\odot}$ MS primary star is likely to be
found first (but see also Shikauchi et al., 2023). This is because such MS
stars are bright, and can be observed even if they are distant. Moreover, they
are longer-lived than RGB stars with similar masses. However, GDR3 does not
present orbital parameters of binary stars with $\gtrsim 10$ $M_{\odot}$ MS
primary stars in AstrospectroSB1 nor Orbital according to the GDR3
binary_masses table obtained with the ADQL query in section 2.2. We do not
know the reason for the absence of such binary stars in GDR3. Nevertheless,
when there are no such binary stars, it may be natural that a BH binary with a
RGB star is first discovered.

We need two types of follow-up observations in order to assess if the BH
binary candidate is true or not. The first type should be spectroscopic
observations to verify the GDR3 spectroscopic data, and to perform spectral
disentangling of the BH binary candidate similar to El-Badry & Rix (2022). The
second type should be deep photometric observations. Such observations could
constrain whether the secondary star is a BH, or consists of multiple stars.
We remark that El-Badry et al. (2023b) have carried out these follow-up
observations, and confirmed it as a genuine BH binary, Gaia BH2. This
demonstrates that these follow-up observations would be important for
confirming or refuting future BH candidates, which may be discovered by our
search methodology in upcoming Gaia data.

##  5 Summary

We first search for BH binary candidates from astrometric binary stars with
spectroscopic data in GDR3. From the sample of 64108 binary stars, we find one
BH binary candidate. The GDR3 source ID is 5870569352746779008\. Since its
primary star is a RGB star, we cannot estimate the mass of the primary RGB
star. However, because of its high astrometric and spectroscopic mass function
($f_{\rm m,astro}>5.68$ $M_{\odot}$ and $f_{\rm m,spectro}>6.75$ $M_{\odot}$
at a probability of 99 %), the secondary star should have more than $5$
$M_{\odot}$, and is likely to be a BH, regardless of the primary mass. Unless
the secondary star is a BH, it must be quadruple or higher-order multiple star
systems with the total mass of $5.68$ $M_{\odot}$. To rule out the possibility
of multiple star systems, we need deep photometric observations. Rather, if it
is quadruple or higher-order multiple star systems, long-term observation may
find modulation of the primary’s orbit (e.g. Hayashi et al., 2020; Hayashi &
Suto, 2020; Liu et al., 2022).

Figure 9: Predicted radial velocities of the BH binary candidate (GDR3 source
ID 5136025521527939072), based on GDR3. The black curve can be obtained from
typical values in GDR3. We generate cyan curves from $10^{2}$ Monte Carlo
random draws, using the covariance matrix. The radial velocity amplitude is
larger and the turnover time spreads more widely than those in the upper panel
of figure 3 in (El-Badry et al., 2023b), because we do not adopt any explicit
link between the spectroscopic and astrometric Thiele-Innes parameters.

The weakness of this paper is that our conclusion entirely relies on the Gaia
DR3. In particular, our BH binary candidate has a period of $\sim 1300$ days,
more than a period of 34 months of the Gaia DR3 data collection. Our
conclusion have to be confirmed by follow-up observations. For example, we
need the time evolution of the radial velocity of our BH binary candidate
similarly to that of Gaia BH1 obtained by El-Badry et al. (2023a). Figure 9
shows the predicted radial velocities of the BH binary candidate. Eventually,
El-Badry et al. (2023b) have confirmed the radial velocity variability, by
observing it $\sim 30$ times from the last half of 2022 to the beginning of
2023, when the radial velocities had steeply decreased.

Previously, RGB stars harboring BHs have not been searched for because of the
difficulty of the mass estimation of the RGB stars (and thus BHs). However,
our tentative discovery in this paper encourages us to explore not only BHs
orbiting around MS stars but also BHs orbiting around RGB stars in the future
data releases of Gaia.

## Acknowledgments

We thank the anonymous referee for many fruitful advices. This research could
not been accomplished without the support by Grants-in-Aid for Scientific
Research (17H06360, 19K03907) from the Japan Society for the Promotion of
Science. KH is supported by JSPS KAKENHI Grant Numbers JP21K13965 and
JP21H00053. NK is supported by JSPS KAKENHI Grant Numbers JP22K03686. TK is
supported by JSPD KAKENHI Grant Number JP21K13915, and JP22K03630. MS
acknowledges support by Research Fellowships of Japan Society for the
Promotion of Science for Young Scientists, by Forefront Physics and
Mathematics Program to Drive Transformation (FoPM), a World-leading Innovative
Graduate Study (WINGS) Program, the University of Tokyo, and by JSPS Overseas
Challenge Program for Young Researchers.

This work presents results from the European Space Agency (ESA) space mission
Gaia. Gaia data are being processed by the Gaia Data Processing and Analysis
Consortium (DPAC). Funding for the DPAC is provided by national institutions,
in particular the institutions participating in the Gaia MultiLateral
Agreement (MLA). The Gaia mission website is https://www.cosmos.esa.int/gaia.
The Gaia archive website is https://archives.esac.esa.int/gaia.

## References

  * Abbott et al. (2019) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2019, Physical Review X, 9, 031040, doi: 10.1103/PhysRevX.9.031040
  * Abbott et al. (2021) Abbott, R., Abbott, T. D., Abraham, S., et al. 2021, Physical Review X, 11, 021053, doi: 10.1103/PhysRevX.11.021053
  * Abdul-Masih et al. (2020) Abdul-Masih, M., Banyard, G., Bodensteiner, J., et al. 2020, Nature, 580, E11, doi: 10.1038/s41586-020-2216-x
  * Andrews (2022) Andrews, J. J. 2022, arXiv e-prints, arXiv:2206.04648.  https://arxiv.org/abs/2206.04648
  * Andrews et al. (2019) Andrews, J. J., Breivik, K., & Chatterjee, S. 2019, ApJ, 886, 68, doi: 10.3847/1538-4357/ab441f
  * Andrews et al. (2021) Andrews, J. J., Breivik, K., Chawla, C., Rodriguez, C., & Chatterjee, S. 2021, arXiv e-prints, arXiv:2110.05549.  https://arxiv.org/abs/2110.05549
  * Andrews et al. (2022) Andrews, J. J., Taggart, K., & Foley, R. 2022, arXiv e-prints, arXiv:2207.00680.  https://arxiv.org/abs/2207.00680
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33, doi: 10.1051/0004-6361/201322068
  * Atwood et al. (2009) Atwood, W. B., Abdo, A. A., Ackermann, M., et al. 2009, ApJ, 697, 1071, doi: 10.1088/0004-637X/697/2/1071
  * Barthelmy et al. (2005) Barthelmy, S. D., Barbier, L. M., Cummings, J. R., et al. 2005, Space Sci. Rev., 120, 143, doi: 10.1007/s11214-005-5096-3
  * Bashi et al. (2022) Bashi, D., Shahaf, S., Mazeh, T., et al. 2022, MNRAS, 517, 3888, doi: 10.1093/mnras/stac2928
  * Binnendijk (1960) Binnendijk, L. 1960, Properties of double stars; a survey of parallaxes and orbits.
  * Bodensteiner et al. (2020) Bodensteiner, J., Shenar, T., Mahy, L., et al. 2020, A&A, 641, A43, doi: 10.1051/0004-6361/202038682
  * Breivik et al. (2017) Breivik, K., Chatterjee, S., & Larson, S. L. 2017, ApJ, 850, L13, doi: 10.3847/2041-8213/aa97d5
  * Bressan et al. (2012) Bressan, A., Marigo, P., Girardi, L., et al. 2012, MNRAS, 427, 127, doi: 10.1111/j.1365-2966.2012.21948.x
  * Buder et al. (2021) Buder, S., Sharma, S., Kos, J., et al. 2021, MNRAS, 506, 150, doi: 10.1093/mnras/stab1242
  * Casares et al. (2017) Casares, J., Jonker, P. G., & Israelian, G. 2017, X-Ray Binaries, ed. A. W. Alsabti & P. Murdin, 1499, doi: 10.1007/978-3-319-21846-5_111
  * Chakrabarti et al. (2022) Chakrabarti, S., Simon, J. D., Craig, P. A., et al. 2022, arXiv e-prints, arXiv:2210.05003.  https://arxiv.org/abs/2210.05003
  * Chawla et al. (2022) Chawla, C., Chatterjee, S., Breivik, K., et al. 2022, ApJ, 931, 107, doi: 10.3847/1538-4357/ac60a5
  * Corral-Santana et al. (2016) Corral-Santana, J. M., Casares, J., Muñoz-Darias, T., et al. 2016, A&A, 587, A61, doi: 10.1051/0004-6361/201527130
  * Cui et al. (2012) Cui, X.-Q., Zhao, Y.-H., Chu, Y.-Q., et al. 2012, Research in Astronomy and Astrophysics, 12, 1197, doi: 10.1088/1674-4527/12/9/003
  * El-Badry & Burdge (2022) El-Badry, K., & Burdge, K. B. 2022, MNRAS, 511, 24, doi: 10.1093/mnrasl/slab135
  * El-Badry et al. (2022a) El-Badry, K., Burdge, K. B., & Mróz, P. 2022a, MNRAS, 511, 3089, doi: 10.1093/mnras/stac274
  * El-Badry & Quataert (2020) El-Badry, K., & Quataert, E. 2020, MNRAS, 493, L22, doi: 10.1093/mnrasl/slaa004
  * El-Badry & Quataert (2021) —. 2021, MNRAS, 502, 3436, doi: 10.1093/mnras/stab285
  * El-Badry & Rix (2022) El-Badry, K., & Rix, H.-W. 2022, MNRAS, 515, 1266, doi: 10.1093/mnras/stac1797
  * El-Badry et al. (2022b) El-Badry, K., Seeburger, R., Jayasinghe, T., et al. 2022b, MNRAS, 512, 5620, doi: 10.1093/mnras/stac815
  * El-Badry et al. (2023a) El-Badry, K., Rix, H.-W., Quataert, E., et al. 2023a, MNRAS, 518, 1057, doi: 10.1093/mnras/stac3140
  * El-Badry et al. (2023b) El-Badry, K., Rix, H.-W., Cendes, Y., et al. 2023b, arXiv e-prints, arXiv:2302.07880.  https://arxiv.org/abs/2302.07880
  * Eldridge et al. (2020) Eldridge, J. J., Stanway, E. R., Breivik, K., et al. 2020, MNRAS, 495, 2786, doi: 10.1093/mnras/staa1324
  * Eyer et al. (2022) Eyer, L., Audard, M., Holl, B., et al. 2022, arXiv e-prints, arXiv:2206.06416.  https://arxiv.org/abs/2206.06416
  * Fu et al. (2022) Fu, J.-B., Gu, W.-M., Zhang, Z.-X., et al. 2022, ApJ, 940, 126, doi: 10.3847/1538-4357/ac9b4c
  * Gaia Collaboration et al. (2016) Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016, A&A, 595, A1, doi: 10.1051/0004-6361/201629272
  * Gaia Collaboration et al. (2018a) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018a, A&A, 616, A1, doi: 10.1051/0004-6361/201833051
  * Gaia Collaboration et al. (2018b) Gaia Collaboration, Babusiaux, C., van Leeuwen, F., et al. 2018b, A&A, 616, A10, doi: 10.1051/0004-6361/201832843
  * Gaia Collaboration et al. (2021) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2021, A&A, 650, C3, doi: 10.1051/0004-6361/202039657e
  * Gaia Collaboration et al. (2022a) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2022a, arXiv e-prints, arXiv:2208.00211.  https://arxiv.org/abs/2208.00211
  * Gaia Collaboration et al. (2022b) Gaia Collaboration, Arenou, F., Babusiaux, C., et al. 2022b, arXiv e-prints, arXiv:2206.05595.  https://arxiv.org/abs/2206.05595
  * Giesers et al. (2018) Giesers, B., Dreizler, S., Husser, T.-O., et al. 2018, MNRAS, 475, L15, doi: 10.1093/mnrasl/slx203
  * Gomel et al. (2022) Gomel, R., Mazeh, T., Faigler, S., et al. 2022, arXiv e-prints, arXiv:2206.06032.  https://arxiv.org/abs/2206.06032
  * Halbwachs et al. (2022) Halbwachs, J.-L., Pourbaix, D., Arenou, F., et al. 2022, arXiv e-prints, arXiv:2206.05726.  https://arxiv.org/abs/2206.05726
  * Hayashi & Suto (2020) Hayashi, T., & Suto, Y. 2020, ApJ, 897, 29, doi: 10.3847/1538-4357/ab97ad
  * Hayashi et al. (2020) Hayashi, T., Wang, S., & Suto, Y. 2020, ApJ, 890, 112, doi: 10.3847/1538-4357/ab6de6
  * Heintz (1978) Heintz, W. D. 1978, Double stars, Vol. 15
  * Holl et al. (2022a) Holl, B., Sozzetti, A., Sahlmann, J., et al. 2022a, arXiv e-prints, arXiv:2206.05439.  https://arxiv.org/abs/2206.05439
  * Holl et al. (2022b) Holl, B., Fabricius, C., Portell, J., et al. 2022b, arXiv e-prints, arXiv:2212.11971.  https://arxiv.org/abs/2212.11971
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Irrgang et al. (2020) Irrgang, A., Geier, S., Kreuzer, S., Pelisoli, I., & Heber, U. 2020, A&A, 633, L5, doi: 10.1051/0004-6361/201937343
  * Janssens et al. (2022) Janssens, S., Shenar, T., Sana, H., et al. 2022, A&A, 658, A129, doi: 10.1051/0004-6361/202141866
  * Jayasinghe et al. (2022a) Jayasinghe, T., Rowan, D. M., Thompson, T. A., Kochanek, C. S., & Stanek, K. Z. 2022a, arXiv e-prints, arXiv:2207.05086.  https://arxiv.org/abs/2207.05086
  * Jayasinghe et al. (2019) Jayasinghe, T., Stanek, K. Z., Kochanek, C. S., et al. 2019, MNRAS, 485, 961, doi: 10.1093/mnras/stz444
  * Jayasinghe et al. (2021) Jayasinghe, T., Stanek, K. Z., Thompson, T. A., et al. 2021, MNRAS, 504, 2577, doi: 10.1093/mnras/stab907
  * Jayasinghe et al. (2022b) Jayasinghe, T., Thompson, T. A., Kochanek, C. S., et al. 2022b, MNRAS, doi: 10.1093/mnras/stac2187
  * Kalogera & Baym (1996) Kalogera, V., & Baym, G. 1996, ApJ, 470, L61, doi: 10.1086/310296
  * Kinugawa & Yamaguchi (2018) Kinugawa, T., & Yamaguchi, M. S. 2018, arXiv e-prints.  https://arxiv.org/abs/1810.09721
  * Kochanek et al. (2017) Kochanek, C. S., Shappee, B. J., Stanek, K. Z., et al. 2017, PASP, 129, 104502, doi: 10.1088/1538-3873/aa80d9
  * Kunder et al. (2017) Kunder, A., Kordopatis, G., Steinmetz, M., et al. 2017, AJ, 153, 75, doi: 10.3847/1538-3881/153/2/75
  * Lallement et al. (2022) Lallement, R., Vergely, J. L., Babusiaux, C., & Cox, N. L. J. 2022, A&A, 661, A147, doi: 10.1051/0004-6361/202142846
  * Lennon et al. (2022) Lennon, D. J., Dufton, P. L., Villaseñor, J. I., et al. 2022, A&A, 665, A180, doi: 10.1051/0004-6361/202142413
  * Liu et al. (2022) Liu, B., D’Orazio, D. J., Vigna-Gómez, A., & Samsing, J. 2022, Phys. Rev. D, 106, 123010, doi: 10.1103/PhysRevD.106.123010
  * Liu et al. (2019) Liu, J., Zhang, H., Howard, A. W., et al. 2019, Nature, 575, 618, doi: 10.1038/s41586-019-1766-2
  * Majewski et al. (2017) Majewski, S. R., Schiavon, R. P., Frinchaboy, P. M., et al. 2017, AJ, 154, 94, doi: 10.3847/1538-3881/aa784d
  * Martin et al. (2005) Martin, D. C., Fanson, J., Schiminovich, D., et al. 2005, ApJ, 619, L1, doi: 10.1086/426387
  * Mashian & Loeb (2017) Mashian, N., & Loeb, A. 2017, MNRAS, 470, 2611, doi: 10.1093/mnras/stx1410
  * Pourbaix et al. (2004) Pourbaix, D., Tokovinin, A. A., Batten, A. H., et al. 2004, A&A, 424, 727, doi: 10.1051/0004-6361:20041213
  * Pourbaix et al. (2022) Pourbaix, D., Arenou, F., Gavras, P., et al. 2022, Gaia DR3 documentation Chapter 7: Non-single stars
  * Ricker et al. (2015) Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2015, Journal of Astronomical Telescopes, Instruments, and Systems, 1, 014003, doi: 10.1117/1.JATIS.1.1.014003
  * Rivinius et al. (2020) Rivinius, T., Baade, D., Hadrava, P., Heida, M., & Klement, R. 2020, A&A, 637, L3, doi: 10.1051/0004-6361/202038020
  * Safarzadeh et al. (2020) Safarzadeh, M., Ramirez-Ruiz, E., & Kilpatrick, C. 2020, ApJ, 901, 116, doi: 10.3847/1538-4357/abb0e8
  * Sana et al. (2012) Sana, H., de Mink, S. E., de Koter, A., et al. 2012, Science, 337, 444, doi: 10.1126/science.1223344
  * Saracino et al. (2022) Saracino, S., Kamann, S., Guarcello, M. G., et al. 2022, MNRAS, 511, 2914, doi: 10.1093/mnras/stab3159
  * Scott (1992) Scott, D. W. 1992, Multivariate Density Estimation
  * Shahaf et al. (2023) Shahaf, S., Bashi, D., Mazeh, T., et al. 2023, MNRAS, 518, 2991, doi: 10.1093/mnras/stac3290
  * Shahaf et al. (2019) Shahaf, S., Mazeh, T., Faigler, S., & Holl, B. 2019, MNRAS, 487, 5610, doi: 10.1093/mnras/stz1636
  * Shao & Li (2019) Shao, Y., & Li, X.-D. 2019, ApJ, 885, 151, doi: 10.3847/1538-4357/ab4816
  * Shapiro & Teukolsky (1983) Shapiro, S. L., & Teukolsky, S. A. 1983, Black holes, white dwarfs, and neutron stars : the physics of compact objects
  * Shappee et al. (2014) Shappee, B. J., Prieto, J. L., Grupe, D., et al. 2014, ApJ, 788, 48, doi: 10.1088/0004-637X/788/1/48
  * Shenar et al. (2020) Shenar, T., Bodensteiner, J., Abdul-Masih, M., et al. 2020, A&A, 639, L6, doi: 10.1051/0004-6361/202038275
  * Shenar et al. (2022) Shenar, T., Sana, H., Mahy, L., et al. 2022, Nature Astronomy, 6, 1085, doi: 10.1038/s41550-022-01730-y
  * Shikauchi et al. (2020) Shikauchi, M., Kumamoto, J., Tanikawa, A., & Fujii, M. S. 2020, PASJ, 72, 45, doi: 10.1093/pasj/psaa030
  * Shikauchi et al. (2022) Shikauchi, M., Tanikawa, A., & Kawanaka, N. 2022, ApJ, 928, 13, doi: 10.3847/1538-4357/ac5329
  * Shikauchi et al. (2023) Shikauchi, M., Tsuna, D., Tanikawa, A., & Kawanaka, N. 2023, arXiv e-prints, arXiv:2301.07207, doi: 10.48550/arXiv.2301.07207
  * Strüder et al. (2001) Strüder, L., Briel, U., Dennerl, K., et al. 2001, A&A, 365, L18, doi: 10.1051/0004-6361:20000066
  * Tanikawa et al. (2020) Tanikawa, A., Kinugawa, T., Kumamoto, J., & Fujii, M. S. 2020, PASJ, 72, 39, doi: 10.1093/pasj/psaa021
  * The LIGO Scientific Collaboration et al. (2021) The LIGO Scientific Collaboration, The Virgo Collaboration, & The KAGRA Scientific Collaboration. 2021, arXiv e-prints, arXiv:2111.03634.  https://arxiv.org/abs/2111.03634
  * Thompson et al. (2019) Thompson, T. A., Kochanek, C. S., Stanek, K. Z., et al. 2019, Science, 366, 637, doi: 10.1126/science.aau4005
  * van den Heuvel et al. (1992) van den Heuvel, E. P. J., Bhattacharya, D., Nomoto, K., & Rappaport, S. A. 1992, A&A, 262, 97
  * van den Heuvel & Tauris (2020) van den Heuvel, E. P. J., & Tauris, T. M. 2020, Science, 368, eaba3282, doi: 10.1126/science.aba3282
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science and Engineering, 13, 22, doi: 10.1109/MCSE.2011.37
  * Vergely et al. (2022) Vergely, J. L., Lallement, R., & Cox, N. L. J. 2022, A&A, 664, A174, doi: 10.1051/0004-6361/202243319
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Weisskopf et al. (2000) Weisskopf, M. C., Tananbaum, H. D., Van Speybroeck, L. P., & O’Dell, S. L. 2000, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 4012, X-Ray Optics, Instruments, and Missions III, ed. J. E. Truemper & B. Aschenbach, 2–16, doi: 10.1117/12.391545
  * Woosley et al. (2002) Woosley, S. E., Heger, A., & Weaver, T. A. 2002, Reviews of Modern Physics, 74, 1015, doi: 10.1103/RevModPhys.74.1015
  * Wright et al. (2010) Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868, doi: 10.1088/0004-6256/140/6/1868
  * Yalinewich et al. (2018) Yalinewich, A., Beniamini, P., Hotokezaka, K., & Zhu, W. 2018, MNRAS, 481, 930, doi: 10.1093/mnras/sty2327
  * Yamaguchi et al. (2018) Yamaguchi, M. S., Kawanaka, N., Bulik, T., & Piran, T. 2018, ApJ, 861, 21, doi: 10.3847/1538-4357/aac5ec