Gaia Data Release 3: Ellipsoidal Variables with Possible Black-Hole or Neutron Star secondaries

As part of Gaia Data Release 3, supervised classification identified a large number of ellipsoidal variables, for which the periodic variability is presumably induced by tidal interaction with a companion in a close binary system. In this paper, we present 6306 short-period probable ellipsoidal variables with relatively large-amplitude Gaia G-band photometric modulations, indicating a possible massive, unseen secondary. In case of a main-sequence primary, the more massive secondary is probably a compact object -- either a black hole or a neutron star, and sometimes a white dwarf. The identification is based on a robust modified minimum mass ratio (mMMR) suggested recently by Gomel, Faigler and Mazeh (2021), derived from the observed ellipsoidal amplitude only, without the use of the primary mass or radius. We also list a subset of 262 systems with mMMR larger than unity, for which the compact-secondary probability is higher. Follow-up observations are needed to verify the true nature of these variables.


Introduction
About 20 dynamically confirmed Galactic stellar black holes (BH) are known to reside in close binary systems. All known stellar mass BHs have been discovered through their X-ray emission, which is fuelled by mass transfer from their non-compact stellar companions (e.g. Fabian et al. 1989;Remillard & Mc-Clintock 2006;Orosz et al. 2007;Ziolkowski 2014). However, a large fraction of BHs in binaries have probably not been detected so far because their optical counterparts are well within their Roche lobes, so that mass is not transferred and X-rays are not generated, making these systems dormant BHs (see the discussion about the frequency of such systems by Breivik et al. 2017;Mashian & Loeb 2017;Yamaguchi et al. 2018;Shao & Li 2019;Yi et al. 2019;Wiktorowicz et al. 2019;Shikauchi et al. 2020). Similar arguments apply to dormant neutron stars (NSs). Only very few known X-ray binaries harbour an NS companion, while most of similar binaries are probably dormant.
Such a dormant binary can be identified either by a large radial-velocity modulation (e.g. Thompson et al. 2019;Zheng et al. 2019;Clavel et al. 2021;Jayasinghe et al. 2022) or by the large stellar ellipsoidal modulation that is induced by tidal interaction with its compact companion (e.g. Rowan et al. 2021). We aim to identify some of the short-period dormant compact systems through this ellipsoidal effect as revealed by the Gaia photometry. A similar study was performed by Gomel et al. (2021c) for the ellipsoidals identified by the OGLE team towards the Galactic bulge ).
Based on the observed ellipsoidal amplitude and the estimated mass and radius of the primary star, one can derive a minimum secondary-to-primary mass ratio, defined as the mass ratio obtained for an inclination of 90 • , provided most of the light is coming from the primary star (e.g. Faigler & Mazeh 2011;Faigler et al. 2015). A binary with a minimum mass ratio significantly higher than unity might have a dormant compact-object companion, which might be a BH, an NS, or even a white dwarf, depending on the estimated minimum mass of the secondary.
Unfortunately, the primary mass and radius are not well known in many cases. Therefore, Gomel, Faigler, & Mazeh (2021a) presented a simple approach that circumvents this problem by suggesting a robust modified minimum mass ratio (mMMR), assuming the primary fills its Roche lobe. The newly defined mMMR depends on the ellipsoidal amplitude and to some extent (by about a few percent) on the effective temperature of the primary, but it does not depend on the mass or radius of the primary star.
The mMMR is always lower than the minimum mass ratio, which in turn is lower than the actual mass ratio. Therefore, binaries with a high mMMR are likely to host a compact object secondary, even if we cannot reliably constrain their primary mass and radius. This is especially true for binaries with an mMMR higher than unity when the periodic modulation is due to the ellipsoidal effect.
We analyse tens of millions of systems that were classified as possible ellipsoidals with varying degrees of certainty by the Gaia DR3 photometric pipeline. We apply our compact com-panion software to analyse the Gaia light curves of these stars to search for systems that might have compact companions.
Our analysis is focused on the short-period binaries and tries to avoid systems with primaries that are not on the main sequence (MS). This is because a high mass ratio is not necessarily an indication of a compact companion for giant stars, for example. As discussed by Gomel, Faigler, & Mazeh (2021a), Algol-type binaries, with sub-giant or giant primaries, are famous counterexamples. The mass ratio of these systems, which probably went through a mass-transfer phase during their evolution (e.g. Fedurco & Parimucha 2018;Chen et al. 2020), can be higher than unity and they can still have an MS secondary (e.g. Negu & Tessema 2018;Samadi Ghadim et al. 2018).
In these binaries, a giant, sub-giant, or stripped giant primary is the lower-mass component, but nevertheless is the brighter star of the system (e.g. Nelson & Eggleton 2001;Budding et al. 2004;Mennekens & Vanbeveren 2017). The recently suggested systems consisting of an evolved primary and a dormant compact object secondary (e.g. Thompson et al. 2019;Liu et al. 2019;Rivinius et al. 2020;Jayasinghe et al. 2021) could be Algoltype binaries or even a binary with stripped giant companion (e.g. van den Heuvel & Tauris 2020; Irrgang et al. 2020;Shenar et al. 2020;Bodensteiner et al. 2020;Mazeh & Faigler 2020;El-Badry & Quataert 2020;El-Badry & Burdge 2022). Therefore, a more restricted list of candidates should only include stars that lie on or near the main sequence, with minimum mass ratios higher than unity. This applies to spectroscopic and photometric candidates alike. The analysis resulted in a catalogue of short-period binary candidates that might have compact companions.
Section 2 lists the constraints we used to define our catalogue of 6306 candidate binaries that have mMMRs higher than 0.5 and therefore might have compact object companions, provided their primaries are on the MS and their modulation is induced by the ellipsoidal effect. Section 3 describes the characteristics of the catalogue, Section 4 examines the CMD location of a subsample of our candidates, and Section 5 cross-matches the catalogue with the available information found in Gaia Data Release 2 (DR2), Simbad astronomical database, the International Variable Star Index (VSX) catalogue and in the Chandra catalogue of X-ray sources. In Section 6, we narrow the catalogue and identify 262 binaries with an mMMR significantly higher than unity, and we present their folded light curves and fitted ellipsoidal modulations. Section 7 presents some details of three stars from the catalogue as examples by comparing the Gaia light curves with OGLE, ASAS-SN, and ZTF photometry. We plot these stars in red in the diagrams of the paper. Finally, Section 8 discusses and summarises our findings.

