Introduction

The mass distribution of compact objects is closely linked to the evolution of massive stars and supernovae. The final fate of massive stars depends on a number of factors including the progenitor star's composition, mass-loss rates, and past binary interactions (e.g., T. Sukhbold et al. 2016). There are estimated to be ∼10 8 stellar-mass black holes (BHs) and ∼10 9 neutron stars (NSs) in the Milky Way (F. X. Timmes et al. 1996), but only a tiny fraction have been observed.

The majority of the observed compact objects are in X-ray binaries (e.g., D. Mata Sanchez et al. 2025) or were discovered in gravitational-wave mergers (e.g., B. P. Abbott et al. 2016). There are expected to be only a few thousand X-ray binaries in the Milky Way, since most must be in tight orbits to produce observable X-ray emission (J. M. Corral-Santana et al. 2016). The gravitational-wave merger sources are all extragalactic, and detections are biased toward more-massive BHs (M ∼ 20-30 M e ; B. P. Abbott et al. 2019), although systems with lower-mass components have also been identified (R. Abbott et al. 2023). Only a tiny fraction of the BHs and NSs in the Milky Way are expected to evolve into gravitational-wave merger sources.

Instead, the vast majority of BHs and NSs are expected to be either free-floating, isolated systems, or in wide, non-interacting binaries. Isolated BHs can only be detected via microlensing surveys (e.g., C. Y. Lam et al. 2022;K. C. Sahu et al. 2022; C. Y. Lam & J. R. Lu 2023), and the upcoming timedomain surveys with the Roman Space Telescope (D. Spergel et al. 2015) and the Vera C. Rubin Observatory (Ž. Ivezić et al. 2019) could detect tens to hundreds of isolated black hole microlensing events, depending on the final survey configurations (N. S. Abrams et al. 2023;C. Y. Lam et al. 2023).

There are multiple observational methods that can be used to search for compact objects in binary systems. The third Gaia data release (Gaia DR3; Gaia Collaboration et al. 2023c) includes >130,000 astrometric orbits. By combining the Gaia astrometric orbit and radial velocities from Gaia and other ground-based facilities, three BHs have been confidently detected. Gaia-BH 1 (S. Chakrabarti et al. 2023;K. El-Badry et al. 2023b) is binary with a 9.6 M e BH and a Sun-like star, with an orbital period of 185 days. The system is at 480 pc, making it the nearest known BH. Gaia-BH 2 (K. El-Badry et al. 2023a; A. Tanikawa et al. 2023) is a binary with a 8.9 M e BH and a red giant in a longer 1277 days orbit at a distance of 1.16 kpc. Standard binary evolution formation channels struggle to explain the observed orbits of these systems, and alternative formation scenarios from dynamical interactions in clusters have been proposed (A. Tanikawa et al. 2024). Finally, Gaia-BH 3 is the most massive (M = 33 M e ) Galactic BH reported to date, with an 11.6 yr orbit and a metal-poor giant companion in the Milky Way halo at a distance of 590 pc (Gaia Collaboration et al. 2024). This system could have formed as a primordial binary (e.g., K. El-Badry 2024) or through dynamical interactions in a stellar stream (e.g., D. Marìn Pina et al. 2024). The three Gaia BHs are all closer than the majority of X-ray binaries (∼2-10 kpc), suggesting many more systems must exist. Gaia astrometry has also been used to detect 21 NS candidates (K. El-Badry et al. 2024). The systems span a wide range of orbital periods (P ∼ 180-1000 days) and distances (d ∼ 240-1600 kpc).

While future Gaia data releases are expected to continue to expand the population of non-interacting compact object binaries (e.g., K. Breivik et al. 2017), there are other avenues to detect and characterize systems with dark companions. For short-period systems, photometric variability from ellipsoidal variations can be used to identify stars gravitationally distorted by massive companions (e.g., R. Gomel et al. 2021; D. M. Rowan et al. 2021; M. J. Green et al. 2023). The orbital inclination of the system can be measured by modeling the ellipsoidal modulations, so a companion mass can be precisely determined given an estimate of the photometric primary mass. However, this method is prone to false positives from contact binaries (e.g., P. Nagarajan et al. 2023;T. N. O'Doherty et al. 2024) or stripped-star imposters (e.g., T. A. Thompson et al. 2019;K. El-Badry et al. 2022;T. Jayasinghe et al. 2022; D. M. Rowan et al. 2024), and an unambiguous identification of non-interacting BHs or NSs starting from a photometric survey has yet to occur.

Spectroscopic surveys can also be used to identify stars with large-amplitude, periodic radial velocity (RV) variability. For example, B. Giesers et al. (2018) identified a stellar mass black hole in NGC 3201 with mass M  4.4 M e . S. Wang et al. (2024) recently identified a 3.6 M e BH candidate in the Gaia DR3 single-lined spectroscopic binary catalog (SB1), but additional RVs are needed to confirm the orbit. There have also been a number of searches (e.g., T. Yi et al. 2022;H. Yuan et al. 2022) for compact objects in the Large Sky Area Multi-Object Fiber Spectroscopic Telescope (LAMOST; X.-Q. Cui et al. 2012), but no strong BH candidates have been identified.

The primary limitation of spectroscopic-based searches, as compared to astrometric and photometric methods, is that the orbital inclination cannot be constrained by RVs alone. Even with an estimate of the mass of the photometric primary, only the minimum mass of the dark companion can be determined. Systems at lower inclinations with black hole companions produce lower RV amplitudes and could be missed when selecting candidates (e.g., T. Jayasinghe et al. 2023).

The observational methods used to search for non-interacting black hole companions parallel those used to identify and characterize exoplanets. Thousands of exoplanet candidates have been characterized using transits (e.g., L. Cacciapuoti et al. 2022) and RV searches (e.g., R. P. Butler et al. 2017). The Roman Space Telescope is expected to expand the population of exoplanets detected with microlensing considerably (M. T. Penny et al. 2019; S. A. Johnson et al. 2020), and thousands of planets are predicted to be detectable using multiepoch astrometry in Gaia DR4 (M. Perryman et al. 2014). While the photometric, spectroscopic, and astrometric variability signals are orders of magnitude larger for systems with black holes than for exoplanets, stars with black hole companions are expected to be much rarer.

One technique used to search for exoplanets that has not been applied to black hole searches is direct imaging (e.g., C. Marois et al. 2008). Here, we report results from a pilot survey combining RV selection with direct imaging using SHARK-VIS (F. Pedichini et al. 2022Pedichini et al. , 2024) ) on the Large Binocular Telescope (LBT; J. M. Hill et al. 2010). Methods such as speckle imaging (e.g., J. W. Davidson et al. 2024), interferometry (e.g., A. C. Rizzuto et al. 2013;R. Klement et al. 2022;K. Deshmukh et al. 2024), and high-contrast imaging (e.g., N. Shatsky & A. Tokovinin 2002;T. Pauwels et al. 2023) have long been used to identify and characterize binary systems and to search for wide companions to close binaries (e.g., S. Isobe et al. 1992;A. Tokovinin et al. 2006). These systems are important for understanding binary population statistics, which informs binary formation theories (e.g., M. Moe & R. Di Stefano 2017).

Direct imaging can be adapted to search for stars with black hole companions. If we resolve a system with large-amplitude RV variability indicating a massive companion, but find no luminous star at the expected separation, the companion must be a compact object. Repeated observations can then measure the orbit of the photometric primary to determine the inclination. If a luminous star is identified, the component masses can be measured with two direct-imaging measurements at different orbital phases. Many of the binaries most suited for this observing strategy are bright (V  10), nearby stars. Even for systems where the companion is determined to be a star, rather than a compact object, direct imaging can provide the first characterization of the stellar companion, contributing to long-term efforts to understand stellar multiplicity. Section 2 describes the search strategy for our program, and Sections 3 and 4 describe the properties of the instruments and the setup for our observations. We fit RV orbits using archival measurements in Section 5 and show how the combination of high-contrast imaging and RVs can be used to directly measure the binary component masses. We characterize the binaries using broadband photometry in Section 6, and we report our first results on four systems in Section 7, and discuss how this search strategy could be applied to more targets in Section 8.

