© The Authors 2025

1 Introduction

Massive stars are well known to display a higher multiplicity than lower-mass stars (e.g., Moe & Di Stefano 2017; Offner et al. 2023). In addition, the components of a massive binary more often display similar masses (mass ratio q = M₂/M₁ ∼ 0.5-1, with M₁ the mass of the most massive primary star, and M₂ the mass of the least massive secondary star). Nevertheless, systems with extremely low mass ratios, q < 0.2, do exist. They are typically expected in two cases: (1) They are born as such, that is, they simply represent the tail of the mass ratio distribution. (2) They evolve into a system with an extremely low mass ratio, that is, they are the products of mass-transfer events.

Both types of systems have been observed. Moe & Di Stefano (2015) reported on the discovery of short-period binaries in the Large Magellanic Cloud (LMC) with extreme mass ratios (16 cases with q < 0.2). Although radial velocity data were missing, the authors considered the low-mass secondaries to be cool pre-main-sequence (PMS) stars and called the systems “nascent eclipsing binaries with extreme mass ratios”. Several Galactic objects have been found to belong to this category (Jerzykiewicz et al. 2013, 2015, 2021; Stevens et al. 2020; Johnston et al. 2021; Stassun et al. 2021; Nazé et al. 2023; Pigulski 2024). As these Galactic studies combined spectroscopy and photometry, it was demonstrated beyond doubt that the low-mass secondaries were also the coolest stars in the systems, which confirmed their PMS status. The typical PMS signature for the Ca triplet (at 8498, 8542, and 8662 Å) was even detected in HD 25631 thanks to a large high-quality spectroscopic dataset (Nazé et al. 2023).

Massive post-interaction systems with low-mass secondaries have also been identified. In this case, the secondary star was initially the most massive object of the binary. It then lost its envelope in a mass-transfer event that left behind a hot low-mass stripped star. Several examples were identified in the UV range because the contrast at these wavelengths is better (Wang et al. 2021; Götberg et al. 2023). Interferometric and spectroscopic data were also used to characterize them (e.g., Mourard et al. 2015; Klement et al. 2024). Most of these stripped stars were found on long-period (P ~ 100 d) orbits around Be stars, but at least one star is known in the Galaxy to form a short-period pair with an O star (HD 96670; Nazé & Rauw 2025). In this context, it may be noted that systems in intermediate stages (i.e., with an ongoing interaction or just after it came to an end) were also detected (e.g., Harmanec et al. 2015; Shenar et al. 2020; for a review on Be binaries, see Labadie-Bartz et al. 2025).

While systems with these extreme mass ratios may be rather rare, they are of the utmost importance, and it is essential to assess their evolutionary stage. For example, they can provide information on the origin of the rotation of massive stars. The primordial binaries inform us on the initial distribution of rotational properties, whereas the post-interaction cases reveal the distribution changes caused by angular momentum transfer.

While the tail of massive fast rotators is often considered to be the product of binary interactions (e.g., de Mink et al. 2013), the detection of primordial fast rotators indicates that at least some of them are born as such (Britavskiy et al. 2024).

Observationally, both types of systems display a large contrast between the binary components, in mass and in effective temperature (hence, also in brightness). This usually makes the secondaries difficult to detect. Only the lines from the massive primary star are usually detected in the spectra, and the Doppler shifts caused by the orbital motion of the primary have a low amplitude because the companion mass is low. Large differences in the effective temperature might lead to strong reflection effects when periods are short, however, which might provide a way to identify them in photometric light curves (Pigulski 2024). We identify in this paper a new set of systems with these effects. All cases also display eclipses, which allows us to lift the uncertainty caused by inclination and to further determine the true characteristics of the stars. Section 2 provides the information on the sample creation and on the data we used in the analysis. Section 3 details the results, and Section 4 discusses their interpretation.

2 Sample and data

2.1 Sample

Since 2018, the Transiting Exoplanet Survey Satellite (TESS) mission (Ricker et al. 2015) provides high-quality photometry for nearly the whole sky. This makes it the largest provider of space-based light curves. At the start of this research work, several catalogs of eclipsing binaries using TESS data were already available (IJspeert et al. 2021; Prša et al. 2022). We extracted from them the list of stars with spectral type O or B (down to B3). In parallel, we also examined the objects that passed the so-called TESS threshold crossing event (TCE)¹, i.e. objects identified as displaying “transits” in the automatic processing of TESS light curves, up to Sector 71. Using Simbad², we identified the TCE cases associated with stars with spectral types earlier than B3 (inclusive). Finally, we compiled a list of OB stars from the general catalogs of Maíz Apellániz et al. (2013), Reed (2003), and Skiff (2009). Stars in these catalogs with a Gaia magnitude R_p > 13 mag were discarded, as were sources in crowded regions (i.e., those not providing at least three quarters of the total flux in the R_p band of all Gaia DR3 sources within one arcminute around the target). These three lists were merged, and the TESS light curves were then extracted and examined for each star one by one. Our sample was built from cases with short periods (P < 10 d) and clearly pronounced reflection effects in addition to eclipses. Cases that were already known and analyzed in the literature were discarded.

After we identified the targets, we retrieved their V magnitudes from Simbad and photogeometric distances from Bailer-Jones et al. (2021). For all stars, the parallaxes appear to be reliable. Their values in the Gaia DR3 catalog are higher than five times the parallax uncertainties (RPlx > 5) and the renormalized unit weight error (RUWE) parameter is lower than 1.5. The reddenings were estimated from dust extinction maps (Lallement et al. 2022). The online tool G-TOMO³ provides A_V extinctions at the positions (i.e., RA, DEC coordinates and distances from the Gaia catalog) of the targets, which were then converted into reddenings E(B - V) using the extinction factor R_V = 3.1. At first, the bolometric corrections were evaluated using the Simbad spectral types and a standard calibration⁴ . Initial L_BOL estimates were then produced and compared to those expected for the Simbad spectral types in the same calibration. Adjustments of the spectral types could then start because they appear clearly too late or too early in some cases. Further refinement was performed considering the effective temperatures derived from the average spectra (see below) or as available in the literature. The bolometric corrections for these temperatures were estimated from the formula (Eq. (6) and Table D.13) of Pedersen et al. (2020). The final spectral types and bolometric luminosities are provided in Table A.1. This table also notes the references reporting the eclipsing nature of the targets.

Some additional information from the literature must be mentioned. Radial velocity (RV) variations were reported for V1061Cen by Feast et al. (1957). AN Dor was reported as a SB2-type binary by Gebruers et al. (2022), and Southworth & Bowman (2022) further found an inclination i = 80.8°, an eccentricity e = 0.054, and a periastron argument ω = 322°. Tkachenko et al. (2024) characterized HD 350685 as a SB1 with e = 0 and a velocity amplitude K₁ = 34.1 kms⁻¹, and HD309317 as SB2 with e = 0, K₁ = 34.5 kms⁻¹ and K₂ = 78.7 kms⁻¹. LS I+61 275 was mentioned as type SB1 by Huang & Gies (2006a), but the authors did not derive an orbital solution because only two spectra were available to them. We discuss these cases in more detail below. None of our targets belongs to a cluster, except for LS I +61 275, which is part of IC 1805 (whose age was estimated to be 3.5 Myr by Sung et al. 2017).