Defining the catalogue
The compact companion software, which is part of the variability analysis software (VariPipe, Eyer et al. 2022), identified ∼ 20 million possible ellipsoidal systems, obtained by combining multiple classifiers (Rimoldini et al. 2022a), with different ellipsoidal probability selection cuts. We searched this large sample for ellipsoidals with possible compact object companions.
In our analysis we used mainly the G light curve (Rimoldini et al. 2022b), as only a small fraction of systems had enough G BP or G RP measurements to derive significant amplitudes. For these, we were unable to detect any colour-based clustering, and therefore we used the colour information only to remove three outlier systems.
We used the cleaned G light curve, folded with P, twice the period found by a generalised Lomb-Scargle (GLS) period search (Zechmeister & Kürster 2009), as we assumed that the Gaia pipeline detected the predominant second-harmonic periodicity of the ellipsoidal. The period uncertainty, P err , was derived by the pipeline.
We included only systems with more than 25 cleaned G fieldof-view (FoV) transits, with a period within 0.25 < P < 2.5 d, P/P err > 10, and with a frequencygram peak 1 > 12 to ensure that the adopted period was genuine and significant. The upper limit of the second condition was used to avoid systems with non-MS primaries.
As shown by Gomel et al. (2021c), the period of ellipsoidal variables with MS primaries in most cases is shorter than 2.5 days. At this stage, we were left with 112 591 candidates.
Next, a three-harmonic model was fitted to each light curve, with seven free parameters, G, a ic , a is , {i = 1, 2, 3}, that characterise the presumed periodic modulation. Each parameter was derived with its corresponding uncertainty. T 0 was chosen so that a 2s = 0. We defined the amplitude of each harmonic as In the next stage of searching for candidates with compact object companions, we chose only systems with 0.33 < A 2 /(range of G) < 0.6, A 2 /A 2,err > 10, A 1 /A 1,err > 3 or A 3 /A 3,err > 3, and A 1 /A 2 < 1 together with A 3 /A 2 < 0.3.
The first two conditions were set to obtain a reliable and significant A 2 . Third condition reflected our expectation for an ellipsoidal light curve to display non-equal minima, and the harmonic ratio of the fourth condition is typical for ellipsoidals. We were left with a sample of 22 914 systems that are highly probable short-period ellipsoidal variables.
The G-mag histogram of the 22 914 systems is shown in Fig. 1. This magnitude was derived with equation (1) and is slightly different from a simple arithmetic mean. It is evident that most systems are relatively faint, covering a G-mag range of 13-20. In the figure we added the 19th mag mark to denote the observational limit of a few existing multi-object spectrographs that can be used to follow the best candidates (see discussion in Section 8).
The final stage of constructing the catalogue of candidates was to identify ellipsoidals that might have compact companions. This was done using equation (1) of Gomel, Faigler, & Mazeh (2021a), which estimates the leading amplitude A 2 of the ellipsoidal as a function of the fillout factor f . This equation combines the primary volume-averaged radius divided by its Roche-lobe volume-averaged radius (Kopal 1959;Paczyński 1971), the orbital inclination i, and the secondary-to-primary mass ratio q, where L is the average luminosity of the star, L 0 is the stellar brightness without any secondary, and E(q) is the Eggleton (1983) approximation for the volume-averaged Roche-lobe radius in binary semi-major axis units. The ellipsoidal coefficient R. Gomel et al.: Ellipsoidal variables with possible compact companions α 2 depends on the linear limb-and gravity-darkening coefficients of the primary and is expected to be in the 1-2 range. The correction coefficient C(q, f ) starts at 1 for f = 0 (no correction), as expected, and rises monotonically as f → 1, obtaining a value of ∼ 1.5 at f 0.9 (Gomel, Faigler, & Mazeh 2021b). Assuming a fillout factor of f = 0.95, an inclination of 90 • , and a typical α 2 of 1.3 for the G band (Claret 2019), we solved for the mMMR,q min , (Gomel, Faigler, & Mazeh 2021a), based on the observed second harmonic amplitude A 2 . We set f to be 0.95 to be able to use the analytical approximation of Gomel, Faigler, & Mazeh (2021b). The mMMR value obtained for a fillout factor of unity would be lower than the value we derived.
The uncertainty ofq min is inherently high because A 2 increases asymptotically as a function of q and because of the uncertainty in the ellipsoidal amplitude approximation we use. In particular, α 2 is not well known, therefore we adopted its uncertainty to be 0.1. We then derived the uncertainty of A 2 /α 2 , and assuming a Gaussian distribution for this ratio, using equation (3), we obtainedq −1σ min , which corresponds to the 15.9 percentile of A 2 /α 2 . In a similar way, we derivedq −3σ min , which corresponds to the 0.135 percentile of A 2 /α 2 .
As pointed out by Gomel, Faigler, & Mazeh (2021a), the mMMR is expected to be lower than the actual mass ratio of the system because it assumes a fillout factor close to unity and an inclination of 90 o . Therefore, we somewhat arbitrarily opted to include systems withq min > 0.5 in our catalogue of ellipsoidals that might have compact companions. This resulted in 6336 candidates. Different thresholds would have yielded different catalogues. In particular, adopting a different fillout factor, for instance, of f = 0.98 to solve equation (3) would yield a smaller catalogue. The adopted catalogue included 27 sources that were identified as RR Lyrae variables, 25 of which are known RR Lyrae stars in the Milky Way bulge, 24 sources from the OGLE survey (Soszyński et al. 2014) and one source from the VVV survey (Contreras Ramos et al. 2018). Two additional sources are not known in the literature but were confirmed as RR Lyrae stars by inspecting their light curves. Therefore we exlcuded the 27 presumably RR Lyrae variables from the final list. Analysis of the G RP and G BP light curves yielded only a small number of systems with periodic colour information significant enough for our use. We derived the second-harmonic amplitude A 2,RP by using equation (1) for the G RP light curve and inspected the A 2,RP /A 2 distribution for the 815 systems with A 2,RP /A 2,RP,err > 10. We found three outliers with amplitude ratios lower than 0.7 that we excluded from the catalogue. This left 6306 candidate ellipsoidal variables with compact companion secondaries.
The resulting catalogue of 6306 short-period ellipsoidals is given online in the DR3 archival data in table gaia_dr3.vari_compact_companion. The first 15 candidate light curves and parameters are given in Fig. 9 and Table 3.