Target Search

We start by selecting candidates from two catalogs of spectroscopic binaries. The first catalog is the SB9 catalog (D. Pourbaix et al. 2004), which has spectroscopic orbits for more than 2300 systems observed as early as 1901 (R. T. Crawford 1901). We only consider single-lined spectroscopic binaries (SB1s) since BH or NS companions will not produce a second set of spectral lines, and remove systems with a reported velocity semiamplitude for the secondary. Many of the binaries in the SB9 catalog are naked-eye stars, and the median V-band magnitude of the catalog is 8.3.

The second is the SB1 catalog in Gaia DR3, which includes >180,000 systems (Gaia Collaboration et al. 2023a;E. Gosset et al. 2025). The Gaia RVS spectrometer operates in the near-IR (845-872 nm) around the CaII triplet and has R ≈ 11,500 (M. Cropper et al. 2018). The Gaia SB1 catalog has a median G-band magnitude of 11.9. Only the RV orbit model is included in Gaia DR3 (the epoch RV measurements are expected to be released in DR4). Comparisons between the Gaia RV orbit model and archival spectroscopy (D. Bashi et al. 2022) and photometry (D. M. Rowan et al. 2023) have found that roughly half of the solutions are unreliable. The majority of the spurious orbits are at short periods and high eccentricities (see Figure 11 of D. Bashi et al. 2022). The problem is likely due to the large number of possible solutions at short periods when fitting well-separated, sparse RVs. Here, we start with the entire Gaia SB1 catalog and vet individual candidates on a case-by-case basis.

Figure 1 shows both catalogs on a Gaia color-magnitude diagram (CMD). We compute extinctions using the mwdust Combined19 3D dust map (J. Bovy et al. 2016, based on R. Drimmel et al. 2003;D. J. Marshall et al. 2006;G. M. Green et al. 2019). For targets in the SB9 catalogs without Gaia photometry, we estimate their Gaia apparent magnitudes using the V and I or Two Micron All Sky Survey (2MASS) H and K s magnitudes. 11 We remove targets below the main sequence where one or both of the components are likely white dwarfs or subdwarfs. We also only include systems where the parallax ϖ to parallax uncertainty σ ϖ ratio is ϖ/σ ϖ > 20, since accurate distances are required to estimate the projected physical separation.

We compute the binary mass function using the SB9 or Gaia RV orbit solution as

where P is the orbital period, K is the velocity semiamplitude, e is the orbital eccentricity, i is the orbital inclination, and M 1 and M 2 are the masses. The binary mass function is an absolute lower limit on the mass of the unseen companion, M 2 , which corresponds to the limit with a zero-mass photometric primary and an edge-on inclination. We can place more meaningful minimum companion mass limits by estimating the mass of the photometric primary.

We use the StarHorse2 (SH2) catalog (F. Anders et al. 2022) to estimate the mass of the photometric primary. SH2 combines broadband photometry and Gaia astrometry to estimate stellar parameters, such as effective temperature, metallicity, and mass. The SH2 model assumes that the targets are single stars, but since these systems are observed as SB1s, we assume that the flux ratio is low enough that SH2 can provide a reasonable mass estimate. Only ∼63% of the stars in the SB9 catalog have SH2 mass estimates, as compared to 96% of the Gaia SB1 catalog. We predict the masses of stars without SH2 mass measurements using their CMD position. We split the sample of stars with SH2 masses into 80% for the training set and 20% for the validation set. We stratify over the evolutionary state so that the fraction of main-sequence, subgiant, and giant binaries in the training and validation sets is the same as in the total data set. We fit a k-nearest neighbors model (J. Goldberger et al. 2004) with scikit-learn (F. Pedregosa et al. 2011) using the extinction-corrected Gaia color and absolute magnitude. As expected, this model is most effective on the main sequence where stars of different masses are more easily separated in color and absolute magnitude. In the validation set, the rms errors for the recovered masses are 0.11 M e , 0.26 M e , and 0.30 M e for main sequence, subgiant, and giant stars, respectively. The binaries in the validation set with large mass errors are found at higher distances and extinctions compared to the full sample and are more concentrated in the Galactic plane, so the extinction-corrected absolute magnitude is likely more uncertain.

Next, we determine the minimum companion mass, M 2,min assuming an edge-on inclination and the SH2 mass of the photometric primary, M 1 , with Equation (1). We determine the maximum angular separation by combining M 1 , M 2,min , P, and the Gaia distance estimate, d, as

Figure 2 shows the companion mass and the maximum angular separation of all of the targets. Targets with large angular separations and higher minimum companion masses are the most promising for searching for non-interacting black holes with SHARK-VIS. We selected four targets that are observable in the Northern hemisphere, bright (R  8), with

and maximum angular separations of >30 mas for this pilot study. They are labeled in Figure 2, and Figure 1 shows their position on the Gaia CMD.

SHARK-VIS Observations

SHARK-VIS (F. Pedichini et al. 2022Pedichini et al. , 2024) is a 400-1000 nm high-contrast optical imager using the right side SOUL Extreme Adaptive Optics (ExAO) system of the LBT (E. Pinna et al. 2016). With an angular resolution of ∼15 mas at 550 nm (M. Mattioli et al. 2019), it is the ideal instrument to search for binary companions to nearby SB1s. SHARK-VIS adopts a fast imaging approach with short integrations of a few milliseconds so that atmospheric turbulence distortions are roughly constant in each exposure (M. Stangalini et al. 2016). Faint companions can then be detected with post-processing by recentering and stacking the frames to avoid smearing the point-spread function (PSF). The frame sequence is cleaned of detector signatures by a standard dark subtraction followed by frame-by frame removal of the dynamical column-amplifier noise measured using the masked pixels surrounding the field of view. The frames are then registered to the Gaussian fit peak of the PSF cores, numerically de-rotated to the sky orientation, and finally averaged, usually after selecting the best frames based on the peak amplitude and sharpness.

We observed the four targets during the nights of 2024 May 17 and 2024 May 18. Table 1 reports the details of the SHARK-VIS configuration for each target. For these observations, we used the "LMIRCam" entrance dichroic. All of the light below 700 nm and half of the light between 700 nm and 1000 nm is sent to the SHARK-VIS instrument, and the rest is sent to the AO wave front sensor (WFS). The beam passes through the atmospheric dispersion compensator and finally an optical filter. We used narrowband filters for HD 137909 and HD 104438, and R-and V-band filters for HD 117044 and HD 176695. Since the WFS does not measure optical abberations along the optical path between the dichroic and the SHARK-VIS detector, a previously calibrated noncommon path aberration correction must be inserted into the AO control loop. Unfortunately, this correction was not perfect for the observations presented here due to an issue with the AO setup, producing trefoil abberations in the image PSFs. This does not affect the first two targets (HD 137909, Section 7.1 and HD 104438, Section 7.2), where the separation of the components is large (100 mas). In the end, the trefoil abberation also did not hinder our analysis of the companions to HD 117044 (Section 7.3) or HD 176695 (Section 7.4).

In all four systems, we find evidence of a luminous companion. Table 2 reports the properties of the detected companions. For three of our targets (HD 137909, HD 117044, and HD 176695), the companion is clearly visible in the stacked image even before selecting best frames. For the fourth target, HD 104438, angular differential imaging based on principal component analysis (PCA-ADI; A. Amara & S. P. Quanz 2012) is required to detect the companion. We discuss the analysis of the SHARK-VIS images for each target individually in Section 7. Only two of the targets are consistent with a simple stellar binary (HD 137909, Section 7.1 and HD 117044, Section 7.3), and the other two are likely hierarchical triple systems (HD 104438, Section 7.2 and HD 176695, Section 7.4).

Little ILocater (Lili)

