© The Authors 2026

1. Introduction

Most hot subdwarf stars (sdO/Bs) are evolved core-helium-burning objects that populate the extreme horizontal branch of the Hertzsprung–Russell diagram (Heber 2016). For sdO/Bs to be formed, their progenitors need to evolve up to the red giant phase and subsequently lose all or most of their outer layers on or close to the tip of the red giant branch, a process driven by binary interactions (Han et al. 2002). Consequently, about one third of these stars are found in short-period binary systems, with mostly brown dwarf, main-sequence red dwarf, or white dwarf (WD) partners (Schaffenroth et al. 2022).

The present-day masses of sdB stars retain information about these binary interactions. Most observed sdB stars have masses clustered around the canonical value for degenerate helium ignition, M ≃ 0.47 M_⊙, but both lower and higher masses are possible if the progenitor exceeded ∼2 M_⊙ and developed a nondegenerate core (Arancibia-Rojas et al. 2024).

Massive short-period sdO/B + WD binaries are particularly interesting because some of them are regarded as potential progenitor systems of type Ia supernovae (SNe Ia; Liu et al. 2023). These powerful transient events are triggered when runaway carbon burning ignites in a C/O white dwarf star, and they are very relevant for distance determinations at the largest scales and therefore cosmology. Hot subdwarf binaries qualify for the double-degenerate (DD) SN Ia channel if their total mass is close to or exceeds the Chandrasekhar limit (M_ch ≈ 1.4 M_⊙) and if they are close enough for gravitational radiation to merge them on the scale of a Hubble time. Because an sdO/B turns into a white dwarf on a timescale of only ∼100 Myr, such systems become WD + WD pairs well before coalescence, given that the post-common envelope period is sufficiently long. To date, three sdO/B + WD systems are known candidates for the super-Chandrasekhar DD merger scenario¹ (Geier et al. 2007; Pelisoli et al. 2021; Luo et al. 2025).

Binary-evolution models also predict sdO/Bs with neutron star (NS) or black hole (BH) companions (Wu et al. 2018), and sufficiently short period sdO/B + NS binaries are predicted to be progenitors of low-mass X-ray binaries (Tauris & van den Heuvel 2023). Despite theoretical predictions, no sdO/B in a confirmed short-period orbit with an NS or BH companion has been found to date, and existing candidates are either wide binaries or systems with poorly constrained parameters (Maitra et al. 2011; Geier et al. 2015, 2023).

In this paper we report the discovery of the massive short-period sdOB + WD binary LAMOST J065816.72+094343.1 (hereafter J0658). J0658 was first observed by the Large Sky Area Multi-Object Fiber Spectroscopic Telescope (LAMOST) and classified as a sdOB type hot subdwarf by Lei et al. (2018). In Section 2 we present our spectroscopic observations and data reduction. Atmospheric and orbital parameters as well as component masses are derived in Sect. 3. We discuss the possible evolutionary future of the system in Sect. 4, and we summarize our conclusions in Sect. 5.

2. Observations

Five archival LAMOST spectra were available for J0658, with a resolution of ∼3.05 Å and covering a wavelength range of 3700 − 9000 Å. Spectral follow-up observations for J0658 were then conducted using the Southern Astrophysical Research (SOAR) telescope’s Goodman High Throughput Spectrograph (Clemens et al. 2004). In total, 42 spectra with a resolution of ∼2.7 Å and covering a wavelength range of 3600 − 5300 Å as well as three spectra with a resolution of ∼0.58 Å and covering 3740 − 4350 Å were obtained across five nights between September 2023 and March 2025. All SOAR spectra were taken with an exposure time of 300 seconds and were reduced and barycentrically corrected using the Multi-Instrument Data Input Reducer (MIDIR²) pipeline. A further nine spectra were taken using the UV-Visual Echelle Spectrograph (UVES) at the ESO Very Large Telescope UT2 (VLT-UT2; Dekker et al. 2000) on the night of 1 March, 2025. These spectra were exposed for 300 seconds, yielding a wavelength coverage of 3300 − 4500 Å in the blue and 5700 − 6650 Å in the red arm at a resolution of R = 40 000. The UVES spectra were reduced using the EsoReflex pipeline (Freudling et al. 2013). Finally, seven additional spectra were taken with the Alhambra Faint Object Spectrograph and Camera (ALFOSC) at the Nordic Optical Telescope (NOT) on 21 April 2025 as well as on the following night. The ALFOSC spectra covered a wavelength range from 3400 − 5350 Å with a resolution of 2.2 Å and were also reduced using the MIDIR pipeline. In combination with the archival LAMOST data, these spectra cover a very long baseline of 4801 days (see Fig. 1).