Catalogue
This section presents some characteristics of the catalogue of the compact companion candidates. Fig. 2 shows the G-mag histogram of the candidates as derived by equation (1). The figure shows that most of the compact companion candidates are relatively faint, in a G-mag range of 15-20. Two brighter sources lie within 13-14 mag, and one of these sources is discussed in Section 7.
A density distribution of the 6306 catalogue stars in Galactic coordinates, coloured according to the orbital period, is presented in Fig. 3. The candidates are spread thorughout the Galactic disk, and most of them are located towards the bulge. Interestingly, the 59 known X-ray binaries (e.g. Corral-Santana et al. 2016) are also concentrated in the Galactic disk, as shown in the figure.
Our ability to identify ellipsoidals, and compact object candidates in particular, clearly depends on the number of FoV transits. This is reflected in Fig. 4, which displays the distribution of the number of G-band FoV transits for the 6306 candidates. Our imposed lower limit of N = 25 is evident. The number of Gband FoV transits varies between 25 to 130, with a median value of 46.
The amplitude of the second harmonic A 2 and the mMMR, q min , are plotted in Fig. 5 as a function of the orbital period P. Each value of A 2 corresponds to a unique value ofq min as ob-  tained from equation (3). The distributions of A 2 ,q min and P are plotted alongside their axes.
As mentioned above, we focused on short-period systems with orbital periods between 0.25-2.5 d, while the lower limit of A 2 ∼ 0.66 mag corresponds to a minimum value ofq min = 0.5, which defines our catalogue of compact companion candidates. Since we chose our catalogue to include only systems with relatively large A 2 and a significant frequencygram peak, we expect a minor contamination, if any, of spurious periods caused by the Gaia scanning law. A slight deficit of low-amplitude systems might be noted in Fig. 5 around P ∼ 0.5 d. This could be a result of our pipeline (Section 2), which was sensitive to the Gaia 6-hour orbital spin.