During our SHARK-VIS observations, the iLocater pathfinder instrument, Little iLocater (Lili), was deployed on the left side of the LBT. Lili is a precursor to the full iLocater spectrograph that is currently under construction. The primary science objective of iLocater is to characterize exoplanets using extreme precision radial velocities (J. R. Crepp et al. 2016;J. Crass et al. 2022). iLocater will be diffraction limited and take advantage of the LBT AO system to get spatially resolved spectroscopy, enabling characterization of binary systems (e.g., J. R. Crepp et al. 2025) and exoplanet hosts. Both Lili and iLocater utilize the iLocater fiber-injection system (acquisition camera), which was installed at the LBT in 2019 allowing for on-sky fiber coupling operations. Lili is a volume phase holographic grating spectrograph designed to validate the on-sky fiber-injection system for iLocater. Lili operates in the near-IR (0.97-1.31 μm) and has a resolving power of R ≈ 1500 recorded at speeds of up to 400 Hz. We observed three Note. The DIT column reports the single-frame detector integration time.

of our targets (HD 137909, HD 104438, and HD 176695) with Lili. For HD 137909, the companion was clearly detected in the iLocater acquisition camera, so we obtained spectra of both components. For HD 104438 and HD 176695, we only obtained spectra of the primary star. The Lili data was reduced using the standard methodology described in R. J. Harris et al. (2025). In short, after standard image corrections, we optimally extracted the spectra using a halogen lamp flat field to define the spectral profile. We then wavelengthcalibrated the spectra using an Ar lamp, and continuum normalized it using a spline fit. We also used the flat field lamp to remove a component of high-frequency striping noise inherent to the InGaAs detector. For each target, we compare the observed Lili spectra to BT-Settl model spectra (F. Allard et al. 2011) matching the stellar parameters of each system. The model spectra are continuum normalized with a spline fit and convolved with a Gaussian function approximating the Lili line spread function to better match the data. We discuss the Lili spectra for each target in Section 7.

RV Orbits

For each target, we fit the archival RVs using a Keplerian orbit model of the form

where γ is the center-of-mass velocity, K is the radial velocity semiamplitudes, f is the true anomaly, and ω is the argument of periastron. The true anomaly, f, is defined based on the eccentric anomaly, E, and the eccentricity, e, by

and the eccentric anomaly is

where P is the period, and T 0 is the time of periastron. We use emcee (D. Foreman-Mackey et al. 2013) to sample over the orbital parameters with Markov Chain Monte Carlo (MCMC) methods. For targets with RVs from multiple instruments, we include an RV offset, δ RV , to fit for differences in RV zero-points.

We include a Gaussian prior on the RV offsets with μ = 0 and σ = 1.5 km s -1 . We also fit for a stellar jitter term, s, which encompasses additional variation from stellar activity and the effects of underestimated RV uncertainties (A. M. Price-Whelan et al. 2017). We run the MCMC chains for 50,000 iterations with 2n walkers, where n is the number of fit parameters. We discard the first 10,000 iterations as burn-in. Table 3 reports the new estimate of the mass function f (M) and median values for all parameters with 1σ uncertainties from the MCMC posteriors.

Figure 3 shows the MCMC posteriors for HD 137909. We discuss 1. The flux ratio is the ratio between the companion and photometric primary flux. For HD 137909 and HD 104438, the flux ratios are measured in narrowband filters. For HD 117044 and HD 176695, we use R-and V-band filters.

Table 3 RV Orbit Fits to the Four Targets HD 137909 HD 104438 HD 117044 HD 176695 Period [days] -+ 3851 3 3 -+ 13000 90 80 -+ 1890 10 10 -+ 7710 50 40 K [km s -1 ] -+ 9.34 0.06 0.06 -+ 6.7 0.6 0.6 -+ 30 10 40 -+ 11.6 0.2 0.2 Ecc -+ 0.529 0.005 0.005 -+ 0.58 0.05 0.04 -+ 0.7 0.1 0.1 -+ 0.662 0.007 0.007 f (M)[M e ] -+ 0.199 0.004 0.004 -+ 0.22 0.05 0.08 -+ 1.3 1.0 7 -+ 0.52 0.02 0.02 γ [km s -1 ] --+ 21.47 0.05 0.05 -+ 24.2 0.3 0.3 -- + 12 1 2 --+ 1.81 0.09 0.09 ω[°] -+ 179.5 0.6 0.4 -- + 80 10 10 -+ 5 3 7 -+ 50 1 1 T 0 -2.42 × 10 6 [days] -+ 9000 10 10 -+ 1300 200 200 -+ 28560 20 30 -+ 13200 80 100

Note. T 0 is the periastron time. We fit for RV offsets when multiple data sets are used. For HD 137909, δ RV,1 , δ RV,2 , and δ RV,3 correspond to the Lick, CASPEC, and DDO observations, respectively. For HD 104438, δ RV,1 , δ RV,2 , δ RV,3 , and δ RV,4 correspond to APF, CORAVEL, the Cambridge Spectrometer, and the Mount Wilson 60in observations, respectively. Finally, for HD 176695, δ RV,1 corresponds to the DDO observations. The orbital solution for HD 137909 in the SB9 catalog only includes the RVs from F. J. Neubauer (1944). Our updated RV orbit model, which includes more recent measurements, does not agree with the SB9 solution at the 1σ level (see Section 7.1).

the RV orbit for each target in Section 7 and compare it to the original solution reported in the SB9 catalog. With a single SHARK-VIS measurement of the angular separation and an estimate of the photometric primary mass, we can measure the orbital inclination and the companion mass. A Keplerian orbit is defined by seven orbital elements, P, T 0 , e, ω, i, a (semimajor axis), and Ω (longitude of ascending node). This first four, P, T 0 , e, and ω, are measured directly from the RV orbit. The semimajor axis, a, is related to the masses M 1 and M 2 and the orbital period through Kepler's third law. If M 1 and i are known, M 2 is determined by combining K, P, and e (Equation ( 1)). The sky-projected angular separation is

where d is the distance (S. Tremaine 2023). The sky-projected angular separation does not depend on Ω, although Ω can be determined from the position angle of the separation. By combining the RV orbit and one measurement of the angular separation, we can determine i and M 2 for an assumed value of M 1 . With two measurements of the angular separation at Table 3 reports the median posteriors and 1σ uncertainties.

different orbital phases, we can also measure M 1 . We show an example of combining multiple angular separation measurements below in Section 7.1.

Spectral Energy Distributions

We jointly fit the broadband photometry and flux ratio from the SHARK-VIS observations to characterize the binaries. We retrieve Galaxy Evolution Explorer (GALEX; L. Bianchi et al. 2011), Gaia (Gaia Collaboration et al. 2023c), 2MASS (R. M. Cutri et al. 2003), and Wide-field Infrared Survey Explorer (WISE; R. M. Cutri et al. 2012) photometry for our targets, when available. We also use Gaia synthetic photometry (GSPC, Gaia Collaboration et al. 2023b), which is determined from the low-resolution Gaia XP spectra. The GALEX and GSPC u, g, r, i, and z photometry use the AB magnitude system, and the 2MASS, WISE, and GSPC U, B, V, R, and I photometry are in Vega magnitudes. Since our targets are bright, some of the photometry is saturated. We discuss the quality of the individual photometry measurements for each target in Section 7.

We fit the spectral energy distribution (SED) defined by the available photometry using the F. Castelli & R. L. Kurucz (2003) model atmospheres included in pystellibs. 12 We use emcee (D. Foreman-Mackey et al. 2013) to sample over the stellar parameters of both stars. We assume a fixed distance from C. A. L. Bailer-Jones et al. (2021) and extinction from mwdust. We also keep the metallicity fixed and use literature metallicity measurements when available. We fit for two stellar components simultaneously, sampling over the effective temperatures, T eff,1 and T eff,2 , surface gravities, g log 1 and g log 2 , total luminosity, L T = L 1 + L 2 , and luminosity ratio L 2 /L 1 . Finally, we include a prior on the flux ratio based on the SHARK-VIS observations in the filter used for each target. For HD 117044 and HD 176695, we use a prior on both the R-and V-band flux ratios. We run the MCMC for 15,000 iterations with 12 walkers and discard the first 5000 iterations as burn-in. Table 4 reports the results of the two-star SED fits. We discuss the SED fits below.