Fig. 1. Phased RV curve and phased and binned TESS light curve for J0658. Both curves are phased to the same period, P, and zero point, t0 (see Table B.1). For the RV curve, the plots show individual RV measurements from the LAMOST survey (salmon), the SOAR Goodman spectrograph (green; Clemens et al. 2004), the ALFOSC spectrograph (light blue), and the UVES spectrograph (ocher). The light curve shows 118287 TESS data points from TESS sectors 33 and 87 that were phased and then binned into 175 bins. Both curves show their respective best-fitting models. For the RV curve, the fit was obtained via the MCMC method, and for the TESS light curve it was obtained using the LCURVE code.

3. Analysis

3.1. Atmospheric parameters and radial velocities

We derived the atmospheric parameters by fitting synthetic spectra to the UVES data via χ² minimization (Irrgang et al. 2014). Synthetic spectra were calculated using combined LTE and non-LTE model grids computed with the ATLAS/DETAIL/SURFACE (ADS) stellar atmosphere and spectral synthesis codes (Irrgang et al. 2018 and references therein; Appendix C.2). Atmospheric parameters from the model fits (Table B.1, Fig. D.1) identify the visible star as a helium-poor sdOB (T_eff = 35 800 ± 750 K, log g = 5.37 ± 0.07, log(n_He) = − 1.76 ± 0.05). The spectral fit also yielded a projected rotational velocity of v sin i = 37.0 ± 1.0 km s⁻¹, which is well constrained due to the metal lines visible in the high-resolution UVES spectra that are modeled in addition to the Balmer and helium lines. A systematic uncertainty of 1.0 km s⁻¹ was added in quadrature to the formal uncertainties from the fit.

Radial velocities (RVs) were obtained by fitting the synthetic model template to all the barycentrically corrected spectra. Typical formal fitting uncertainties are ≃0.8 km s⁻¹ (UVES) and ≃5 km s⁻¹ (SOAR, ALFOSC, LAMOST). To correct for typical systematic offsets we added 3 km s⁻¹ in quadrature to the uncertainties of the UVES RVs and 15 km s⁻¹ to the others.

3.2. Spectral energy distribution

Fitting the observed spectral energy distribution (SED) using atmospheric priors yields the angular diameter, which we combined with the parallax to derive the stellar radius and mass. The χ²-minimization method used to construct and fit the SED is described in detail in Heber et al. (2018). For J0658, we adopted the Gaia DR3 parallax and uncertainty, incorporating the zero-point correction of Lindegren et al. (2021) and the uncertainty adjustment of El-Badry et al. (2021), and we obtained a precise parallax for the star of ϖ = 0.79 ± 0.05 mas. Interstellar extinction was accounted for using the reddening law of Fitzpatrick et al. (2019) with E(44 − 55) as a free parameter and R(55) = 3.02 (a galactic average). The SED (Fig. D.2) does not show an infrared excess. When fitting the SED with the same models as described in Sect. 3.1, we found the sdOB to have a large radius, R_sd = 0.309 ± 0.02 R_⊙, and a high mass, M_sd = 0.82 ± 0.17 M_⊙, compared to other He-poor hot subdwarfs (Latour et al. 2026). The color excess measured in the fit (E(44 − 55) = 0.1233 ± 0.0025 mag) is comparable to the value given by current reddening maps at the position of J0658 (E(44 − 55) = 0.1042 mag; Zhang & Green 2025).

3.3. Combined spectroscopic and photometric orbital solution

J0658 was observed in sector 33 and sector 87 of the Transiting Exoplanet Survey Satellite (TESS). We used the PDCSAP flux, which accounts for some crowding present in the TESS aperture (CROWDSAP = 0.64). The resulting light curve with a baseline of 1433 days between the 27-day observations revealed ellipsoidal variations, implying that the hot subdwarf is deformed by the gravitational influence of a compact companion.