One star was added to our sample, that is, the binary V1208 Sco, which belongs to the NGC 6231 cluster (3.5-5 Myr, Banyard et al. 2022). It has long been known to be an eclipsing and spectroscopic binary (Levato & Morrell 1983; Balona & Engelbrecht 1985). Only spectroscopic data were analyzed in detail, however. The most recent SB1 solution was provided by Banyard et al. (2022), and the SB2 nature was reported by Banyard et al. (2023). No photometric modeling was done, and we therefore decided to include it here. The extinction of this system could not be derived from the dust maps, and we therefore evaluated it by comparing the observed and expected B - V color for the reported spectral type (B0.5V).

Finally, it is important to note that two stars were excluded from our sample: LsV+20 28 and HD 298425. The first was reported to be eclipsing by IJspeert et al. (2024), but the Kepler light curves clearly revealed that the eclipsing object is the neighboring star LS V +20 27 (LaCourse et al. 2015 even further suggested that the eclipsing system is in fact a very close neighbor of LS V +20 27, not LS V +20 27 itself). It only appears to be eclipsing in TESS data because the two objects are close (32″) and the light from LS V +20 27 contaminates the TESS data of LS V +20 28 (but not the Kepler data because the Kepler PSF is better than that of TESS). Gaia data indicated for the second star that the eclipsing binary is a neighbor at 19.5" of HD 298425 (Mowlavi et al. 2023). Both cases can thus be securely discarded from the category of eclipsing binaries.

2.2 Photometry

The TESS satellite observes sectors of 24° × 96° on the sky for about one month before it shifts to a neighboring region. All our targets were observed in several such sectors (see the last column of Table A.1), but not always with the same mode. When they were preselected, the targets were observed with a 2 min cadence (sectors marked with an asterisk in this table), while other targets were imaged with a cadence of 30 min in Sectors 1-26, 10 min for Sectors 27-55, and 200 s since Sector 56.

A pipeline that was specifically designed for TESS data (and was inspired by the pipeline designed for the Kepler mission) performed the main reduction steps (pixel-level calibration, background subtraction, flat-fielding, and bias subtraction). For data with a 2 min cadence (i.e., preselected targets), the MAST portal⁵ also provides so-called conditioned light curves (PDC-SAP) in addition to the simple aperture photometry (SAP) time series. These PDCSAP curves are further corrected for crowding, the limited size of the aperture, and instrumental systematics. For each target, we kept only the best-quality (quality flag = 0) data points, and we compared the two datasets. They were often indistinguishable, or the PDCSAP curves appeared to be better (fewer long-term trends).

For the cases that were not preselected, individual light curves were derived from image cutouts of 51×51 pixels using aperture photometry performed with the Python package Lightkurve⁶. Only data with a high quality were selected (quality_bitmask=‘hard’ in task SEARCH.TESSCUT). The source mask was defined from pixels with fluxes above a few times the median absolute deviation over the median flux (the exact threshold value was adapted for each target, depending on crowding). The background mask was defined by pixels with fluxes below the median flux (i.e., below the null threshold). To estimate the background, we used two methods. First, a principal component analysis, with five components, and second, a simple median. In each case, the two background-subtracted light curves were compared. Most of the time, the light curves that were corrected with the principal component analysis were similar to or better than the median-corrected ones (fewer long-term trends) and were therefore kept.

The TESS fluxes were then converted into magnitudes, and the average magnitude was subtracted from each individual light curve. A sample of individual light curves are shown in Fig. B.1 and figures showing all light curves are available at Zenodo⁷. As the TESS pixels are large (21″), the time series might be exposed to crowding issues. We examined the close (within 1′) environment of our targets in the Gaia DR3 catalog. There are no bright (∆G < 2.5 mag) neighbors in that area, except for HD 154407C⁸ and the additional target, V1208Sco. The spectra of HD 154407C (see below) were used to verify that the RV variations were coherent with the photometric variations, and because they were, the target was kept. To ensure that the analyzed curves suffer as little as possible from crowding, we relied on the PDCSAP curves because they are supposed to be corrected for these effects. Three targets lack 2 min data, however: HD350685, LSI+61 275, and V1208 Sco. For these targets, we extracted the TESS light curves with the TESS-GAIA light curve (TGLC) package⁹ (Han & Brandt 2023). This package uses neighboring sources in the Gaia catalog, assuming they are constant. We choose the aperture extraction¹⁰ to derive the cleaned photometry.

Kepler data are also available from the MAST portal for TYC 1881-933-1 (see ‘K’ in the last column of Table A.1). As for TESS data, we kept only the high-quality data points of the PDCSAP curve and then converted its fluxes into magnitudes. A long-term trend was fit by a straight line and then taken out, so that the Kepler and TESS curves display the same base level.

2.3 Spectroscopy

While eclipses with reflection provide much information on the systems, complementary spectroscopy is required to assess the true nature of the systems (pre/post-interaction). We therefore gathered spectroscopic data for all the systems.

For three systems, the velocities are available in the literature. LS I +61 275 was studied by the Georgia state team (Huang & Gies 2006a,b) and by the Apache Point Observatory Galactic Evolution Experiment (APOGEE, Jönsson et al. 2020). The former references provided effective temperature, projected rotational velocity, and two radial velocity (RV) measurements. The latter reference provides 8 RV points. While both indicate that the RVs decrease around the main eclipse, they suggest very different RV amplitudes: about 100 km s⁻¹ for the former ones, and four times lower for the latter (Fig. B.2). This strong discrepancy on the RV amplitudes, coupled with the impossibility to derive a good photometric fit when using a mass function derived from the more numerous APOGEE data, prevented us from performing a detailed analysis of this star before a new high-quality spectroscopic campaign is performed. TYC 1881933-1 was studied in the framework of the Large Sky Area Multi-Object Fibre Spectroscopic Telescope (LAMOST) survey. Corrected RVs, determined from the blue and red channels, were published by Zhang et al. (2022). The blue-channel RVs behave similarly to the red channel RVs, but have more scatter. We therefore only kept the latter to calculate an orbital solution. One value, at BJD = 2 458 122.20438, appears to be discrepant and was therefore taken out, which left 43 RV points. Effective temperature, gravity, and projected rotational velocity were estimated for each LAMOST spectrum by Xiang et al. (2022). The average values and their dispersions are reported in Table A.1. Finally, the velocities for the primary of V1208 Sco were provided by Banyard et al. (2022), and its effective temperature, gravity, and projected rotational velocity were evaluated by McSwain et al. (2009). These literature velocities are recalled in the appendix (Table B.1).