Discussion of Individual Targets

HD 137909

HD 137909 (β CrB) is a bright (V = 3.7) Ap main-sequence star (Figure 1) with an overabundance of rare-earth elements (D. W. Kurtz et al. 2007) and a complex magnetic field (J. M. Seach et al. 2022). HD 137909 was first identified as a spectroscopic binary in 1907, and the first orbit was published in J. B. Cannon (1912) using 153 spectra taken at the Dominion Observatory between 1910 and 1912.foot_3 F. J. Neubauer (1944) reported 341 RV measurements taken between 1931 and 1942 with the Mills spectrograph at the Lick Observatory (W. W. Campbell 1898). F. J. Neubauer (1944) found an orbital solution of P = 10.496 yr and K = 9.19 km s -1 , and this is the orbital solution reported in the SB9 catalog (D. Pourbaix et al. 2004).

HD 137909 has since been targeted by a variety of spectrographs. K. W. Kamper et al. (1990) collected 121 RVs from the David Dunlap Observatory (DDO). P. North et al. (1998) took 78 RV measurements of HD 137909 with the CORAVEL spectrograph on the 1 m telescope at the Olservatoire de Haute Provence. HD 137909 was later observed by G. Mathys (2017) as part of a program to measure magnetically split lines of Ap stars using CASPEC on the ESO 3.6 m telescope, the AURELIE spectrograph on the 1.5 m telescope at Observatoire de Haute Provence (OHP), and the Kitt Peak National Observator (KPNO) Coudé spectrograph.

We include all 67 measurements in our RV fit and use a single RV offset for the G. Mathys (2017) RVs, since they do not report which observations correspond to which instrument. There are also 93 HARPS RVs of this target (T. Trifonov et al. 2020) from 2007, but they are not zero-point corrected and were all taken in the same night, so we do not include them in our orbit fit.

Figure 4(b) shows the RVs and random samples from the MCMC posteriors, and Table 3 reports the RV orbit fit from the combined RV data sets. Our RV orbit fit prefers a slightly longer orbital period, 3851 ± 3 days (10.55 yr), than the SB9 catalog orbit P = 3833.58 days, which is only based on the F. J. Neubauer (1944) RVs. We also measure a higher orbital eccentricity, 0.529 ± 0.005 versus 0.406 ± 0.025. Our updated RV model leads to a lower binary mass function of f (M) = 0.199 ± 0.004 M e instead of the f (M) = 0.24 M e reported in F. J. Neubauer (1944).

HD 137909 is not included in the StarHorse2 catalog, which uses Gaia DR3 parallaxes and magnitudes, but it was included in the original StarHorse catalog (F. Anders et al. 2019) based on Gaia DR2. The StarHorse mass is

. The Gaia DR3 mass estimate from the FLAME model pipeline is significantly higher, with M 1 = 2.14 ± 0.05 M e .

The primary is a main-sequence star, so the effective temperature should be a good proxy to estimate the stellar mass. However, there are a number of disparate effective temperature estimates in the literature. R. O. Gray et al. (2003) measured T eff = 7624 K using spectra from the 0.8 m telescope at the Dark Sky Observatory. The ELODIE spectra library reports T eff = 8763 K (P. Prugniel et al. 2007). A reanalysis of the archival ELODIE spectra using a principal component analysis based on a set of spectroscopic standard stars found a lower temperature from the combined set of ELODIE spectra A9/F0 star (M ∼ 1.7 M e ) has T eff = 7200-7400 K (M. J. Pecaut & E. E. Mamajek 2013), which is cooler than any of the spectroscopic effective temperature measurements. Since the companion was not detected in the spectra used to derive the RVs, we observed this target with SHARK-VIS to search for evidence of a luminous companion. However, we later found that H. Bruntt et al. (2010) had previously discovered the companion using Very Large Telescope (VLT)/NACO adaptive optics imaging.

We observed HD 137909 with SHARK-VIS on 2024 May 17 with a plate scale of 6.43 mas pix -1 and a narrowband 647 nm filter. We used 1.2 ms detector integration times for a total exposure of 600 s. Figure 4(a) shows the sky-aligned, stacked image. The companion is clearly visible without using PCA-ADI to subtract the photometric primary. In post-ExAO high-resolution imaging, the PSF is characterized by a very complex shape where the diffraction core is surrounded by a number of variable and quasi-static speckles (e.g., M. Stangalini et al. 2016). The lack of a smooth background makes reliable photometric and astrometric measurements within the AO control radius (∼250-300 mas; Figure 4(a)) challenging. We characterize the companion by taking a Bayesian modeling approach based on numerical information field theory (NIFTy; P. Arras et al. 2019). We maximize the posterior distribution for the unknown PSF, A(x, y), and companion location, (dx, dy), and the flux ratio, α,

where I(x, y) is the stacked image after selecting the best 1.5% of frames. We find a companion flux ratio of 0.182 ± 0.001 and a separation 130.56 ± 0.05 mas (Table 2). The separation is less than the predicted maximum angular separation used to select this target as a candidate (205 mas). The detection of a stellar companion means a non-interacting companion can be ruled out. By combining our measurement of the SHARK-VIS angular separation with the earlier observation from VLT/NACO (H. Bruntt et al. 2010) with the RVs, we can break the degeneracy between the primary mass and the inclination and measure the companion mass. Figure 4(c) shows the constraints on M 2 as a function of i from the RV orbit, SHARK-VIS observation, and VLT/NACO observation. We find that M 1 = 2.03 ± 0.08 M e and i = 65°± 3°, and the resulting companion mass is M 2 = 1.49 ± 0.3 M e .

Figure 4(d) shows the SED of HD 137909. The 2MASS IR photometry is saturated, so we instead use J-and K-band magnitudes from M. Morel & P. Magnenat (1978). The WISE W1 photometry is also saturated. Figure 4(d) shows the SED model from H. Bruntt et al. (2010). Since they use the Hipparcos parallax ϖ = 28.60 ± 0.69, we use the same value here rather than the zero-point-corrected Gaia DR3 parallax ϖ = 27.96 ± 0.97. The model under predicts the flux in the Gaia synthetic photometry, but the apparent G-band magnitude is 3.6, and the XP spectra for stars G  5 could be affected by saturation (Gaia Collaboration et al. 2023b). The flux ratio of the model at 647 nm is α = 0.19 ± 0.03, which is consistent with the measured SHARK-VIS flux ratio. The SED parameters are reported in Table 4.

Figure 4(e) shows the Lili spectrum of each component. BT-Settl model atmospheres (F. Allard et al. 2011) are shown for comparison using the effective temperature, surface gravity, and metallicity reported in Table 4. Both components have Paschen lines, and there is evidence for SI and CI lines in the companion, but these features are absent in the spectrum of the primary.

In summary, HD 137909 is a long-period (10.55 yr) binary consisting of two main-sequence stars with spectral types of A5Vp and F2V. By combining our SHARK-VIS measurement with the VLT/NACO observation from H. Bruntt et al. (2010), we measure the component masses to be 2.03 ± 0.08 M e and 1.49 ± 0.03 M e (Figure 4(c)).

HD 104438