We determined the orbital period by combining information from the light curve and RV data (Table B.2). We calculated the Lomb-Scargle Periodogram (Zechmeister & Kürster 2009) for the TESS light curve and used it as a prior for the period when fitting the RV curve via the Markov chain Monte Carlo (MCMC) approach. For a short-period hot subdwarf binary, tidal forces are expected to efficiently circularize the orbit. We limited the eccentricity to e < 0.04 in the RV curve and thus adopted e = 0 for our model, fitting the RVs with v_rad(t) = Ksin(2π(t − t₀)/P) + γ, where v_rad is the RV, K is the RV half-amplitude of the sdOB component, γ is the systemic velocity, t₀ is the zero point of the ephemeris, and P is the orbital period.

This approach yielded a precise period measurement of the system at P = 0.31955193 ± 0.00000028 d (7.669246 ± 0.0000067 h) and a RV half-amplitude of 186.2 ± 3.3 km s⁻¹ among the other orbital parameters listed in Table B.1.

3.4. Companion mass

We determined the inclination angle and used the binary mass function to derive the mass of the unseen companion. The sdOB in the binary, at its high mass and with its short orbital period, is expected to be tidally synchronized by theory (Ma & Fuller 2024). We also detected clear light curve variations due to ellipsoidal deformation and Doppler beaming. When comparing the phasing of the light curve with that of the RV curve, we found no phase shift in the tidal bulges, which indicates no deviation from synchronous rotation. We therefore conclude that the rotation of the sdOB in J0658 is synchronized to its orbital motion.

We fit the phased TESS data with the LCURVE code (Copperwheat et al. 2010) using MCMC sampling. Limb and gravity darkening coefficients used in the fit were sourced from Claret et al. (2020). Free parameters were the mass ratio, q = M_comp/M_sd; inclination, i; velocity scale, v_s = (1 + 1/q)K_sd/sin i; and scaled sdB radius, r_sd = R_sd/a = 2πR_sd/Pv_s. Priors for all but i were obtained through Monte Carlo sampling of the spectroscopic solution and SED, i.e., independent of the light curve. The results of the MCMC fit to the light curve can be seen in Fig. D.3. Due to the very small amplitude of the ellipsoidal modulation and the subsequently low signal-to-noise level of the TESS light curve, the fit only constrains the inclination to the wide range of $i_{\mathrm{LC}}={50}_{-12}^{+18}\mathrm{deg}$.

To obtain a constraint on the inclination angle under the assumption of tidal synchronization, we recognized that the rotational velocity of the hot subdwarf in such a system is v_rot = 2πR_sd/P, where R_sd is the radius of the hot subdwarf and P is the orbital period. Since we had measured v_rot sin i, we could then obtain the inclination angle as i = arcsin(P v_rot sin i/(2πR_sd)). Using this expression with the highly accurate v_rot sin i obtained with UVES, we obtained a value of $i_{v_{\mathrm{rot}}}=49.6_{-4.1}^{+5.2}\mathrm{deg}$. This inclination angle is consistent with the one obtained from the light curve, but more precise, and we used it for further calculations. Solving the binary mass function numerically with all determined values, we calculated the companion mass to be $M_{\mathrm{comp}}=1.{30}_{-0.26}^{+0.31}{\mathrm{M}}_{\odot }$.

4. Evolutionary future

The unusually massive yet unseen companion of J0658 allows several possible evolutionary pathways, depending on its nature. Its mass is consistent with either an ultramassive CO/ONe WD (Camisassa 2025) or an NS companion. To explore these possibilities, we modeled the WD companion cases with Eggleton’s stellar evolution code (Eggleton 1971), adopting the same assumptions for metallicity, angular-momentum loss, and the Roche-lobe overflow (RLOF) mass transfer as in Liu et al. (2018a). We computed one reference track with parameters close to the adopted values of the observed system (Fig. D.4) (M_i, WD,  M_i, comp,  log P_i) = (1.23 M_⊙,  0.79 M_⊙,   − 0.4955 log(d)).

The calculation started with a CO/ONe WD + He–MS binary. After ≃3.3 × 10⁷ yr, the He star fills its Roche lobe on the subgiant branch, initiating RLOF at $\dot{M}\sim {10}^{-7}{M}_{\odot }{\mathrm{yr}}^{-1}$. The accreted material triggers weak He shell flashes, steadily increasing the WD mass. After 2.3 × 10⁵ yr, the transfer rate falls below 6 × 10⁻⁸ M_⊙ yr⁻¹, and the flashes become strong, stopping further mass accumulation. About 3.3 × 10⁴ yr later, the donor’s envelope is exhausted, and it turns into a second WD, leaving a detached double-WD system. At the end of the simulation (1.9 × 10⁵ yr after RLOF), we obtained M_f, WD = 1.2576 M_⊙, M_f, comp = 0.75 M_⊙, and log(P_f/day) = − 0.4594, values comparable to those of the recently discovered double-WD SN Ia progenitor WD J181058.67+311940.94 (Munday et al. 2025). If the accretor is a CO WD, a SN Ia could follow via the DD channel in ∼30 Gyr, which is much longer than the Hubble time.