Location of the candidates in the CMD
As mentioned above, we are interested in ellipsoidal variables with primaries that are not giants or sub-giants. To determine whether this is indeed the case, we plot in Fig. 6 the position of a sub-sample of our candidates in the Gaia colour-magnitudediagram (CMD). Only 513 candidates were included. The relative parallax precision for them was better than 20% ( / σ > 5), the parallax was higher than 0.5 mas ( > 0.5 mas), and a Gaia extinction estimate was available.
The location of the sub-sample in the CMD is shown after correcting for their extinction, taken from Gaia DR3 data (Fouesneau 2022). As a background we plot the Gaia stars that have the same characteristics as the ellipsoidal sub-sample, but have a relative parallax precision better than 10%. The background stars mark the location of the MS, and even the red clump is visible in the upper part of the figure. We also added two theoretical (green) curves that mark the expected edges of the MS, as derived with the ISOCHRONE 2 software package (Morton 2015). Fig. 6 clearly shows that at least for the plotted sub-sample, we succeeded to avoid the giant and the sub-giant branch of the CMD, probably due to the period limit of 2.5 days. Therefore, if the candidates are indeed ellipsoidals and if their mass ratios are higher than unity, they might have compact companions.
We note, however, that the minimum mass ratio of our sample is 0.5 because we chose candidates withq min > 0.5, and therefore the secondary of many of our candidates can still be a MS component. This is probably reflected in the figure, which shows the location of quite a few of the candidates above the MS, possibly due to the light contribution of the MS secondary. To show this point, we plot in the figure (dashed blue curve) the theoretical upper G value that is expected for binaries with two MS stars, which was obtained by shifting the upper bound of the MS by 0.75 mag.
However, we emphasise that CMD locations can contain large uncertainties in absolute G due to errors in parallax, photometry, and extinction estimates, which may be highly inaccurate (Gaia Collaboration et al. 2018;Andrae et al. 2018;Anders et al. 2019). In addition, the extinction uncertainty can shift the G BP − G RP colour index. We demonstrate this point by plotting nine candidates for which the minimum mass ratio is higher than unity, as discussed in Section 6. Their companions are less likely to be on the MS, but the locations of six of them nevertheless appear above the MS. For these, we plot the G 1σ uncertainty, showing that all six systems might still be, within 1 or 2σ, on the MS strip.
We note that the figure suggests that many of the candidates in our sample might have a mass lower than 0.8 M , for example. Consequently, their compact companions might have a relatively low mass in the white dwarf range.