HD 104438 is an R ≈ 5.0 red giant (Figure 1). The Gaia photogeometric distance is 95.5 pc (C. A. L. Bailer-Jones et al. 2021). HD 104438 was first characterized as a possible RV variable in the Bright Star Catalog (D. Hoffleit & C. Jaschek 1991) based on RV measurements from W. Harper (1923). 14The SB9 catalog includes 75 RVs, the majority of which are from the CORAVEL spectrograph at Observatoire de Haute Provence or the Cambridge CORAVEL spectrograph obtained between 1992 and 2008 by R. F. Griffin (2008). Since R. F. Griffin (2008) did not report RV uncertainties, we assume σ RV = 1.0 km s -1 for the CORAVEL observations, σ RV = 2.5 km s -1 for observations taken with the Cambridge Spectrometer, and σ RV = 5.0 km s -1 for observations taken with the Mount Wilson 60in (W. S. Adams & A. H. Joy 1923;W. Harper 1923;H. A. Abt 1970). As described above, we include a term for stellar jitter in our RV orbit fits, which can account for underestimated RV uncertainties. HD 104438 was also observed by APOGEE in 2013, extending the baseline of the RV observations. We obtained three high-resolution (R ≈ 80,000) spectra for HD 104438 with the Automated Planet Finder (APF) Levy spectrograph on the Lick Observatory 2.4 m telescope (S. S. Vogt et al. 2014) to refine the orbital ephemeris and improve constraints on the predicted angular separation. The observations used the  ´ 0 . 5 8 . 0 Decker-N slit and cover a wavelength range of 3730-10206 Å. We used the California Planet Survey pipeline (A. W. Howard et al. 2010) to create 1D spectra, and the spectra were then blaze-corrected by fitting polynomials to the continuum in each order. We used 33 orders spanning 4600-7813 Å for the RV analysis that were chosen to avoid wavelengths affected by telluric lines.

We derive RVs using the cross-correlation method with a PHOENIX stellar atmosphere model (T. O. Husser et al. 2013) chosen to match the spectral type of the photometric primary (Table 4). We measure the RV for each echelle order independently and combine them following the method described in S. Zucker (2003). Table 5 reports the three RV measurements. We simultaneously fit the Keplerian RV orbit

Table 5 APF Radial Velocity Measurements for HD 104438 Julian Date RV σ -2,460,000 (km s -1 ) ( km s -1 ) 648.89 29.51 0.2 663.85 29.67 0.2 653.90 29.64 0.2

model to the archival RVs and our APF measurements. The Gaussian priors on the archival RV offsets have μ = 0 and σ = 1.5 km s -1 (Section 5). We use σ = 0.5 km s -1 for the APF offset prior, but we find that the choice does not significantly change the posteriors on the orbital parameters. (2016) also used the APOGEE spectrum of HD 104438 to predict its mass based on its C and N abundances and found M = 1.22 M e . Figure 5(c) shows the companion mass M 2 for HD 104438 as a function of the orbital inclination for three different values of M 1 . The unseen companion would only be massive enough to be a black hole with M 2 > 3 M e if the orbital inclination was 30°. In order for the companion to be a neutron star with M 2 > 1.4 M e , the inclination would have to be <55°for M 1 = 1.2 M e , which corresponds to a probability P(i < 55°) = 0.42 assuming an isotropic distribution of inclination angles. Based on the RV orbit, primary mass estimate, and Gaia distance, we predict that the maximum angular separation is 148 mas (Figure 2), making this an ideal target for SHARK-VIS.

HD 104438 was observed with SHARK-VIS on 2024 May 8 in an Hα filter with a 5 nm bandwidth. We used a larger 1.5× magnification lens to get a plate scale of 4.32 mas pix -1 . The total exposure time was 300 s with a detector integration time of 10 ms and a field rotation of 8 deg, allowing for the use of PCA-ADI to search for a faint companion. The combination of fast imaging with a narrowband filter and the focal extender was possible because the source is bright (G = 5.29).

After subtracting a dark image and dynamical bias correction, all 30,000 frames were coregistered and processed using the PCA-ADI implementation in scikit-learn (F. Pedregosa et al. 2011). This revealed a faint companion at a separation of ∼100 mas (Figure 5(a)). Forty principal components were chosen to maximize the signal-to-noise ratio (SNR) of the companion. We use the NFC technique (C. Marois et al. 2010), the standard for PCA-ADI photometry, to accurately measure the contrast and position of the companion by injecting negative scaled copies of the primary until the companion is canceled. We found a contrast ratio of (1.90 ± 0.02) × 10 -4 and a separation of 107.3 ± 0.4 mas (Table 2). To confirm that the detection is from a stellar companion rather than an instrumental artifact, we inject fake companions at the same angular separation and different position angles with contrast ratios 2 × 10 -4 . The injected sources have similar shapes and brightness to the detected companion. We also run PCA-ADI on the first and second half of the images independently and find the companion at the same position in both data segments.

Figure 5(b) shows the RV orbit and marks the epoch of the SHARK-VIS observation. HD 104438 was observed shortly after periastron, so the separation should be much less than the maximum separation predicted for an edge-on orbit (148 mas). Even if the orbit had i ∼ 30°, we would only expect the companion to be separated by ∼70 mas at the epoch of the SHARK-VIS observation. Figure 5(c) shows the constraints on the companion mass using the RV orbit and the observed SHARK-VIS projected separation of 107 mas. The two curves intersect at ∼17°, which would require M 2 ; 10 M e assuming M 1 = 1.22 M e . This solution is unphysical because a 10 M e star would dominate the SED and be inconsistent with the observed flux ratio. One explanation for the system could be that the faint companion seen with SHARK-VIS is a high-mass ratio binary of a black hole and a low-luminosity M-dwarf or white dwarf. However, it is more likely that the object detected in SHARK-VIS is not the source of the observed RV variability, but is instead a wide tertiary or a chance alignment. We discuss both of these options below.

The source of the RV variability in HD 104438 could either be a stellar companion that is too close to the giant to be resolved or a compact object. Figure 6 shows the expected orbital separation as a function of time for an edge-on orbit and M 1 = 1.22 M e . For an edge-on orbit, we predict an angular separation of 30 ± 10 mas at the time of our observation. The large uncertainty is driven by the lack of RVs near RV minima (Figure 5(b)), which contributes to the uncertainty on K and e and consequently Δθ. Since the nominal separation detection limit is 30 mas for our SHARK-VIS observations, we cannot rule out a luminous companion at smaller separations. The orbital separation will continue to increase (Figure 6), so additional SHARK-VIS observations will constrain the companion that is the source of the long-period RV variability. For example, the expected separation 1 yr after this observation is 70 ± 10 mas.

The low flux-ratio source separated by 107 mas could be a chance alignment of a fainter star or a wide tertiary companion. It is not intrinsically faint with V ≈ 15 mag. First, we fit a single-star SED to the red giant primary to estimate the total luminosity and estimate the possible spectral type of the faint source. We sample over the effective temperature, surface gravity, and luminosity of the giant, and the distance, extinction, and metallicity. We use a Gaussian prior on the distance based on the parallax uncertainty. We also include priors on the T eff and g log of the giant from the APOGEE spectrum. As with the two-star SED models described in Section 6, we use emcee and run the MCMC for 15,000 iterations with 12 walkers and discard the first 5000 iterations as burn-in. The 2MASS photometry is saturated in the J, H, and K s bands and has a "D" quality flag. Similarly, the WISE W1 and W2 photometry is saturated with SNR < 2. The W3 and W4 bands have higher SNR (70.5 and 54.6, respectively) and are tagged as "A" quality photometry. The GALEX NUV is saturated at m NUV = 13.7, but we include the FUV point in our analysis. We also include the Tycho B-and V-band measurements (C. O. Wright et al. 2003) in our fit. Figure 5(d) shows the UV, optical, an IR photometry and the SED fit, and Table 4 reports the parameters and 1σ uncertainties.

Figure 5(d) also shows the SED of a 0.37 M e M-dwarf and WD cooling track models at the distance of HD 104438. We use D. Koester (2010) WD model atmospheres and cooling tracks from A. Bédard et al. (2020). Based on the low flux ratio of the SHARK-VIS detection, the wide companion is most likely an M-dwarf. We assess the orbital stability of such a triple using the criterion from P. Eggleton & L. Kiseleva (1995) 1)) and the red curves show the constraint from the projected angular separation (Equation ( 6)). The two curves intersect for a high companion mass, which is inconsistent with the observed flux ratio, so the observed wide companion is likely a tertiary companion. For comparison, the blue and purple curves show the same curves if we had detected a companion at a separation of 30, 40, or 50 mas. Panel (d): single-star SED fit to the optical and IR photometry of the red giant. The blue square shows the measured SHARK-VIS flux ratio. The red curve shows the SED corresponding to a M-dwarf, and the blue curves show WD model atmospheres for a range of cooling ages. Panel (e): Lili spectrum of the photometric primary (blue) and BT-Settl model atmosphere (black). In addition to the Paschen lines, there are some metal lines between 1000 and 1100 nm.