We collected spectroscopy for the other targets using several dedicated campaigns. Four spectra were collected for HD 350685 with the Hermes echelle spectrograph installed on the Mercator telescope (La Palma Observatory). In the framework of the Optical Infrared Coordination Network (OPTICON) program, spectra were collected for us with the SOPHIE¹¹ echelle spectrograph at the Haute-Provence Observatory for HD 254346 (four spectra; program 23B004, PI Nazé) and with the CARMENES¹² echelle spectrograph at the Calar-Alto Observatory for LSI+61 145 and TYC741-1565-1 (two spectra each; program 24B010, PI Nazé). Finally, between five to seven spectra were collected at the European Southern Observatory (ESO) with the XShooter spectrograph for AN Dor, TYC741-1565-1, HD112485, V1061Cen, and HD154407C under run ID=114.26ZK.001 (PI: Nazé). These spectra were complemented by archival spectra available in the ESO archives¹³ for HD 309317 (eight FEROS¹⁴ spectra) and ANDor (six XShooter and three FEROS spectra). The XShooter spectra were taken with a medium spectral resolution (R between 5500 and 10000), and all other data relied on high-resolution spectrographs (R between 40 000 and 95 000). All spectra were provided in reduced form by the observatories. The CARMENES spectra were calibrated for vacuum wavelengths, however, while the other data used an air wavelength calibration. The CARMENES wavelengths were therefore converted using λ_air = λ_vac/1.000272. Cosmic rays were eliminated using a median filter of size 2 (MIDAS FILTER/MEDIAN command) for all but XShooter and CARMENES spectra (which had very short exposures and hence no cosmic-ray hits).

3 Results

3.1 Light-curve analyses

First, we performed a period search on all individual light curves up to 25 d⁻¹ to assess their frequency content (Heck et al. 1985). The orbital signal clearly dominated in all cases. All TESS curves of the same target were then merged, and the period search was repeated around the main signal, with a smaller frequency step (adapted to the longer time interval covered by the data). This allowed us to obtain preliminary values for the orbital period and reference time of the main deeper eclipse, which were then refined by examining the folded light curves. The final best-fit ephemerides can be found in Table A.2. An average phase-folded light curve (with 200 phase bins for TESS data and 100 for Kepler data) was then calculated from all crowding-corrected curves, except when there was too much noise or some obvious change between individual light curves (the data we used are identified by boldface in Table A.1). These average curves were used to remove the binary signal from the data, and the residuals were examined in search of pulsations. The appendix provides the results of this search.