Cross-match with Gaia DR2, Simbad, Chandra, and VSX
In this section we cross-match our catalogue with the Simbad astronomical database, 4 with the Gaia DR2 variability data set (Gaia Collaboration et al. 2018), with the Chandra source catalogue, 5 release 2.0, and with the VSX catalogue. 6 Fig. 5. Second-harmonic amplitude (A 2 ) and mMMR (q min ) as a function of orbital period (P) for the 6306 compact companion candidates. The mMMR was derived using a typical α 2 of 1.3 for the G band and assuming a fillout factor of 0.95. The three sources discussed in Section 7 are shown in red.

Cross-match with the Simbad catalogue
The cross-match with Simbad used a search radius of 3 and resulted in 528 systems. They are summarised in Table 1. The cross-match yielded 522 stars classified as variables and 6 additional systems without a main identifier of variability. Two out of the 6 systems were detected by the Chandra mission (Evans et al. 2010) as having significant X-ray emission. Most systems that our catalogue has in common with the Simbad catalogue were classified by the OGLE team as eclipsing binaries. We inspected a few doezn of them and found in all cases that the OGLE I and Gaia G light curves did not present clear narrow or flat-bottom eclipses. They showed possible very broad two V-shape eclipses, if any, centred around phases 0 and 0.5, suggesting they might either be contact binaries or ellipsoidal variables. It is difficult to differentiate between the two modulations because they are both spread over the entire binaryperiod phase and have similar shapes. Nevertheless, we suggest that the contamination of the catalogue by contact binaries is small, as discussed in Section 8.
We visually inspected the modulation of the two stars identified by OGLE as RR Lyrae variables (Gaia DR3 4116610292178919936, 4107297257038868736) and found that the OGLE classification is more probable, at least for the second object. The fact that only two cases were found as RR Lyrae stars out of 528 cross-matched systems suggests that about half a percent of the stars of our catalogue might be RR Lyrae variables.

Cross-match with Gaia DR2 variables
A cross-match of the candidate list with the Gaia DR2 variability data set (Gaia Collaboration et al. 2018), resulted in two systems, 2024410711652386560 and 241721038597143296. The two sources were classified as rotational variables, with a rotational period that is half the value presented in our catalogue. Our analysis suggests that they are really ellipsoidal variables, however. None of the sources appears in the Simbad database.

Cross-match with the Chandra source catalogue
A cross-match of the candidate list with the Chandra source catalogue (Evans et al. 2010), using a search radius of 1 , resulted in Article number, page 5 of 13 A&A proofs: manuscript no. output Fig. 6. Compact companion candidates in the extinction-corrected Gaia CMD. Only 513 candidates with a minimum parallax of 0.5 mas and parallax-over-error > 5 are plotted. The background grey-scale density map of stars with a minimum parallax of 0.5 mas and a parallax-over-error > 10 serves as a proxy for the expected CMD. Magnitudes and extinction corrections for the candidates and the background stars are from Gaia DR3 data (Fouesneau 2022). Nine systems withq min > 1 are plotted in red (see text). The two green curves mark the theoretical expected edges of the MS. The dashed blue curve is the theoretical upper G value expected for binaries with two identical MS stars, obtained by shifting the upper bound of the MS by 0.75 mag. We show the evolutionary tracks of 0.8 and 1 M stars within their MS phase in black. They were obtained with the ISOCHRONE 3 software. The two sources DR3 4042390512917208960 (upper point) and DR3 4056017172771375616 (lower point) are discussed in Section 7 and are shown in red with a black edge.  three systems that are presented in Table 2. They are all included in Table 1. The third object, Gaia DR3 5966509571940818048, is in the category "Other stars", classified by Simbad as a candidate young stellar object. Its X-ray flux is ∼ 2 × 10 −15 erg/cm 2 /s, with a relatively soft spectrum, consistent with the X-ray luminosity (derived from its Gaia EDR3 parallax) coming from the stellar surface. The other two X-ray sources are too faint for a significant flux measurement.