that the orbit is stable if r > Y EK , where

and a in and a out are the semimajor axes, and e in and e out are the orbital eccentricities. Figure 7 shows the ratio r/Y EK for a range of e out and a out assuming M 1 = 1.22 M e , M 3 = 0.37 M e , and that the inclination of the inner binary is edge-on. The detected tertiary is consistent with a stable hierarchical triple for wide orbits with a out  60 au. We observe the system at a projected separation of ∼10 au, which suggests that the orbit is inclined.

The probability of the detected source being a chance alignment of a background star is near zero. HD 104438 is at l = 169  .7, b = 76  .2, so the local stellar density is low. We calculate the probability of spurious association of a background star in two ways. First, we use Gaia DR3 to determine the probability P s that a star of comparable brightness is found in a random 0.1 mas radius patch corresponding to the angular separation of the M-dwarf. We draw 10,000 random circles with radius 0.1 mas on the sky within 5°of HD 104438 and count the number of stars in Gaia DR3. We find no matches even without imposing an apparent magnitude cutoff, so we place a limit of P s < 1 × 10 -4 . 22 M e and i = 90°. The vertical black lines mark the median and 16th and 84th quantiles, and the red line shows the nominal detection limit for our SHARK-VIS observation. Since the projected angular separation will be increasing until ∼2029, an additional SHARK-VIS observation could be used to rule out a luminous companion responsible for the RV orbit.

Figure 7. Orbital stability criterion for the HD 104438 triple system assuming M 1 = 1.22 M e and i = 90°for the inner binary. Orbits with r/Y EK > 1 are stable, and the red line shows where r/Y EK = 1. The faint companion seen in SHARK-VIS at 107 mas separation is consistent with a stable triple system if it is in a wide (a out  60 au) orbit. We also assess the probability of spurious association using the Besançon model (M. A. Czekaj et al. 2014), which is a stellar population synthesis model of the Milky Way that can produce an estimate of the star counts along the line of sight. The model has 131 stars with V < 15 within the 1 deg patch toward HD 104438, so the probability that there is a random star in a given 0.1 mas region is P s ∼ 1 × 10 -7 . It is therefore unlikely that the detected star is a chance alignment, but additional SHARK-VIS observations are needed to confirm the tertiary companion.

ROSAT identified HD 104438 as an X-ray source with a separation of  16 . 6 (S. Freund et al. 2022) from the optical source. The X-ray luminosity at the Gaia distance is L x = 2.6 × 10 29 erg s -1 . Figure 8 shows the distribution of eROSITA coronal line emitter X-ray luminosities (S. Freund et al. 2024) compared to L x for HD 104438. The X-ray emission could come from either the giant primary or the M-dwarf tertiary, but it is more consistent with the M-dwarf population, and M-dwarf X-ray emitters are more common than red giant X-ray sources.

Finally, we obtained a spectrum of the giant primary using Lili. The wide companion was not detected in the iLocater acquisition camera, and we did not obtain a Lili spectrum of the tertiary. Figure 5(e) shows the Lili spectrum of the photometric primary compared to a BT-Settl model spectra for a giant using the APOGEE effective temperature, surface gravity, and metallicity. The spectrum has a strong H Pa β line and some possible metal lines between 1050 and 1100 nm.

We conclude that HD 104438 is most likely a hierarchical triple system. The inner binary contains the red giant photometric primary, and the secondary is undetected by SHARK-VIS. Additional RVs are needed to refine the ephemeris and place limits on a luminous companion. The SHARK-VIS observations reveal a low flux ratio, 1.9 × 10 -4 (Table 2, Figure 5) star that must be the tertiary based on the phase of the inner binary RV orbit. This tertiary is most likely an M-dwarf and is unlikely to be the chance alignment of a background star. Additional SHARK-VIS observations are necessary to search for the RV counterpart of the inner binary and confirm the tertiary. The RVs are consistent with a long-period eccentric orbit, but the velocity semiamplitude is uncertain because of the lack of observations near periastron. We fit the RV observations with a broad Gaussian prior on K with a mean of 14.6 km s -1 , corresponding to half the range of the observed RV time series, 2). Panel (b): archival RVs from J. L. Halbwachs et al. (2012). As compared to the other three targets, the orbit is less well constrained, especially near periastron, so the posteriors on K and f (M) are large (Table 3). Panel (c): constraints on the companion mass for different values of f (M). Panel (d): two-star SED fit to the multiband photometry with a prior on the R-and V-band flux ratios from SHARK-VIS. and a standard deviation of 10 km s -1 . As expected, the RV posteriors are asymmetric, and RV orbits with high K cannot be ruled out. Our RV orbit solution agrees with the solution reported in the SB9 catalog. It is possible that this system has a large f (M) > 1 M e , but additional RV observations are needed to confirm the orbit and better constrain K. Unfortunately, the next periastron passage is not expected until 2028 January.

The Starhorse2 catalog (F. Anders et al. 2022) reports

. Figure 9(c) shows the constraint on M 2 as a function of the orbital inclination for three different values of f (M). If f (M) = 0.3 M e , the companion could have M 2 > 3.0 M e for inclinations i < 37  .8. The secondary would be more massive than the primary unless the inclination is i < 66  .7. At higher mass functions, the companion mass could be very large, making this an ideal target for SHARK-VIS to search for a luminous companion. Based on the RV orbit and primary mass estimate, the maximum angular separation is 32 mas.

We observed HD 117044 with SHARK-VIS on 2024 May 18 in the R and V bands. The total exposure times were 630 s and 300 s in the R and V bands, respectively, and both used a 2.3 ms detector integration time. Since the companion is expected to be separated by ∼30 mas (Figure 2), we used the 1.5× magnification lens in order to maximize the plate scale (4.32 mas pix -1 ). The companion is clearly visible in the stacked image made by averaging the 3% of frames selected to have the largest Strehl ratios (Figure 9(a)). The trefoil abberation produces three lobes around both stars, which complicates the flux ratio determination, so we applied the Bayesian modeling approach used for HD 137909 (Section 7.1) to characterize the companion. The R-band flux ratio is 0.226 ± 0.002 and the separation is 30.42 ± 0.50 mas (Table 2). The V-band flux ratio is 0.2119 ± 0.0005, and the uncertainties include both the optimization uncertainties and the uncertainties from the choice of box position. Because the RV orbit is poorly constrained, and the projected separation depends on ω and e (Equation ( 6)), we cannot place meaningful constraints on the inclination in Figure 9.

HD 117044 is the faintest of the four targets in the R band, and the broadband photometry is less affected by saturation in the near-IR. The 2MASS J, H, and K s magnitudes all have an "A" quality flag, as do the WISE W1-W4 magnitudes. HD 117044 has also been observed by GALEX and has NUV = 13.0, but since GALEX photometry is saturated and nonlinear for NUV < 14.0, we exclude it from our SED fit. Figure 9(d) shows the optical to IR photometry of HD 117044. We assume the metallicity is solar since there is no previously reported spectroscopic metallicity measurement. We fit a twostar SED model with a prior on the R-band and V-band flux ratios measured by SHARK-VIS. Our SED fit (Table 4) prefers a solution with two main-sequence stars. The primary has an effective temperature =

-+ T 7010 eff,1 10 20 K and radius = R 1  -+ R 1.880 0.007 0.006 , and the secondary has = -+ T 6750 eff,2 60 50 K and  = -+ R R 0.96 2 0.01 0.02

, suggesting that the primary is near the end of its main-sequence lifetime.