The average binned light curves were fit using the model by Moe & Di Stefano (2015). This analytical model uses two Gaussians to represent the eclipses and a cosine function to reproduce the reflection effect: $⟨I⟩+\mathrm{\Delta }{I}_1[{\exp }^{-{\varphi }^2/2{\theta }_1^2)}+{\exp }^{(\varphi -1{)}^2/2{\theta }_1^2}]+\mathrm{\Delta }{I}_2{\exp }^{-(\varphi -{\varphi }_2{)}^2/2{\theta }_2^2}+0.5\mathrm{\Delta }{I}_{refl}\ast [\cos (2\pi (\varphi -1/2))+1$, where the numbers i = 1 and 2 correspond to the deeper and secondary eclipses, respectively, ∆I_i and θ_i are the eclipse depths and widths, respectively, ∆I_refl is the amplitude of the reflection, φ₂ is the phase of the second eclipse, and φ is the orbital phase derived from the average ephemeris (see Cols. 2 and 3 of Table A.2). The best-fit parameters, found through a LevenbergMarquardt algorithm, are listed for each target in Table A.2. Despite its simplicity, this model is very useful to catch the main features of the light curves: the eclipse depth, the reflection strength, and the potential eccentricity. For example, the phase φ₂ of the shallower eclipse is compatible with 0.5 for all but three cases: ANDor, HD254346, and V1061 Cen. This indicates that these three systems have eccentric orbits. Their lowest eccentricities, e_min = π|φ₂ - 0.5|/2 (Moe & Di Stefano 2015), are 0.049, 0.010, and 0.008, respectively (and the typical values should be about 1.6 times e_min for a uniform distribution of the periastron angle ω).

The model of Moe & Di Stefano (2015) was then used to precisely derive the individual eclipse times. For each orbital cycle covered by the data, the times of the main and secondary eclipses were predicted from the average ephemerides and the phase of secondary eclipse (Table A.2). Then, the eclipse models were shifted within the range ±2 × θ_i (with i = 1,2 for the main and secondary eclipses, respectively), and the χ² between light curve and model was evaluated for each shift within a phase interval of ±3 × θ_i around each eclipse. We then examined not only the differences between the best-fit eclipse times and their predictions but also the interval between primary and secondary eclipses of the same cycle. This allowed us to detect apsidal motion in AND or and global eclipse shifts in TYC 1881-9331, HD 309317, and V1061 Cen. Dedicated analyses of these stars are provided in Sections 3.4 and 3.5, respectively.

3.2 Spectral analyses

Two independent derivations of the RVs were made for the spectra. In the first derivation, normalization was achieved by spline fitting through dedicated spectral windows, and then RVs of Hβ and He Iλ5876Å were evaluated independently by crosscorrelation with a 20kK model spectrum from the BSTAR2006 database (Lanz & Hubeny 2007). Finally, an average of the two RVs was calculated. In the second derivation, normalizations were made for the whole wavelength range using the tool called spectrum normalization neural network (SUPPnet) by RóZański et al. (2022)¹⁵. The RVs were then derived using the iSpec tool (Blanco-Cuaresma 2019)¹⁶. To do this, iSpec cross-correlated the observed spectra with a generic TLUSTY model for an earlyB star (T_eff = 20kK, log(g) = 4, and v sin(i) = 150 km s⁻¹), taken from the BSTAR2006 database¹⁷. The wavelength ranges were 3650-5600 Å (XShooter UV subspectra), 3800-6800 Å (FEROS and Hermes spectra), 3900-6800Å (SOPHIE spectra), and 5700-6800 Å (CARMENES spectra). The derived RVs, provided in Table B.2, agreed well with those derived from the previous method. Finally, we also made a trial with templates of lower temperatures (8 kK and 14 kK, which bracket the derived secondary temperatures; see below), but to no avail. The results for the 14 kK template were similar to those of the 20 kK template, and the 8 kK template resulted in velocities with extremely large errors and/or that were constant (or affected by the primary).

The individual spectra of each target were shifted by minus the determined RVs to place them in the same zero velocity reference frame, and they were then averaged. A fitting of these average spectra was then performed within iSpec to determine the closest BSTAR2006 TLUSTY models and thereby obtain the best physical parameters of the (primary) star. This was done first with the metallicity set free to vary and over the whole wavelength range. The fitting was then repeated with the metallicity fixed to solar and considering only the ranges 4000-5100 Å, 5860-5885 Å, and 6540-6700 Å, in which the most interesting and most numerous spectral lines appear. The best-fit values for the effective temperature T_eff, gravity log(g), and projected rotational velocity v sin(i) are provided in Table A.1, except for LS I+61 145, for which the CARMENES red spectra contained too few stellar lines for a meaningful estimate of the effective temperature and gravity. The best-fit temperatures were used, as mentioned above, to improve the evaluations of the spectral types and bolometric corrections. Finally, the bolometric luminosities, effective temperatures, gravities, and projected rotational velocities were entered in BONN Stellar Astrophysics Interface (BonnSAI)¹⁸ (Schneider et al. 2014) to derive estimates of the evolutionary masses, radii, and ages of the massive primaries (Table A.1).

At first sight, the spectra only showed one stellar signature. To assess whether a weak second component existed, we shifted the spectra with minimum and maximum RVs to the zerovelocity reference frame. We plotted the superimposed spectral lines and also calculated their difference. Except for the case of V1061 Cen (which is further analyzed in Sect. 3.6), the line profile appeared to be constant, without any weak additional component. The targets thus appear to be of SB1 type in our dataset. Most probably, more spectra of the highest quality are needed to detect the secondary signature, as for HD 25631 (Nazé et al. 2023). Tkachenko et al. (2024) concluded from the same FEROS dataset as ours, however, that HD 309317 was SB2. Figure 1 compares the line profiles of AN Dor, HD 309317, and V1061 Cen after the spectra were shifted in a zero-velocity reference frame. Only the latter target clearly shows the signature of a second component in the wings of the main line. Other lines and other phases were also checked, but we were still unable to detect any signature typical of a SB2 nature for HD 309317 (e.g., a second component or line doubling, a shallower and broader line at quadratures or a deeper and narrower line at conjunctions). In this context, we also note that the K₁ and K₂ mentioned by Tkachenko et al. (2024) would lead to a mass ratio of 0.44, implying a late-B star nature for the secondary, whose lines should then be rather obvious in the data (but even a disentangling trial failed to unveil a meaningful secondary spectrum). Therefore, within the limitation of current data, the SB2 classification must be revised to SB1.

Fig. 1 Line profiles of He I λ4388 Å in the most blue- and redshifted spectra (blue and red lines, respectively) for AN Dor, HD 309317, and V1061 Cen. The spectra were corrected for the measured RVs, and the black line provides their difference (multiplied by four and then shifted by +1 to facilitate the comparison with the lines).

3.3 System characterization

After applying the methods outlined in Sections 3.1 and 3.2, we present the folded light curves and RVs in Fig. B.2 for all targets. One common characteristic then appears obvious. Around the time of the deepest eclipse, the RV of the massive star is found to decrease. This implies that the deepest eclipse corresponds to the eclipse of the massive star by its companion (i.e., a conjunction with the massive star behind its companion), and this massive star is the hotter component of the system. Our targets thus possess the same nature. They are young pre-interacting systems with fast-rotating B-type primaries and very low-mass secondaries.

With light curves and RVs at hand, the physical parameters of the systems were constrained in detail. For V1208 Sco, a SB2 fitting was performed by Banyard et al. (2023). The velocity amplitudes amount to K₁ = 44.7 ± 0.4 kms⁻¹ and K₂ = 259.6 ± 0.8 kms⁻¹ (L. Mahy, private communication), which leads to a mass ratio of q = 0.172 ± 0.002 and to mass estimates of M₁ sin³(i) = 12.99 ± 0.11 M_⊙ and M₂ sin³(i) = 2.24 ± 0.03 M_⊙. For the other targets, we started by fitting the RVs with sinusoids, except for AN Dor (see the subsection below). Only three systems were found to be eccentric, and two of them (HD 254346 and V1061 Cen) have a low e_min, together with only a few spectra. In these conditions, the eccentricity cannot be constrained precisely enough, and a noneccentric solution provides a sufficiently good fit for our purposes. Table A.3 yields the best-fit parameters, and we note that our RV amplitude of HD 350685 agrees well with that published by Tkachenko et al. (2024). Finally, the RV solution for LS I +61 145 should be considered preliminary because only two RV points are available.

From the RV amplitude and the period, we calculated the mass function. Using the primary masses from Table A.1, we derived first estimates of the secondary mass and the mass ratio q for an inclination i = 70-90°. The average light curves were then fit using Nightfall v2.0.2 (Wichmann 2011), assuming I band for the data (the closest to TESS bandpass), detailed reflection, and a quadratic limb darkening. The eccentricity, mass ratio, plus the mass and effective temperature of the primary star were fixed to the values appropriate for the target under consideration. A preliminary fitting provided a first inclination value, which was iteratively used, in conjunction with the mass function, to refine the mass ratio. The final best-fit parameters are provided in Table A.3. One caveat must be mentioned. The errors in this table correspond to formal errors around the provided solution. Exploration of the parameter space reveals that solutions of similar quality but with different parameters can be found. As in Nazé et al. (2023), the true uncertainty is thus larger than that around a specific individual solution. The realistic error values are a few degrees for inclinations, about 10% on radii, and several kK for the temperature of the secondaries.

We used the effective temperatures and radii to evaluate the bolometric luminosities. In all cases, the secondary appears to be much fainter because their luminosities are <2% of the primary luminosities. Furthermore, the derived total luminosities agree well with those derived from the V magnitude, distance, and extinction (Table A.1). The radii of the primaries also agree with those found from evolutionary tracks (Table A.1). Finally, the gravities of the primaries were computed from the (assumed) mass and the derived radius. A good agreement was again found with those estimated from the spectral fitting by atmosphere models (Table A.1). With values around log(g) = 4, the gravities also confirm the main-sequence status of the primaries.

Fig. 2 Light curve observed with ASAS (four epochs; see text) and TESS for AN Dor (three sectors shown). The separation between the primary and secondary eclipses evolves.

3.4 Apsidal motion in AN Dor

The intervals between the main and secondary eclipses remain compatible with a constant for all targets but one: AN Dor. Increasing shifts between predicted and observed eclipse times were also observed for this system. All this is typical of apsidal motion (Giménez & Bastero 1995; Claret et al. 2021). To obtain a longer temporal baseline, that is, to improve the characterization of this effect, we also downloaded All-Sky Automated Survey (ASAS ¹⁹, Pojmanski 1997) data of this system and kept only data points with the highest quality (flag=‘A’). We then cut the ASAS curve into four parts (HJD-2 450 000 < 2350, 2350-3200, 3200-4000, and > 4000). For each interval, we first estimated the shift between the folded ASAS data and the average binned TESS curve within ∆(0) = ±0.15 of the main eclipse. After we determined the best-fit time of the main eclipse, the Moe & Di Stefano model was fit to the folded ASAS data to derive the phase (hence, the time) of the secondary eclipse. The evolving light curve of AN Dor is shown in Fig. 2, and the best-fit model parameters are provided at the bottom of Table A.2, except for the first interval, which lacks enough data points to provide meaningful results.

Newtonian apsidal motion is a secular change in the longitude of the periastron ω that is caused by the nonspherical shape of the gravitational potential in an eccentric close binary (e.g., Giménez & Bastero 1995; Schmitt et al. 2016; Rosu et al. 2020, and references therein). Giménez & Bastero (1995) provided the mathematical expressions of the times of photometric minima as a function of eⁿ cos n ω and eⁿ sin n ω, where e is the eccentricity, and n takes integer values between 1 and 6. These quantities are multiplied by functions of cot i and csc i, where i is the orbital inclination.

We used 469 times of minima inferred from the TESS data of AN Dor along with 6 times of minima derived from the ASAS observations to search for apsidal motion. For this purpose, we adopted relation (15) of Giménez & Bastero (1995) to compute the theoretical times of minima and compared them to the observed values of the difference with respect to best-fit linear ephemerides (O - C_linear). In addition to the orbital inclination i, the theoretical expression of O - C depends on four parameters: the anomalistic orbital period P^an, the eccentricity e, the rate of apsidal motion ω, and the longitude of periastron ω₀ at a reference time. For values of i between 70° and 90°, by steps of 5°, we built a four-dimensional grid to sample the parameter space. For each combination of parameters of our grid, we adjusted the linear ephemerides in such a way as to achieve a zero mean residual. When we identified a combination of parameters that yielded the lowest χ² residuals between the theoretical and observed O - C_linear values, we refined the grid to achieve a higher resolution. As a general result, we found that the dependence of the solutions on i is very weak in the range of inclinations we considered. Only the best-fit e and ω₀ parameters changed slightly as a function of i. The largest differences were found in the best-fit eccentricity, which increased from 0.0496 to 0.0528 as i changed from 70° to 90°. Table 1 lists the best-fit parameters and their 1σ error bars for i = 80°, and the corresponding best-fit adjustment of time series of O - C_linear values is displayed in Fig. 3.

Given the high number of TESS data points and their excellent accuracy, our best fit is clearly dominated by these points. The inclusion of the ASAS data still provided a valuable confirmation of the inferred parameters. The times of primary eclipses derived from the ASAS data agree very well with our best solution. The secondary eclipses are subject to much higher uncertainties in the ASAS observations, but even there, the trend is consistent with our solution. From what we find here, the rate of apsidal motion corresponds to a period U = (55.4 ± 4.3) yr.

Keeping the eccentricity, rate of apsidal motion and anomalistic orbital period fixed to their best-fit values derived from our analysis of the photometric data, we then adjusted a RV curve that explicitly accounted for apsidal motion (see Eqs. (1) and (2) of Rauw et al. 2016) to determine the time of periastron passage and the amplitude of the SB1 RV curve. The results are given in Table 1 and illustrated in Fig. 4.

Because of the apsidal motion, we did not combine the whole set of TESS light curves. We instead averaged only a subset of curves from Sectors 63, 64, and 65, which cover HJD= 2 460 014 - 97. At the time (Table 1), the periastron argument should be around ω = 340.76°. We therefore calculated a photometric solution as for other stars, but fixing the periastron argument to this value, e = 0.052 (Table 1), $T_{eff}^1=21529$ K, M₁ = 8.2M_⊙ (Table A.1), and q = 0.205 (found from the mass function and the inclination of a preliminary fitting run). The best-fit parameters are provided in Table A.3.

Fig. 3 Best-fit adjustment to the differences between the observed times of minima of AN Dor and the linear ephemerides O − Clinear for i = 80°. The orange and pink symbols show the eclipses by the primary and secondary stars, respectively. Data prior to 2010 are taken from the ASAS survey, and data after 2018 are taken from TESS photometry. The long and short dashed lines yield our best-fit theoretical relations for the times of minima. In this plot, the linear ephemerides are T0 = 2 460 204.6031 and Psd = 2.0326782 d.

Fig. 4 Comparison between the measured RVs of AN Dor (blue symbols with error bars) and the best-fit orbital solution explicitly accounting for apsidal motion.

3.5 Eclipse time shifts

Even when the interval between the two eclipses remained stable, evolving shifts of the light curves can be detected. This was so for three cases: TYC 1881-933-1, HD 309317, and V1061 Cen (Fig. 5). The most likely explanation for this change is a lighttravel time effect (LTTE). In a triple system, the motion of the inner eclipsing binary about the center of mass of the triple system leads to a periodic change in the eclipse times because the velocity of light is finite. This is expressed as $$(O-C)(t)=\frac{a_{bin}\sin i}{c}\frac{1-e^2}{1+e\cos \psi (t)}\sin (\psi (t)+\omega )$$(1)

where a_bin, e, i, and ω stand for the semimajor axis, eccentricity, inclination, and argument of periastron of the orbit of the inner binary around the center of mass of the triple system. For each epoch t, the true anomaly ψ(t) of the inner binary was computed solving Kepler’s equation, which involves e, the time of periastron passage $T_0^p$, and the orbital period P_out of the outer orbit. The comparison of the observed O - C values with Eq. (1) then allowed us to constrain several parameters of the outer orbit. These are $\frac{a_{bin}\sin i}{c}$, e, ω, $T_0^p$, and P_out.

For each of the three targets, three epochs of photometric observations exist. For HD 309317 and V1061 Cen, the data span a total duration of 1500 d, and the differences are mostly detected between the second epoch and the other two (which have similar O - C). The O - C changes in that second epoch by about 0.02 d and 0.01 d, respectively. With basically only three points available, the outer orbit is clearly not well constrained. Because the inner binary details remain unknown for V1061 Cen (see Section 3.6), we only provide a preliminary estimate for HD 309317. For the latter system, the orbital solution suggests the binary mass to be about 10 M_⊙ (Table A.3). We further assumed a zero eccentricity for the outer orbit and inclinations between 30° and 90°. We also considered two (simple) possibilities: either the observations show half a cycle or a full cycle, that is, $\frac{a_{bin}\sin i}{c}=0.02$ d and P ~ 3000 d or $\frac{a_{bin}\sin i}{c}=0.01$ d and P ~ 1500 d. Kepler’s law indicates $(a_{bin}\sin (i){)}^3={\left(\frac{P}{2\pi }\right)}^2\frac{G{M}_3^3\sin (i{)}^3}{(M_{bin}+M_3{)}^2}$. Therefore, M₃ is found to be 4-14 M_⊙. While RVs measured on spectra taken over only 60 d should not be affected by the outer orbit, the spectral lines from a star at the higher end of this mass interval should have been detected in the spectra of HD 309317, but they only show a single component (Fig. 1). The coverage of the outer orbit is scarce, however. More observations and therefore a more precise derivation of the outer orbit are now needed to reach a better knowledge of the tertiary component.

Finally, TYC 1881-933-1 displays clear O - C trends during each epoch, not only between observing campaigns (see the insets in the top left panel of Fig. 5). More precise constraints can then be derived in this case because the O - C needs to fold in a coherent way if it is due to LTTE. First, we determined the period of the LTTE effect in this target. The Fourier periodogram of the O - C values of the times of primary minimum of TYC 1881-933-1 displays a large number of peaks at frequencies below 0.01 d⁻¹ (see the right panels of Fig. 5). The frequencies of most of these peaks are integer multiples of the lowest frequency (3.7 10⁻⁴ d⁻¹, corresponding to 2703 d). Others are likely aliases. The complex spectral window has strong peaks (amplitude higher than 0.8) at 1.13 10⁻³ d⁻¹,1.45 10⁻³ d⁻¹, 2.58 10⁻³ d⁻¹ and 4.03 10⁻³ d⁻¹, followed by somewhat lower peaks (amplitude above 0.65) at 2.90 10⁻³ d⁻¹ and 3.2 10⁻⁴ d⁻¹.

We then adjusted the light-travel time relation of Eq. (1) with the periods corresponding to each of these peaks and retained the periods that yielded the best adjustments to the data. They are 169.2, 225.2, 337.8 and 675.2 d. The formally best-quality fit was obtained for the longest of these periods, although the rms of the fit is nearly independent of the period. The three shorter periods are the n^th harmonics of the 675 d period, with n = 2, 3 and 4. Given the nonsinusoidal shape of the light-travel time relation, the 675 d period thus appears to be the most promising candidate. The parameters for the four possible periods are given in Table 2 and are displayed in Fig. 6. Based on the periods and delays of that table and a binary mass of 7.3 M_⊙ (Table A.3), Kepler’s law indicates a tertiary mass of M3 = 11-44 M_⊙ for i = 30-90°. The RV changes due to the gravity of the tertiary component would then amount to $K=a_{bin}\sin (i)\frac{2\pi }{P\sqrt{1-e^2}}\sim 60\mathrm{k}\mathrm{m}\mathrm{s}$. There are two problems with these evaluations. First, lines from a star that is this massive would be detectable in the spectra, which is not the case (only lines from the B2 primary are seen). Second, these strong RV changes would be seen while the RV curve remains similar, regardless of the epoch under consideration (LAMOST RVs were taken over 800 d with eight campaigns with durations of 1-8 d; see Fig. B.2 and Table B.2; and they sample phases near apastron and near periastron). Again, more data are needed to clarify this result.

Fig. 5 Left : differences between observed and predicted eclipse times based on the ephemerides of Table A.2. The data for the deeper eclipse are shown by blue dots, and the data for the shallower eclipse are shown as red circles. Close-ups are provided for TYC 1881-933-1. Right : Fourier spectrum of the O − C values of the times of primary minimum of TYC 1881-933-1 (top panel) and the associated spectral window (bottom panel).

3.6 Spectral disentangling of V1061 Cen

For V1061 Cen, the main spectral component appears to be stationary, while the fainter secondary component seems to shift with orbital phase. Fitting two Gaussians to the line profiles (one for the main component, the other for the secondary component) did not work well. Since lines appear in absorption, a first estimate of the spectra of the main component was built by taking the flux maximum of the five spectra at each wavelength step. This average spectrum was subtracted from all original spectra to isolate the secondary component, and a single Gaussian was then fit to the lines that appeared in this difference. Using these first estimates of the RVs, we performed a spectral disentangling on the whole spectral set using the shift-and-add procedure of González & Levato (2006), as coded by Rauw & Nazé (2016). This was done for three domains, 4000-4300 Å, 4300-4600 Å, and 4600-4950 Å. After fewer than 200 iterations, the results appeared to be stable. The final RVs, averaging the results from each of the three domains, are provided in Table B.2 and shown in Fig. B.2. The spectral appearance (Fig. 7) suggests spectral types of B0V and B2.5V for the main and secondary components, respectively. As the main component is stationary, V1061 Cen clearly is a triple system, but only one component of the eclipsing binary is detected. We attempted a disentangling with three components, but found no significant trace of the missing binary component. With an amplitude of 209.6 ± 5.2kms⁻¹ and a period of 2.21 d (Table A.2), the mass function reaches f (m) = 2.11 ± 0.16 M_⊙. Assuming the detected star has 6 M_⊙ and considering i ~ 90° because of the presence of eclipses, the mass of the other component should be around 7 M_⊙, implying a star more massive than the component detected in the spectra. Since the RV of the detected component decreases around the main eclipse (Fig. B.2), the other component should also be cooler. Without a trace of this component in the spectra, however, we were prevented from deriving a full solution for this system. A more thorough investigation is now needed that probes several epochs with a better phase coverage.

Fig. 6 Light-travel time effect adjustment of the O − C values of the times of primary minimum of the inner eclipsing binary of TYC 1881-933-1. Phase 0.0 corresponds to the periastron passage of the outer orbit.

Fig. 7 Result of the disentangling for V1061 Cen. The main component is shown in blue, and the fainter component is shown in red (continuum shifted by +0.2 for clarity). No adjustment for the light ratio was made. The line depths are those observed in the input spectra, i.e., with respect to the combined light of all stars in the system.

4 Discussion and conclusions

Our combined analysis of photometric and spectroscopic data adds a dozen systems in which a massive star is paired with a luminous companion of much lower mass. The system characteristics can be found in Tables A.1 and A.3. In all cases, this companion is cooler than the massive star, and it is therefore not expected to be a stripped star resulting from a binary interaction.

The location of our targets in Hertzsprung-Russell diagrams and mass-radii diagrams are provided in Fig. B.3. An additional plot provides a comparison between masses and ages for the secondaries alone because all primaries appear to be beyond the zero-age main sequence (ZAMS; Fig. 8). The mass ratios of AN Dor and HD 112485 are slightly higher than q = 0.2, which means that they do not fulfill the criterion of an extreme mass ratio. In addition, considering their ages, their secondaries might be young MS stars with early-F and late-B spectral types, respectively, rather than PMS stars. The other secondaries, however, are clearly less massive and still in a PMS stage.

In Figures B.3 and 8, our targets (shown by red dots) are compared to other young systems with an early-type primary (earlier than B3), a short period (P < 10 d), and an extreme mass ratio (q < 0.2). This set comprises 16 of the nascent systems from Moe & Di Stefano (2015) in the Magellanic Clouds, shown by black circles, which were analyzed on the basis of photometry alone (spectroscopy is still required to confirm their status). In addition, this set also contains several Galactic cases whose analysis usually combined spectroscopy and photometry (Table B.3). Some of them (shown by open orange triangles in the figure) may not be true nascent binaries (i.e., systems with PMS secondaries), but rather young binaries with stars very close to the ZAMS: HD 163892 (Mahy et al. 2022), HD 165246 (Johnston et al. 2021), and HD210478 (V446Cep, Çakirli et al. 2014). Other cases, shown by filled orange triangles, clearly harbor PMS secondaries, as demonstrated by their position in the age-mass diagram: HD 145482 (c02 Sco, Gullikson et al. 2016; Pigulski 2024), HD 149834 (Stassun et al. 2021), HD 120307, HD 138690A (νCen and γLupA, Jerzykiewicz et al. 2021), HD25631, HD 191495, and HD46485 (Nazé et al. 2023). Three additional systems (shown by orange stars in the age-mass diagram) most probably also belong to this group, although only part of their physical parameters have been derived, and more study is required to fully constrain the full set of their properties: HD158926A (λ ScoA, Tango et al. 2006), HD 152200 (Banyard et al. 2023; Britavskiy et al. 2023), and HD62747 (V390Pup, Pigulski 2024).

Binary systems, especially nascent ones, allow us to peer directly into the results of the stellar formation processes. As an example, Smith et al. (2024) used the Gaia DR3 catalog to study close (d < 500 pc) astrometric binaries with periods of 1003000 d and intermediate-mass primaries (spectral types between B8 and F1). They found that systems with the lower mass ratios (q < 0.32) and/or lower eccentricities (e < 0.2) have primaries that rotate more slowly and display a higher spin-orbit alignment than primaries in other systems. This excess of slow rotators supports the disk-fragmentation scenario, in which a very low-mass fragment is formed from a disk instability, then migrates inward, and finally accretes from the (now circumbinary) disk. In this case, higher accretion rates lead to more similar stars in more circular orbits, so that modest accretion is required to keep the mass ratio at low values. The lower eccentricity would then need to naturally arise from the fragmentation+migration process (Smith et al. 2024). The slow rotation of the primaries is then explained by two effects: (1) the circumbinary accretion mostly occurs on the secondary, so that after fragmentation, the primary star gains less in mass and angular momentum; (2) magnetic braking may transfer angular momentum from the primary through the (circumprimary) disk into the orbit.

For systems with shorter periods, as studied here, the situation is expected to be quite different. The proximity of the two stars leads to a truncation of the circumprimary disk, prohibiting magnetic braking (Smith et al. 2024 and references therein). Furthermore, because of their very low masses, the secondaries exert little tidal effects on the massive primaries (Britavskiy et al. 2024), which limits braking of the rotation and therefore synchronization onto the orbital period. All this implies that stars born rotating fast will remain so. In this context, we compared the orbital period P to the rotation period P_rot of the primary derived from the projected rotational velocity, the radius, and the binary inclination. This revealed that only TYC 1881-933-1 and LS I+61 145 are synchronized. The rotation of other primary stars is clearly supersynchronous (Table A.3): P_rot/P ∼ 0.3 for HD254346, HD112485, and HD350685 or 0.5-0.8 for AN Dor, TYC 741-1565-1, HD 309317, HD 154407C, and V1208 Sco. This fast rotation²⁰ was previously found for the three systems that we analyzed in Nazé et al. (2023), and this also applies to the other nascent systems reported in the literature, except for HD 120307 and HD62747 (Table B.3 and Fig. B.4).

Supersynchronous rotation in nascent systems may appear to be unsurprising, but the high incidence of these cases in this group is clearly not expected. For example, Huang et al. (2010) studied rotation in B stars and demonstrated that the higher-mass (M > 8 M_⊙) cases rotate more slowly, and very fast rotation (v/v_crit > 0.7, see their Fig. 6) is even totally absent for them. One may wonder whether the observed large fraction of super-synchronous rotators among the nascent short-period binaries with extreme mass ratios is an observational bias. A priori, stellar rotation plays no role in the appearance of eclipses and reflection effects that were used to detect them. Observations may thus reflect the actual situation. In support of this idea, we may note that the vast majority of young HAeBe stars analyzed by Pinzón et al. (2021) are not slow rotators, but display v sin(i) > 100 km s⁻¹. It is true that the sample of nascent systems currently remains limited, however, with only 18 cases. Increasing the sample would certainly improve the reliability of these conclusions and clarify the rotation distribution. In this context, we may note that disk fragmentation modeling and observations revealed the largest fraction of systems with very low mass ratios to lie in long-period (P > 50 d) binaries (Tokovinin & Moe 2020). While long periods make the spectroscopic detection of the orbital motion of the primaries even more difficult than for the cases studied here, the study of these systems would also help us to reach a broader view of this population, and it would thereby further constrain the star formation scenarios.

Fig. 8 Comparison of mass and age for secondary stars of extreme mass ratio systems. Our targets are shown by red symbols, the MC systems of Moe & Di Stefano (2015) are shown with black circles, and the other Galactic cases are shown with orange triangles or stars (see text for details). The ZAMS track was defined from models by Tognelli et al. (2011), considering a central hydrogen abundance Xc equal to 99.5% of the initial value (which corresponds to the stabilization of the total nuclear-based luminosity in these models).

Data availability

Figures showing all light curves are available on Zenodo (https://doi.org/10.5281/zenodo.17257979).

Acknowledgements

The Liège authors acknowledge support from the Fonds National de la Recherche Scientifique (Belgium) and the University of Liège (through its Special Funds for Research and the IPD-STEMA program). This paper includes data collected by the TESS mission, which are publicly available from the Mikulski Archive for Space Telescopes (MAST). Funding for the TESS mission is provided by NASA’s Science Mission directorate. ADS and CDS were used for preparing this document. This project has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No 101004719 (OPTICON/ORP) and J.L.-B. was co-funded by the European Union (ERC, MAGNIFY, Project 101126182) - this material reflects only the authors views and the Commission is not liable for any use that may be made of the information contained therein.

References

Appendix A Tables presenting the target and their properties

Appendix B Additional information

A sample of individual light curves of each target is shown in Fig. B.1 (the full set is available on Zenodo), while folded ones are shown in Fig. B.2, along with the folded RV curves. Literature velocities are recalled in Table B.1 and our newly determined RVs are provided in Table B.2. Table B.3 lists the whole sample of nascent, short-period massive binaries with extreme mass ratios and Fig. B.4 provides histograms of the rotational velocities, highlighting their high values.

Finally, the photometric data were cleaned from the binary signals using the average folded and binned light curves. A period search was then performed on the residuals. A forest of low-frequency signals with increasing amplitudes (“red noise”) are obvious in these frequency spectra, notably for LS I +61 145, TYC741-1565-1, HD 154407C, and V1208 Sco. While we cannot fully exclude the presence of imperfect calibration effects and some small long-term variations in these curves, stochastic low-frequency signals are well known to be present in massive stars and our targets would thus simply be no exception. In addition, isolated signals are also present. Those with low frequencies (ν < 4 d⁻¹) usually indicate Slowly Pulsating B-stars (SPB) activity. In our sample, low-frequency signals are detected for AN Dor, HD254346, HD309317, HD112485, and V1061Cen. In contrast, isolated signals at high frequencies (ν > 4 d⁻¹) usually indicate the presence of β Cep pulsations. Eze & Handler (2024) already reported the detection of β Cep activity in LS I+61 145, HD 112485, HD254346, and V1061 Cen, plus ANDor (see also Southworth & Bowman 2022) and HD 154407C as candidates. We confirm their detections and add to the list the cases of HD309317, HD350685, and V1208Sco, plus TYC 1881-933-1 as candidate. The cases of HD 154407C and V1208 Sco must of course be taken with caution as the TESS data alone cannot disentangle the source of pulsations which could thus be these stars or their close neighbors. Nevertheless, the presence of both eclipsing and pulsational signals provides a golden opportunity for asteroseismology. The orbital analysis helps pinpoint the stellar properties, which can then serve as input for a detailed asteroseismic modeling. The fact that most of our sources are rotating quite fast constitutes both a challenge and an asset in this domain, as the pulsations of such stars remain difficult to model.

Fig. B.1 Individual light curves of the targets - sample. The full set is provided as an additional pdf on Zenodo.

Fig. B.2 Binned and averaged light curves of all targets (bottom panels - the data used to produce these curves are identified by boldface in Table A.1; residuals are shown if a fitting was achieved), together with RV curves (top panels). The folding was done using the ephemerides of Table A.2. The best-fit solutions from Table A.3 are superimposed to the data using red lines. For LS I +61 275, the black RV points come from APOGEE while the red ones come from Huang & Gies (2006a), showing the discrepancy between them. For V1061 Cen, the black filled dots and blue circles correspond to the main and fainter spectral components, respectively (see Table B.2 and text for details).

Fig. B.3 Hertzsprung-Russell diagrams for primaries (top left, with the post-PMS evolutionary tracks of Brott et al. (2011) without initial rotation in blue and with the fastest initial rotation available in these models in violet lines) and secondaries (top right, with the PMS tracks of Tognelli et al. 2011 in green). Isochrones are shown in lighter color and with dashed lines. Symbols as in Fig. 8. For our targets, the secondary effective temperatures are those without reflection and the bolometric luminosities were calculated from these temperatures and the stellar radii (Table A.3) while the bolometric luminosities of the primaries are those from Table A.1. Bottom panels provide the mass-radius diagrams, using the same conventions. The comparison with theoretical evolutionary tracks is not without limitations: we recall that the errors shown here are most probably underestimated (see Sect. 3.3) and that the tracks correspond to single stars, not to binary systems with strong reflection effects.

Fig. B.4 Histogram of the projected rotational velocities v sin(i) of the 18 nascent short-period binaries with extreme mass ratios (Table B.3) and of the ratios between these velocities to the synchronization velocities (last bin comprises all v sin(i)/vsync > 3.5 cases).

All Tables

All Figures

	Fig. 1 Line profiles of He I λ4388 Å in the most blue- and redshifted spectra (blue and red lines, respectively) for AN Dor, HD 309317, and V1061 Cen. The spectra were corrected for the measured RVs, and the black line provides their difference (multiplied by four and then shifted by +1 to facilitate the comparison with the lines).
In the text

	Fig. 2 Light curve observed with ASAS (four epochs; see text) and TESS for AN Dor (three sectors shown). The separation between the primary and secondary eclipses evolves.
In the text

Fig. 3 Best-fit adjustment to the differences between the observed times of minima of AN Dor and the linear ephemerides O − Clinear for i = 80°. The orange and pink symbols show the eclipses by the primary and secondary stars, respectively. Data prior to 2010 are taken from the ASAS survey, and data after 2018 are taken from TESS photometry. The long and short dashed lines yield our best-fit theoretical relations for the times of minima. In this plot, the linear ephemerides are T0 = 2 460 204.6031 and Psd = 2.0326782 d.

In the text

	Fig. 4 Comparison between the measured RVs of AN Dor (blue symbols with error bars) and the best-fit orbital solution explicitly accounting for apsidal motion.
In the text

Fig. 5 Left : differences between observed and predicted eclipse times based on the ephemerides of Table A.2. The data for the deeper eclipse are shown by blue dots, and the data for the shallower eclipse are shown as red circles. Close-ups are provided for TYC 1881-933-1. Right : Fourier spectrum of the O − C values of the times of primary minimum of TYC 1881-933-1 (top panel) and the associated spectral window (bottom panel).

In the text

	Fig. 6 Light-travel time effect adjustment of the O − C values of the times of primary minimum of the inner eclipsing binary of TYC 1881-933-1. Phase 0.0 corresponds to the periastron passage of the outer orbit.
In the text

	Fig. 7 Result of the disentangling for V1061 Cen. The main component is shown in blue, and the fainter component is shown in red (continuum shifted by +0.2 for clarity). No adjustment for the light ratio was made. The line depths are those observed in the input spectra, i.e., with respect to the combined light of all stars in the system.
In the text

Fig. 8 Comparison of mass and age for secondary stars of extreme mass ratio systems. Our targets are shown by red symbols, the MC systems of Moe & Di Stefano (2015) are shown with black circles, and the other Galactic cases are shown with orange triangles or stars (see text for details). The ZAMS track was defined from models by Tognelli et al. (2011), considering a central hydrogen abundance Xc equal to 99.5% of the initial value (which corresponds to the stabilization of the total nuclear-based luminosity in these models).

In the text

	Fig. B.1 Individual light curves of the targets - sample. The full set is provided as an additional pdf on Zenodo.
In the text

Fig. B.2 Binned and averaged light curves of all targets (bottom panels - the data used to produce these curves are identified by boldface in Table A.1; residuals are shown if a fitting was achieved), together with RV curves (top panels). The folding was done using the ephemerides of Table A.2. The best-fit solutions from Table A.3 are superimposed to the data using red lines. For LS I +61 275, the black RV points come from APOGEE while the red ones come from Huang & Gies (2006a), showing the discrepancy between them. For V1061 Cen, the black filled dots and blue circles correspond to the main and fainter spectral components, respectively (see Table B.2 and text for details).

In the text

Fig. B.3 Hertzsprung-Russell diagrams for primaries (top left, with the post-PMS evolutionary tracks of Brott et al. (2011) without initial rotation in blue and with the fastest initial rotation available in these models in violet lines) and secondaries (top right, with the PMS tracks of Tognelli et al. 2011 in green). Isochrones are shown in lighter color and with dashed lines. Symbols as in Fig. 8. For our targets, the secondary effective temperatures are those without reflection and the bolometric luminosities were calculated from these temperatures and the stellar radii (Table A.3) while the bolometric luminosities of the primaries are those from Table A.1. Bottom panels provide the mass-radius diagrams, using the same conventions. The comparison with theoretical evolutionary tracks is not without limitations: we recall that the errors shown here are most probably underestimated (see Sect. 3.3) and that the tracks correspond to single stars, not to binary systems with strong reflection effects.

In the text

	Fig. B.4 Histogram of the projected rotational velocities v sin(i) of the 18 nascent short-period binaries with extreme mass ratios (Table B.3) and of the ratios between these velocities to the synchronization velocities (last bin comprises all v sin(i)/vsync > 3.5 cases).
In the text