Cross-match with the VSX catalogue
Following the referee's suggestion, we cross-matched our candidates with the VSX catalogue 7 and found 2044 common systems with derived periodicities and available light curves. Figure  7 shows the VSX periods versus those from Gaia. Most systems (1853) have identical periods in the two catalogues. Some systems (175) display a VSX period that is half the Gaia period, consistent with them being ellipsoidal variables for which the dominant frequency of the modulation is the second harmonic. Three sources (shown as blue circles in Figure 7), Gaia DR3 461521311430643456, Gaia DR3 1834964858132757504, and Gaia DR3 1825053924741093376, analysed by the ZTF project (Bellm et al. 2019), display a clear discrepancy between the Gaia and the VSX periods. We re-analysed the ZTF light curves of the three systems 8 by deriving their corresponding power spectra, as shown in Figure 8. The frequencies of the strongest peaks of the three ZTF power spectra are within 1σ of the corresponding Gaia frequencies.
The spectra are dominated by strong side lobes on the two sides of the strongest peak, separated by ∆ f = 1, reflecting the daily window of ground-based observations. The frequencies of the highest peaks of the three systems, which correspond to half the Gaia orbital periods, are close to a whole number, in units of 1/day, and therefore have one side-lobe peak with a low frequency, as shown in the figure. This side lobe corresponds to the period suggested by the ZTF analysis. This exercise shows the potential advantage of space-mission measurements, which are not subject to the daily cycle of ground-based observations.

Candidates with an mMMR higher than unity
We chose 262 ellipsoidals withq −1σ min > 1, which we considered to have an mMMR significantly higher than unity. This makes them more promising candidates for having compact companions. Forty one of these candidates were found in Simbad, 39 of which are classified as eclipsing binaries and 2 as ellipsoidals by the OGLE project . One of the two ellipsoidals is relatively bright, with a G magnitude of ∼ 13.8, and is discussed in detail in Section 7.
Folded light curves and a summary table for the first 15 candidates, in descendingq min order, are presented in Fig. 9 and Table 3. Similar figures and a table for all 262 candidates are given in the online supplementary document. 9 The data are also available at the CDS. The table includes the Gaia DR3 identifier (id), Gaia orbital period, reference time T 0 [BJD-2455197.5], which was chosen so that a 2s = 0, average G magnitude G, and the cosine and sine Fourier coefficients a ic , a is {i = 1, 2, 3} of the three-harmonic model, defined by equation (1), each with its uncertainty. The table also gives the total number of G-band FoV transits, N, the derived mMMR, and the lower-percentile mMMR.
All included parameters of Table 3, except for N, are given online in the DR3 archival data, in the table gaia_dr3.vari_compact_companion. The archival data also include reference times, averaged magnitudes, and Fourier coefficients with their uncertainties for the Gaia G BP and G RP light curves, the 0.135th percentile mMMR,q −3σ min , and the ellipsoidal coefficient α 2 , defined in equation (3).
Out of the 262 light curves, the light curve of Gaia DR3 4068402346632484864 has an untypical shape for ellipsoidal variables, as shown in Fig. 9 (marked in bold), suggesting that this star is an RR Lyrae variable. This is consistent with the above discussion that our catalogue is contaminated by ∼ 1% of RR Lyrae stars.
A close look at Fig. 9 might give the impression that the two maxima of some of the light curves have different heights. This is reflected by the difference between the coefficient of a 1s and a 3s , indicating either that the light curves are not due to an ellipsoidal modulation or that we are witnessing an additional effect. However, this is probably not the case. The Gaia data are quite sparse and in many cases not sufficient to identify such local small features of the modulation. One good example is Gaia DR3 4056017172771375616, which is discussed in the next section. Fig. 12 shows that the Gaia model, which is based on only 44 measurements, does display two different peaks, while the OGLE light curve, which has 70 points, shows equal maxima.