Depending on the specific initial parameter set chosen from the possible parameter space as constrained by the observations, the outcome diverges from the one in this detailed scenario (see Fig. 2). For combined masses M_sys ≲ 2.15 M_⊙, the system does not merge within a Hubble time and forms a wide CO/ONe WD + CO WD binary instead. For CO WD companions at higher masses, the system produces a SN Ia via the DD merger channel at coalescence. For ONe WD companions, the WD undergoes accretion-induced collapse into an NS soon after RLOF, potentially evolving further into an intermediate mass binary pulsar (Liu et al. 2018a). If the companion is already an NS (Fig. D.5), stable RLOF during sdOB expansion leads to a transient phase as an intermediate-mass X-ray binary, after which the system becomes an intermediate-mass binary pulsar, as in the ONe case.

Fig. 2. Results of the stellar evolution and BPS simulations. Blue points mark initial configurations producing wide, detached double-WD binaries; red points reach Mch, triggering an SN Ia or accretion-induced collapse. The background shading shows BPS population densities by companion type, and the orange diamond marks the adopted masses for J0658.

To assess the likelihood of the possible companion types we performed Monte Carlo binary population synthesis (BPS) using the Hurley rapid binary evolution code (Hurley et al. 2002) with methods similar to those in Liu et al. (2018b). We selected J0658-like systems by requiring masses and periods near the observed values (0.72 < M_sd < 0.92 M_⊙, 1.08 < M_comp < 1.4 M_⊙, and −0.8 < log Pⁱ < −0.2 log(d)). Results of this simulation can be seen overlaid in Fig. 2. The BPS results overlap best with the ONe-WD and NS scenarios, and the model explicitly disfavors the SN Ia outcome. For an NS companion, we estimate the maximum current X-ray flux through wind accretion as F_X, max ∼ 10⁻¹³ erg cm⁻² s⁻¹ (Sect. C.1), one magnitude larger than the eROSITA limit of 3.2 ⋅ 10⁻¹⁴ erg cm⁻² s⁻¹ at this position (Merloni et al. 2024). The nondetection in eROSITA does not exclude the NS scenario, however, since the calculated maximum flux is idealized and does not factor in detrimental factors such as absorption losses. Overall, the ONe WD or NS companion scenarios can currently be seen as the most likely based on the BPS results.

5. Conclusions

We have identified and characterized J0658 as a new, highly massive short-period hot subdwarf binary. Atmospheric analysis of the UVES spectra showed the visible component to be a helium-poor sdOB star with T_eff = 35 800 ± 750 K, log g = 5.37 ± 0.07, and a projected rotational velocity of v sin i = 37.0 ± 1.0 km s⁻¹. The likely tidally locked nature of the system allowed for a precise measurement of the inclination angle and a subsequent determination of the companion mass. With a primary mass of M_sd = 0.82 ± 0.17 M_⊙ and a companion mass of $M_{\mathrm{comp}}=1.{30}_{-0.26}^{+0.31}{\mathrm{M}}_{\odot }$, it is the most massive of the now four known potential DD merger SNe Ia progenitor systems.

The companion’s mass range overlaps with those of the most massive CO and ONe white dwarfs and extends into the low-mass neutron-star regime. Notably, the current data do not allow for an unambiguous identification among these possibilities. If the companion is a CO WD, binary-evolution calculations predict stable He accretion and a double-WD merger within ∼30 Gyr, leading to a SN Ia or the formation of a detached double WD. An ONe WD companion would instead lead to accretion-induced collapse and the formation of an intermediate-mass binary pulsar, which is the same outcome as expected for an NS companion, though the system would first go through a transient phase as a intermediate-mass X-ray binary in this case.

References

Appendix A: Acknowledgements