In summary, we find that HD 117044 is a binary with two main-sequence stars. Even though the uncertainty on the binary mass function is large compared to our other three targets, we can characterize the binary and reject this as a black hole candidate with a single-epoch SHARK-VIS observation. Additional RVs could be used to better constrain the RV orbit and estimate the masses and inclinations (as in Figure 4(c)).

The flux ratio is the closest to unity in our sample (Table 2), but the companion would likely still be challenging to detect spectroscopically without techniques like spectral disentangling (e.g., R. Seeburger et al. 2024).

HD 176695

HD 176695 is a red giant star at d = 290 pc (C. A. L. Bailer-Jones et al. 2021). J. F. Heard (1956) measured five RVs between 1948 and 1953 from the DDO and identified it as an RV variable. R. F. Griffin (1995) later took 161 more observations between 1973 and 1995 with CORAVEL, confirming the RV variability and constraining the RV orbit. Figure 10(b) shows the combined set of RV measurements and the joint fit to the data. The orbital period is -+ 7710 50 40 days (21.1 yr), the velocity semiamplitude is 11.6 km s -1 , and this solution agrees with the orbit reported in the SB9 catalog. Despite the CORAVEL observations only covering one orbit, the orbit is well constrained because of the large number of measurements and the dense phase coverage.

The StarHorse2 catalog reports a primary mass of

, and no other mass estimates are available. Figure 10(c) shows the constraints on the companion mass as a function of the orbital inclination. Even if the mass of the primary is lower than predicted, M 1 = 1.1 M e and the orbit is edge-on, the companion must be >1.54 M e . If M 2 > M 1 , the companion would also have to be a giant of greater or comparable luminosity, yet no companion is observed in the spectrum (R. F. Griffin 1995). The mass ratio can only be q 1 if the primary is more massive, 2.1 M e and the orbit is edgeon, which would mean the secondary is an A/F star. Based on the RV orbit and primary mass estimate, the maximum angular separation is 38 mas for an edge-on orbit.

We observed HD 176695 with SHARK-VIS on 2024 May 17 in the R and V bands using a 3 ms integration time for a total exposure time of 300 s with a plate scale of 6.43 mas pix -1 for both filters. Figure 10(a) shows that the companion is visible in the registered and sky-oriented average stack built from the best 3% (47%) of frames in the R band (V band). We use our Bayesian modeling approach (Section 7.1) to determine a companion R-band flux ratio of (9.14 ± 0.02) × 10 -2 and a separation of 42.54 ± 0.02 mas (Table 2). The V-band flux ratio is (10.0 ± 0.6) × 10 -2 . As with HD 117044, there is trefoil abberation around the primary in the SHARK-VIS image. We can distinguish this instrumental artifact from the astrophysical companion because the abberation is fixed with respect to the telescope and has the same orientation in all frames before derotation.

Figure 10(c) shows the companion mass as a function of orbital inclination from the RV orbit and projected separation constraints for different estimates of the primary mass. The measured angular separation (42.54 ± 0.02 mas) requires inclinations of ∼65°-75°and companion masses between 1.8 and 2.0 M e for a 1.1 and 1.8 M e primary. An additional SHARK-VIS observation could measure M 1 as well (See Section 7.1).

We fit the Gaia, 2MASS, and WISE photometry with a twostar SED model with a prior on the R-and V-band flux ratios from SHARK-VIS. There is no spectroscopic metallicity, so we assume [Fe/H] = 0.0 for the two-star SED fit. Figure 10(d) shows the two-star SED fit to HD 176695. The MCMC model finds a total luminosity of 66 L e where the companion has T eff,2 = 5500 K and luminosity L 2 = 5.7 L e . For MIST 2). The bright component to the lower left is an artifact of trefoil abberation (Section 3). Panel (b): RV orbit of HD 176695 based on J. F. Heard (1956) and R. F. Griffin (1995) RV measurements. The vertical line marks the time of the SHARK-VIS observation. Panel (c): constraints on the companion mass for different combinations of M 1 and i from the RV orbit (black) and SHARK-VIS angular separation (red). The two intersect for inclinations 65°-75°, which corresponds to companion masses between 1.8 and 2.0 M e . Panel (d): two-star fit to the multiband photometry with a prior on the R-band flux ratio from SHARK-VIS. Figure 11 shows the MIST evolutionary tracks, and we find that this solution is not compatible with the RV orbit. It is more likely that the HD 176695 is a triple system, and Figure 12 shows an example of a possible SED for a triple system model. Panel (e): Lili spectrum of the photometric primary (blue) and BT-Settl model atmosphere (black). The giant has strong Paschen lines and some metal lines between 1000 and 1100 nm but is otherwise relatively featureless. evolutionary tracks (J. Choi et al. 2016;A. Dotter 2016) of solar-metallicity stars of different masses, we find that this combination of effective temperature and luminosity is consistent with a ∼1.35 M e companion (Figure 11).

However, a companion mass M 2 = 1.35 M e is not consistent with the observed RV orbit and angular separation. From the RV orbit constraint alone, if the companion was 1.35 M e , the primary would have to be 0.85 M e . Such a system can be rejected because the companion would be more evolved than the primary and the total luminosity would be much less than the 66.0 L e needed to match the observed photometry. A system with M 2 = 1.35 M e would also not match the combined projected angular separation, even as M 1 goes to zero.

This means that the most likely companion to HD 176695 is an inner binary. With only optical+IR photometry and the flux ratio in two optical bands, we cannot fit for the parameters of the inner binary components. We instead combine MIST evolutionary tracks and F. Castelli & R. L. Kurucz (2003) model atmospheres to search for systems that match the constraints from the RV orbit, flux ratio, projected separation, and combined SED. Figure 12 shows an example system with M 1 = 1.5 M e , M 2 = 1.37 M e , and M 3 = 0.56 M e . The SED is dominated by the evolved primary, and M 2 is near the end of its main-sequence lifetime, contributing the 9% to the R-band flux necessary to match the SHARK-VIS flux ratios. The low-mass star in the inner binary contributes negligibly to the observed flux, but adds enough mass to make M 2 + M 3 > M 1 , satisfying the RV orbit and observed separation constraints. We emphasize that this is simply a possible solution. More observations are needed to characterize the inner binary. Ultraviolet observations could be used to search for evidence of the inner binary component with M 2 ∼ 1.37 M e since it dominates the flux at these wavelengths. Unfortunately, since HD 176695 is near the Galactic plane (b = 10  .9), it is not included in the GALEX survey, and no other archival UV observations are available.

The tertiary in this proposed model (M 3 = 0.56 M e ) is most likely a main-sequence star rather than a white dwarf. This is because the progenitor of a 0.56 M e WD would have a zeroage main sequence mass of ∼1.1 M e based on the IFMR (M. A. Hollands et al. 2024). Such a progenitor would not have formed a WD before the M 1 = 1.5 M e and M 2 = 1.37 M e stars. However, this triple system model is sensitive to the assumed primary mass, and additional SHARK-VIS observations could be used to determine M 1 and measure M 2 + M 3 .

Finally, we also obtained a spectrum of the photometric primary using Lili. As with HD 104438, the companion was not detected in the acquisition camera, and we did not measure a spectrum of the companion. Figure 10(e) shows the Lili spectrum of the red giant and a model atmosphere for comparison. The spectrum is relatively featureless, aside from the Paschen lines.

In summary, HD 176695 is most likely a hierarchical triple. Based on the RV orbit, the companion must be more massive than red giant primary unless the giant is massive >2 M e . However, such a binary is not consistent with the observed SED or the SHARK-VIS flux ratio. The only way to construct a system that matches the SED, RV, and SHARK-VIS flux ratio constraints is for the companion be a binary with an unequal mass ratio (Figure 12).

Discussion and Conclusions

We observed four bright (R  8) nearby (d  300 pc) SB1s in the SB9 catalog of spectroscopic binaries. We selected these targets based on their mass functions and projected angular separations (Figure 2) as candidates for hosting massive, dark companions. We used SHARK-VIS to directly resolve these systems and look for evidence of a second luminous star in the system and use Lili for spectroscopic characterization.