Three examples of compact companion candidates
This section concentrates on three candidates from the catalogue, examining their photometric modulations as obtained by other surveys, OGLE Pawlak et al. 2016), ASAS-SN (Shappee et al. 2014;Kochanek et al. 2017), and ZTF (Bellm et al. 2019), and compares them with the modulation of Gaia. Then we discuss their candidacy of reliably having compact object companions. This exemplifies the potential of the candidates in the catalogue. The three stars are marked in the pertinent figures above, except for Gaia DR3 4070409432055253760, which does not have an extinction estimate and therefore does not appear in Figure 6. Gaia DR3 4042390512917208960 is the second brightest candidate in our catalogue. It is marked in bold in Table 3 and Fig. 9. Figures 10-12 show the available photometry of the three examples, folded with the Gaia period. It clearly displays the Gaia modulation in the different independent data sets. Table 4 lists some parameters of the photometric analysis of the light curves of the three examples: source identifier, variability classification, second harmonic Fourier coefficient with its uncertainty, and the mMMR and its lower-percentile value.
The first star, Gaia DR3 4042390512917208960, was classified as an ellipsoidal variable by the OGLE team, but as a contact eclipsing binary (EW in their terminology) by the ASASSN-V Gaia Period group (Jayasinghe et al. 2020). Based on the Gaia and OGLE classification, we suggest that the system is an ellipsoidal variable. The relatively small amplitude of the ASAS-SN light curve is probably due to light contamination from a neighbouring star, as expected in the dense Galactic bulge. The second star, Gaia DR3 4070409432055253760, was classified as an ellipsoidal variable by the OGLE team, which reinforces our classification. The two stars have already been analysed by Gomel et al. (2021c) in a search for compact companions. The authors derived similar values forq −1σ min .
The third star, Gaia DR3 4056017172771375616, with the shortest period of the three systems, was classified as a contact binary (EC in their terminology) by the OGLE team. Although it is quite difficult to distinguish between the two types of modulations, the modulation shape is consistent with an ellipsoidal Fig. 9. Folded Gaia light curves in the G band of the first 15 candidates withq −1σ min > 1, in descendingq min order. For each candidate, the Gaia DR3 id is given, together with the Gaia period in days and the value ofq min . The epoch of the second-harmonic minimum corresponds to phases 0 and 0.5, and a three-harmonics model is plotted with a solid red line. All curves are plotted with a mag range of 0.4 for convenience. The residuals are plotted in the lower panels and are shown between -0.1 to 0.1 mag for clarity. Marked identifiers are of Gaia DR3 4068402346632484864, which might be an RR Lyrae variable, discussed in Section 6, and Gaia DR3 4042390512917208960, which is one of our brightest candidates, discussed in Section 7. modulation, and follow-up observations are required to determine the correct variability. Table 6 lists the astrophysical parameters of the three systems, including effective temperature, log gravity, metallicity, radius, mass, and stellar age, as derived by Fouesneau (Gaia DR3, 2022). The last column of the table lists the minimum mass of the secondary, derived from the mMMR of Table 3 and the pri-mary mass. In parentheses we show the 15.9 percentile of the mass, using theq −1σ min value. The second star does not have an estimate of the astrophysical parameters, and therefore no minimum secondary mass was derived.
The two stars with radius and mass in Table 5 display a radius that indicates a slightly evolved state. In principle, we tried to compose our catalogue of candidates with MS primaries alone Article number, page 9 of 13 A&A proofs: manuscript no. output to avoid Algol-type binaries, for which the evolved primary can outshine a more massive companion (see below). Nevertheless, the two examples of Table 5 are still good candidates, as their assumed secondaries, with a minimum mass of 3.2 and 1.9M , are massive enough to be detected in the combined luminosity of the systems if they were MS components.