Acnowledgements. F. M., H. D., and M. P. are supported by the Deutsche Forschungsgemeinschaft through grants GE2506/17-1, GE2506/9-2 and GE2506/18-1, respectively. B.B. and C.D. are supported by faculty research start-up funds kindly provided by the UNC Department of Physics & Astronomy and the UNC College of Arts & Sciences. This research was supported by Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy - EXC 2121 "Quantum Universe" – 390833306. Co-funded by the European Union (ERC, CompactBINARIES, 101078773). Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council. Neither the European Union nor the granting authority can be held responsible for them. This study is supported by the National Natural Science Foundation of China (Nos 12288102, 12225304, 12090040/12090043, 12273105), the National Key R&D Program of China (No. 2021YFA1600404), the Youth Innovation Promotion Association of the Chinese Academy of Sciences (No. 2021058), the Yunnan Revitalization Talent Support Program (Yunling Scholar Project and Young Talent Project), the Yunnan Fundamental Research Project (No 202201BC070003, 202401AV070006, and 202201AW070011), and the International Centre of Supernovae, Yunnan Key Laboratory (No. 202302AN360001). Based on observations obtained at the Southern Astrophysical Research (SOAR) telescope, which is a joint project of the Ministério da Ciência, Tecnologia e Inovações (MCTI/LNA) do Brasil, the US National Science Foundation’s NOIRLab, the University of North Carolina at Chapel Hill (UNC), and Michigan State University (MSU). The data presented here were obtained in part with ALFOSC, which is provided by the Instituto de Astrofisica de Andalucia (IAA) under a joint agreement with the University of Copenhagen and NOT. Based on observations collected at the European Southern Observatory under ESO programme 114.28KP.001. This paper includes data collected with the TESS mission, obtained from the MAST data archive at the Space Telescope Science Institute (STScI). Funding for the TESS mission is provided by the NASA Explorer Program. STScI is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5–26555. The Pan-STARRS1 Surveys (PS1) and the PS1 public science archive have been made possible through contributions by the Institute for Astronomy, the University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society and its participating institutes, the Max Planck Institute for Astronomy, Heidelberg and the Max Planck Institute for Extraterrestrial Physics, Garching, The Johns Hopkins University, Durham University, the University of Edinburgh, the Queen’s University Belfast, the Harvard-Smithsonian Center for Astrophysics, the Las Cumbres Observatory Global Telescope Network Incorporated, the National Central University of Taiwan, the Space Telescope Science Institute, the National Aeronautics and Space Administration under Grant No. NNX08AR22G issued through the Planetary Science Division of the NASA Science Mission Directorate, the National Science Foundation Grant No. AST–1238877, the University of Maryland, Eotvos Lorand University (ELTE), the Los Alamos National Laboratory, and the Gordon and Betty Moore Foundation. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation. The national facility capability for SkyMapper has been funded through ARC LIEF grant LE130100104 from the Australian Research Council, awarded to the University of Sydney, the Australian National University, Swinburne University of Technology, the University of Queensland, the University of Western Australia, the University of Melbourne, Curtin University of Technology, Monash University and the Australian Astronomical Observatory. SkyMapper is owned and operated by The Australian National University’s Research School of Astronomy and Astrophysics. The survey data were processed and provided by the SkyMapper Team at ANU. The SkyMapper node of the All-Sky Virtual Observatory (ASVO) is hosted at the National Computational Infrastructure (NCI). Development and support of the SkyMapper node of the ASVO has been funded in part by Astronomy Australia Limited (AAL) and the Australian Government through the Commonwealth’s Education Investment Fund (EIF) and National Collaborative Research Infrastructure Strategy (NCRIS), particularly the National eResearch Collaboration Tools and Resources (NeCTAR) and the Australian National Data Service Projects (ANDS). This work presents results from the European Space Agency (ESA) space mission Gaia. Gaia data are being processed by the Gaia Data Processing and Analysis Consortium (DPAC). Funding for the DPAC is provided by national institutions, in particular the institutions participating in the Gaia MultiLateral Agreement (MLA). The Gaia mission website is https://www.cosmos.esa.int/gaia. The Gaia archive website is https://archives.esac.esa.int/gaia.

Appendix B: Supplementary tables

Appendix C: Supplementary calculations

C.1. Estimated X-ray luminosity

In order to constrain the nature of the unseen companion, we estimate the current X-ray luminosity of the binary for the scenario of a WD and NS companion, and compare it with the point-source detection limit of eROSITA. We use the best-fit mass-loss rate function determined by Krtička et al. (2016) to estimate the wind mass-loss rate Ṁ_w of the hot subdwarf star as