For all four systems, we identify a luminous companion and rule out a BH/NS for three systems. We model the archival RV observations and the multiband photometry to characterize the companions. Only two of the four systems, HD 137909 and HD 117044, are consistent with simple stellar binaries. We find that HD 104438 and HD 176695 are likely triples.

HD 137909 (Section 7.1) is a massive main-sequence binary. The companion was previously identified with VLT CONICA infrared observations (H. Bruntt et al. 2010). We measure an SHARK-VIS flux ratio of 0.182 ± 0.001, which is consistent with the H. Bruntt et al. (2010) SED model, where the companion is another slightly less-massive main-sequence star (Figure 4(d)). By combining the RV solution and the orbital separation measurements from VLT/NACO and SHARK-VIS, we measure the orbital inclination and both component masses (Figure 4(c)).

HD 104438 (Section 7.2) is a red clump star, and no companion has been previously detected. The RV orbit is long, P = 35.6 yr, and we identify a faint companion separated by 107.3 mas. Based on the RV orbit, this star is not the source of the RV variability and is likely a wide tertiary (Figure 5(c)). We find that it is unlikely to be a chance alignment with a background star. The small flux ratio (1.90 ± 0.02) × 10 -4 (Table 2) implies that it is most likely an M-dwarf (Figure 5). Since the binary companion to the red giant is not detected, we cannot rule out a compact object companion. Additional RV observations near RV minima would refine the orbital solution, and 1.45 M e colored by stellar age. The blue point shows the effective temperature and luminosity derived for the secondary of HD 176695 using the two-star SED fit with a prior on the SHARK-VIS flux ratio (Figure 10(d)). The T eff,2 and L 2 indicate a 1.35 M e companion, but this is not consistent with the RV orbit, so the companion is instead most likely a binary (Figure 12). but the next RV minimum is not until 2056. Since the orbit is now past periastron and the separation is increasing (Figure 6), it is more promising to use additional SHARK-VIS observations to search for a luminous companion.

HD 117044 (Section 7.3) is a main-sequence star, and the RV orbit is poorly constrained, so the mass function is uncertain. We identify a luminous companion in the SHARK-VIS image separated by 30.42 ± 0.50 mas (Figure 9), and it is consistent with a second, lower-mass main-sequence star, suggesting that the mass function is low, f (M) ∼ 0.3 M e . The binary has the smallest projected angular separation in our sample, but the companion is clearly detectable in the SHARK-VIS image even by eye.

Finally, HD 176695 (Section 7.4) is an evolved giant with a large mass function, f (M) = 0.52 M e . Figure 10(c) shows that the combination of the mass function and the measured angular separation requires M 2 > M 1 . The companion should then also be a red giant, but the SHARK-VIS observations reveal a companion with an R-band flux ratio of (9.14 ± 0.02) × 10 -2 (Table 2). While we can find a solution by fitting the SED with a two-star SED model using a prior on the R-and V-band flux ratios, comparisons to MIST evolutionary tracks show that such a system is not consistent with the RV orbit (Figure 11). Instead, we find that it is more likely that HD 176695 is a triple system, with an unequal-mass inner binary made up of an intermediate-mass and low-mass star. Figure 12 shows an example of the SED fit and evolutionary tracks for an example system that meets the constraints of the RV orbit, combined SED, and SHARK-VIS flux ratio. Additional UV photometry and flux ratio measurements at other wavelengths are needed to confirm and characterize the inner binary.

High-contrast imaging is typically used to search for and characterize exoplanets (e.g., G. Chauvin et al. 2004;C. Marois et al. 2008). The planetary companions are orders of magnitude fainter than the host stars, and techniques such as PCA-ADI are needed to carefully model and remove the signal of the host star. The observations of the binary targets described here are much more straightforward. We expect the luminous stellar companions to have a flux ratio of 1% for evolved photometric primaries and 10% for main-sequence primaries. For three out of four systems, the companions are easily identifiable in the stacked image or after a best frames selection. PCA-ADI was only necessary for identifying the likely tertiary companion to HD 104438 (Figure 5). While it is possible that luminous companions to these system could have been identified with multiple high-resolution spectra using techniques like spectral disentangling (e.g., S. Ilijic et al. 2004;R. Seeburger et al. 2024), here we use only a single SHARK-VIS observation with a typical total time on-source of ∼10 minutes to directly detect the luminous companion. This is especially useful for long-period systems or systems with uncertain mass functions, like HD 117044, where multiple years of RV observations would likely be needed to confirm the orbit and search for a black hole companion. Only one SHARK-VIS observation is needed to confidently reject this system as a non-interacting black hole candidate. Since the observation time is short, future SHARK-VIS observations will be done in multiple filters/bands to better constrain the companion SEDs.

Combining astrometric and spectroscopic orbits to constrain binary orbits is a well-established technique (e.g., C. L. Morbey 1975). High-contrast imaging has also been combined with RVs to measure the masses of planets, brown dwarfs, and stars (e.g., J. R. Crepp et al. 2012;A. Boehle et al. 2019;T. D. Brandt et al. 2019). Although the primary objective of our survey is to identify stars with compact object companions, the systems presented here show how SHARK-VIS can be used to make benchmark stellarmass measurements. Combining the SHARK-VIS separation with the radial velocity mass function determines the companion mass M 2 and orbital inclination i given the mass of the primary. With a second measurement of the projected orbital separation, as was available from the earlier VLT/NACO observations of HD 137909 (Section 7.1), the primary mass can be directly determined as well. This demonstrates that it would be straightforward to use SHARK-VIS to fully characterize the Figure 12. Left panel: spectral energy distribution (SED) of HD 176695 compared to a hierarchical triple-star SED model. We use MIST evolutionary tracks to explore triple systems with different component masses and search for systems that match the SED, SHARK-VIS flux ratio, and RV orbit constraints. In this example, the primary star is 1.5 M e , and the inner binary consists of a 1.37 M e and a 0.56 M e star. The colored lines show the contributions of the various components to the inner binary and the combined SED, and the blue square shows the expected flux ratio between the photometric primary and the inner binary based on the SHARK-VIS observation. Right: evolutionary tracks for the same system. The primary is an evolved red giant and the more massive component of the inner binary is near the end of its main-sequence lifetime.

properties of the 100 systems with predicted separations larger than ∼30 mas in Figure 2 using two SHARK-VIS observations. Two of the four systems we observed are most likely triple systems. Triples have been previously suggested as false positives in the search for non-interacting black holes and neutron stars (e.g., E. P. J. van den Heuvel & T. M. Tauris 2020), but by far the most common false positives have been stripped stars (e.g., T. Jayasinghe et al. 2022) and massive white dwarfs (e.g., M. A. Tucker et al. 2024). Direct imaging surveys are less likely to find stripped stars, since post-mass transfer binaries will in general have smaller orbital separations than we can resolve with SHARK-VIS or similar instruments. Massive white dwarfs could be detected through this method, and multiple observations could be used to place constraints on the orbital inclination and measure the companion mass more precisely.

These four systems highlight how SHARK-VIS can efficiently search for stars with compact object companions. One of our targets remains a candidate for hosting a compact object, and we show how two measurements of the projected separation can break the degeneracy between primary mass, M 1 and inclination, i. There are >100 other candidates in the SB9 and Gaia SB1 catalogs (Figure 2) that can be observed by SHARK-VIS. Other instruments, such as the Large Binocular Telescope Interferometer (P. M. Hinz et al. 2016) can characterize systems at angular separations 40 mas, and instruments with smaller angular resolution, such as VLTI GRAVITY (GRAVITY Collaboration et al. 2017), could be used to expand the search to binaries with separations <10 mas. Future Gaia data releases are also expected to expand the sample of binaries with astrometric and spectroscopic orbits. Some of these will likely also be candidates for direct imaging where a luminous companion can be identified with a single observation.