Discussion
We constructed a catalogue of 6306 variable stars that may have a massive, probably compact, companion on a short-period orbit according to the modulation of their ellipsoidal effect. As shown in Fig. 4, the period distribution of the candidates peaks at ∼ 0.3-0.4 days, with a moderate decline towards longer periods. Our search is clearly more sensitive at short periods, and the decline is therefore probably due to selection effects. Out of the whole catalogue, we selected 262 candidates with an mMMR higher than unity, with a higher probability to have a compact companion.
A&A proofs: manuscript no. output As mentioned above, the reality of the candidacy of the binaries in our catalogue and the sample of stars with a high mMMR in particular, depends on two main assumptions: -The observed periodic variability is due to ellipsoidal modulation. -The star is on the MS.
As shown, the two assumptions are not always fulfilled. Some systems are contact binaries (CB), with two components that are not detached. As pointed out above, it is not straightforward to distinguish between CBs and ellipsoidal variables. For example, Gaia DR3 4042390512917208960, discussed in Section 7, was classified by our pipeline and by OGLE as an ellipsoidal variable, but it is listed as a contact eclipsing binary by ASAS-SN. Nevertheless, we suggest that the contamination of the catalogue by contact binaries is not strong because CBs must by nature have short orbital periods, mostly shorter than 0.25 d, as shown by Rucinski (2010). We have avoided periodic variables with such short periods.
In other cases, the modulations might be due to single-star variability, such as stellar pulsation or rotation. One example is probably Gaia DR3 4068402346632484864, whose folded light curve resembles an RR Lyrae ab-type modulation, as shown in Figure 9. Furthermore, RR Lyrae c-type variables may mimic ellipsoidal variability better. However, we suggest that only a few single-star variables are hidden in the catalogue because the CU7 classifier is highly efficient (Rimoldini et al. 2022a).
Finally, some of the stars in the catalogue might be real ellipsoidal variables but have an MS companion. In some cases, the derived mMMR is lower than unity. We nevertheless included them in the catalogue because the mMMR can be substantially lower than the actual minimum mass ratio, depending on the actual stellar mass and radius of the primary. In other cases, the primary star is slightly evolved and therefore might outshine an MS secondary. However, as shown in Fig. 6, a considerable portion of the catalogue stars with extinction estimates are on the main sequence, and even cases with slightly evolved stars, such as the two examples discussed in Section 7, cannot hide an MS companion.
Nevertheless, despite all cautious measures, many of the candidates of the catalogue may not have a compact object companion. Therefore, only spectroscopic radial-velocity (RV) followup observations can validate the high mass ratio of the candidates. However, even before RV resources are devoted to followup observations, more photometric measurements can be added to the light curve to better study the shape of the periodic modulation, in particular at different photometric bands. This is expected in future Gaia data releases, 10 for which the number of FoV transits is going to double, including the G BP and G RP measurements, which were barely used here. The mass, radius, and temperature of the primary star might also be measured. These can be obtained either from the stellar position on the CMD 10 https://www.cosmos.esa.int/web/gaia/release and/or from its spectral energy distribution, when infrared observations and parallax are available. Any deviation from the massradius-temperature relations of MS stars might indicate that either the primary is not on the MS or that the secondary substantially contributes to the brightness of the system, making the star a less attractive candidate for having a compact companion.
Another indication for a compact companion may come from X-ray observations (e.g. Forman et al. 1978;Voges et al. 1999;Evans et al. 2010;Webb et al. 2020). An example is the discussion of Gaia DR3 5966509571940818048 in Section 5.3. Some dormant compact objects might have temporal X-ray emission (e.g. Remillard & McClintock 2006;Belloni 2010) and therefore have not yet been detected as X-ray sources. Therefore, any X-ray survey or pointed observations of the best candidates can be very useful.
When follow-up RV is performed, only a few measurements with medium precision, on the order of 10 km/s, should be enough to establish the binarity and the minimum mass ratio of each system because the presumed orbital period is known and the expected RV amplitude is on the order of 100 km/s. Unfortunately, Gaia RVS measurements (Katz et al. 2019) are not suitable for the candidates of our catalogue because their stellar brightness and/or temperature are not in the RVS effective range. Thus, our candidates necessitate a dedicated program of followup observations.
An RV follow-up project like this can be done, for example, with multi-object spectrographs such as SDSS-V (Kollmeier et al. 2017) and the upcoming 4MOST (de Jong et al. 2019). The magnitude limit of these spectrographs clearly depends on the exposure time, spectral information of the systems, and on the required RV precision. All in all, we assume that for our purpose, the limit is about 19th mag. This limit is marked in Figs. 1 and 2, indicating that most of our candidates can be followed by these spectrographs.
The catalogue presented here, although the level of contamination is unknown, has the capacity of opening a new window for studying short-period binaries with a compact object companion, either a BH or an NS, and sometimes a white dwarf, when some of the candidates are confirmed by additional observations.