$$\begin{array}{c}\hfill {\dot{M}}_w=(3.59\pm 1.72)·{10}^{-12}{\mathrm{M}}_{\odot }{\mathrm{yr}}^{-1}.\end{array}$$(C.1)

Using a terminal wind velocity of the star estimated from the model values in Krtička et al. (2016) as v_w,∞ ≈ 1420 km s⁻¹ together with the Bondi–Hoyle–Lyttleton accretion rate and the best-fit determined system parameters, we obtained

$$\begin{array}{cc}\hfill {\dot{M}}_{\mathrm{acc}}& =\alpha \frac{4\pi {(G{M}_{\mathrm{comp}})}^2\rho }{v_{w,\infty }^3}\approx \frac{G^2{M}_{\mathrm{comp}}^2{\dot{M}}_w}{a^2{v}_{w,\infty }^4}\hfill \end{array}$$(C.2)

$$\begin{array}{cc}\hfill & =(8.58\pm 6.16)·{10}^{-15}{\mathrm{M}}_{\odot }{\mathrm{yr}}^{-1}.\hfill \end{array}$$(C.3)

The corresponding X-ray luminosity can be approximated using the accretion rate-luminosity relationship (Frank et al. 2002; Eq. 1.5) and approximating R_NS ∼ 12 km and R_WD  ∼  3000 km based on neutron star equation of state constraints from Özel & Freire (2016) and the well-known empirical mass-radius relationship for white dwarfs by Nauenberg (1972):

$$\begin{array}{cc}& L_{\mathrm{X}},\mathrm{N}\mathrm{S}\approx G{M}_{\mathrm{comp}}{\dot{M}}_{\mathrm{acc}}/R_{\mathrm{comp}}=(7.7\pm 5.8)·{10}^{31}\mathrm{erg}/\mathrm{s}\hfill \end{array}$$(C.4)

$$\begin{array}{cc}& L_{\mathrm{X}},\mathrm{W}\mathrm{D}\approx (3.1\pm 2.3)·{10}^{29}\mathrm{erg}/\mathrm{s},\hfill \end{array}$$(C.5)

which, at a distance of 1266 ± 80 pc (as calculated from the adopted parallax) corresponds to an ideal flux of

$$\begin{array}{cc}& F_{\mathrm{X}},\mathrm{N}\mathrm{S}=\frac{L}{4\pi d^2}=(4.0\pm 3.1)\times {10}^{-13}\mathrm{erg}{\mathrm{cm}}^{-2}{\mathrm{s}}^{-1}\hfill \end{array}$$(C.6)

$$\begin{array}{cc}& F_{\mathrm{X}},\mathrm{W}\mathrm{D}=(1.60\pm 1.22)\times {10}^{-15}\mathrm{erg}{\mathrm{cm}}^{-2}{\mathrm{s}}^{-1}\hfill \end{array}$$(C.7)

for J0658. This calculation does not include potentially detrimental effects such as absorption through neutral hydrogen, which means that the true flux of the system is likely lower than this value.

C.2. Reliability of the spectroscopic solution

The reliability of the atmospheric parameters we determine through model fitting is central to our analysis of J0658. We thus give some additional details on our method here and compare our results to the literature.

The synthetic spectra were calculated using a state-of-the-art hybrid LTE and non-LTE approach using the ADS code and are described in detail by Möller (2021). In addition to hydrogen and helium lines, 24 additional elements are modeled in these synthetic spectra, though only N II, Si III and Si IV metal lines are actually present in the UVES spectra with equivalent widths of > 10 mÅ. The atmospheric parameters were determined via χ² minimization by fitting these synthetic spectra simultaneously to all individual UVES data, while keeping a common set of atmospheric parameters.

Luo et al. (2021) previously classified J0658 as an sdOB-type hot subdwarf, finding a divergent surface gravity (log g _Luo+ = 5.23). Their analysis is based on a single co-added low-resolution spectrum from LAMOST DR7 with a combined exposure time of 45 minutes. Since the exposure time point of this spectrum coincides with the RV curve’s turnaround point, it is affected by significant smearing due to an RV change of ∼36 kms⁻¹. When fitting the individual LAMOST exposures making up the co-added spectrum Luo et al. (2021) used with our models and method, we get atmospheric parameters that are closely aligned with our UVES-based solution (T_eff = 37000 ± 1000 K, log g = 5.32 ± 0.08, log(n_He) = − 1.76 ± 0.07). The discrepancy in derived parameters can therefore be attributed to the differing models, fitting methods and data quality between the two works.

Appendix D: Supplementary figures

Fig. D.1. Stacked and re-binned (50000 bins) UVES spectra for J0658. The black line represents the spectral data, the red line the best-fitting model spectrum. The gap between the end of the UVES blue arm and red arm (4500Å − 5700Å) and sections of the spectrum without any absorption lines are not shown.

Fig. D.2. Spectral energy distribution of J0658. Valid fluxes from all available photometric surveys are shown. The best-fitting synthetic spectral model is overlaid as a gray line. Horizontal error bars represent the filter width at tenth maximum. (Data publications: Pan-STARRS1 (yellow): Chambers et al. 2016; 2MASS (salmon): Cutri et al. 2003; SkyMapper (green): Onken et al. 2019; and Gaia G/BP/RP (blue) and Gaia XP (light blue): Gaia Collaboration 2022).

Fig. D.3. Posterior distributions of the free parameters in the light curve fit, obtained via MCMC sampling with 10 000 000 samples of the TESS light curve of J0658.

Fig. D.4. Binary evolution model for J0658 assuming a CO–white-dwarf companion. Left: Hertzsprung–Russell track of the helium star (solid line) and simultaneous change of the orbital period (dash-dotted line) as well as the current position of J0658 (orange diamond). Right: Mass-transfer history after the onset of Roche lobe overflow. Solid, dashed and dash-dotted curves denote the donor’s mass-loss rate, the WD accretion rate, and the growing WD mass, respectively. Insets list the adopted initial parameters and the final values reached at the end of the computation.

Fig. D.5. Binary evolution model for J0658 assuming an NS companion. Contents are analogous to those in Fig. D.4.

All Tables

All Figures

Fig. 1. Phased RV curve and phased and binned TESS light curve for J0658. Both curves are phased to the same period, P, and zero point, t0 (see Table B.1). For the RV curve, the plots show individual RV measurements from the LAMOST survey (salmon), the SOAR Goodman spectrograph (green; Clemens et al. 2004), the ALFOSC spectrograph (light blue), and the UVES spectrograph (ocher). The light curve shows 118287 TESS data points from TESS sectors 33 and 87 that were phased and then binned into 175 bins. Both curves show their respective best-fitting models. For the RV curve, the fit was obtained via the MCMC method, and for the TESS light curve it was obtained using the LCURVE code.

In the text

	Fig. 2. Results of the stellar evolution and BPS simulations. Blue points mark initial configurations producing wide, detached double-WD binaries; red points reach Mch, triggering an SN Ia or accretion-induced collapse. The background shading shows BPS population densities by companion type, and the orange diamond marks the adopted masses for J0658.
In the text

	Fig. D.1. Stacked and re-binned (50000 bins) UVES spectra for J0658. The black line represents the spectral data, the red line the best-fitting model spectrum. The gap between the end of the UVES blue arm and red arm (4500Å − 5700Å) and sections of the spectrum without any absorption lines are not shown.
In the text

Fig. D.2. Spectral energy distribution of J0658. Valid fluxes from all available photometric surveys are shown. The best-fitting synthetic spectral model is overlaid as a gray line. Horizontal error bars represent the filter width at tenth maximum. (Data publications: Pan-STARRS1 (yellow): Chambers et al. 2016; 2MASS (salmon): Cutri et al. 2003; SkyMapper (green): Onken et al. 2019; and Gaia G/BP/RP (blue) and Gaia XP (light blue): Gaia Collaboration 2022).

In the text

	Fig. D.3. Posterior distributions of the free parameters in the light curve fit, obtained via MCMC sampling with 10 000 000 samples of the TESS light curve of J0658.
In the text

Fig. D.4. Binary evolution model for J0658 assuming a CO–white-dwarf companion. Left: Hertzsprung–Russell track of the helium star (solid line) and simultaneous change of the orbital period (dash-dotted line) as well as the current position of J0658 (orange diamond). Right: Mass-transfer history after the onset of Roche lobe overflow. Solid, dashed and dash-dotted curves denote the donor’s mass-loss rate, the WD accretion rate, and the growing WD mass, respectively. Insets list the adopted initial parameters and the final values reached at the end of the computation.

In the text

	Fig. D.5. Binary evolution model for J0658 assuming an NS companion. Contents are analogous to those in Fig. D.4.
In the text