© The Authors 2025

1 Introduction

More than 97% of all main-sequence stars in our Galaxy (those with M ≤ 8–11 M_⊙) will become or have already become white dwarfs. With the cessation of nuclear reactions in their interiors, these stellar remnants slowly cool while releasing their thermal energy. This cooling process can extend over more than 10 billion years, which is why white dwarfs provide important information about the history of our Galaxy and its constituents (Althaus et al. 2010; García-Berro & Oswalt 2016; Isern et al. 2022).

White dwarfs are usually found isolated, but they can also be part of binary systems. In progenitor main-sequence binaries, each component evolves as a single star if the system is sufficiently separated (≥ 10 AU; Willems & Kolb 2004; Farihi et al. 2010). The more massive star leaves the main sequence first and eventually loses mass to produce a white dwarf, thus resulting in the expansion of the orbit and the formation of a widely separated white dwarf plus main-sequence star binary (WDMS). Wide WDMSs are excellent tools for tackling, for instance, the age–metallicity (Rebassa-Mansergas et al. 2021a) and age–velocity dispersion relations (Raddi et al. 2022) of the Galactic disk, as well as the age–activity–rotation relation of low-mass main-sequence stars (Rebassa-Mansergas et al. 2023; Chiti et al. 2024) since white dwarf ages can be accurately derived in these systems (Toloza et al. 2023).

If, conversely, the original main-sequence stars are close enough, mass transfer interactions are likely to ensue when the more massive star becomes a giant or an asymptotic giant star. This generally results in an unstable process of mass transfer, which leads to the formation of a common envelope surrounding both stars (Paczynski 1976; Webbink 2008). Due to friction involving both the giant's core (i.e., the future white dwarf) and the main-sequence companion with the material of the envelope, the orbital separation reduces dramatically and the orbital energy released is used to eject the envelope. As a consequence, the newly formed WDMS is in a much closer orbit with periods ranging from a few hours to several days (Nebot Gómez-Morán et al. 2011; Rebassa-Mansergas et al. 2012). These post-common-envelope binaries experience further angular momentum losses due to gravitational wave radiation and/or magnetic braking, which continues to shorten their orbital periods until a second phase of (stable) mass transfer begins producing exotic objects such as cataclysmic variables or super-soft X-ray sources. A second common envelope can begin to form if the secondary star has time to leave the main sequence, which leads to the formation of a double white dwarf binary (DWD) in a tight orbit.

Close WDMSs are excellent tools for constraining a wide range of open problems in astronomy. To begin with, the reconstruction of their past evolution places direct constraints on how efficiently the release of orbital energy is used to expel the common envelope. Thus, since the pioneering work of Zorotovic et al. (2010), this process has consistently been found to be inefficient (Camacho et al. 2014; Cojocaru et al. 2017; Toonen & Nelemans 2013; Iaconi & De Marco 2019; Zorotovic & Schreiber 2022). By predicting the future evolution of close WDMSs, it is also possible to analyze their different evolutionary paths, which, ultimately, can be used to place constraints on the different channels toward type Ia supernovae (Parsons et al. 2016, 2023; Hernandez et al. 2021, 2022). Moreover, eclipsing WDMSs enable us to derive independent masses and radii for both components, which in turn enables us to test the theoretical mass-radius relation of white dwarfs (Parsons et al. 2017; Brown et al. 2023), main-sequence stars (Parsons et al. 2018), and even subdwarf stars (Rebassa-Mansergas et al. 2019) and brown dwarfs (French et al. 2024).

Despite the significance of WDMSs, it is important to bear in mind that the currently known population is severely affected by selection effects. This is particularly the case for unresolved systems since the main-sequence components generally outshine the white dwarfs in the optical (Rebassa-Mansergas et al. 2010). One way to overcome this issue is by making use of ultraviolet spectroscopy and/or photometry, which minimizes the contribution of the companions at those wavelengths (Parsons et al. 2016; Rebassa-Mansergas et al. 2017; Ren et al. 2020; Anguiano et al. 2022; Nayak et al. 2024; Sidharth et al. 2024). An additional problem is that most WDMS catalogs are the result of mining data from magnitude-limited samples, and, to make things worse, in most cases the WDMSs were observed following selection criteria designed for targeting different objects. Volume-limited samples have been possible to achieve thanks to the excellent photometry and astrometry provided by Gaia (Inight et al. 2021). Using Early Data Release 3, RebassaMansergas et al. (2021b) acquired and characterized the most homogeneous volume-limited sample of unresolved WDMSs to date within 100 pc. Although this sample is still far from complete – since it only includes the bridge region in the Gaia color-magnitude diagram between single main-sequence stars and single white dwarfs – it is a good benchmark for testing theoretical models and assumptions because it can be used to generate synthetic populations for comparison.

In this sense, binary population synthesis studies have been widely used, particularly in the study of DWD and WDMS systems (e.g., Nelemans et al. 2001b,a; Toonen et al. 2014, 2017; Camacho et al. 2014; Cojocaru et al. 2017; Canals et al. 2018). More recently, Torres et al. (2022), taking advantage of the accurate astrometric and photometric data provided by the Gaia mission, succeeded in tuning the model parameters so that the synthetic population reproduces the observed percentages of resolved double-degenerate and WDMS systems. Among these parameters, they identified accurate prescriptions for the binary fraction, initial mass ratio distribution, and initial separation distribution, among others.

The population synthesis code used in Torres et al. (2022), largely developed by our group, incorporates a comprehensive description of the most relevant physical processes in stellar binary evolution (mass transfer, tidal interactions, stellar wind losses, magnetic braking, gravitational radiation, etc.), as well as a realistic characterization of the different Galactic components – thin and thick disks and the halo – for a wide range of masses and metallicities. In addition, a detailed set of selection criteria has been implemented to achieve significant statistical comparisons between the observed and synthetic samples.

We used this consolidated binary population synthesis code (Torres et al. 2022) to study the 100 pc unresolved WDMS sample derived from Gaia data (Rebassa-Mansergas et al. 2021a). In particular, we evaluated the degree of completeness of the observed sample and provide insights into the evolution of close binary systems, especially parameters related to the commonenvelope phase.

In Sect. 2 we describe the Gaia WDMS sample. In Sect. 3 we provide details of our Monte Carlo simulator, which produces synthetic samples of WDMSs to be directly compared to the observed ones. In Sect. 4 we present our results, and we conclude in Sect. 5.

2 The Gaia 100 pc WDMS population

To test the models used and the assumptions made in our simulations (details are provided in Sect. 3), we compared the synthetic sample of unresolved WDMSs with the observed population obtained from Gaia within 100 pc. This distance limit was imposed with the aim of dealing with the largest nearly complete observed volume-limited sample (Jiménez-Esteban et al. 2018). To that end, we adopted the unresolved WDMS catalog of Rebassa-Mansergas et al. (2021b), who searched for candidates within the Gaia color-magnitude region between the single white dwarf and the single main-sequence star loci (this is represented by solid blue lines in Fig. 1). The sample initially included all Gaia sources with parallax and the three band relative error fluxes under 10%. Moreover, objects with a large excess factor were excluded from the list following Riello et al. (2021), which resulted in 2001 WDMS candidates within the considered region. To further exclude contaminants, RebassaMansergas et al. (2021b) used the Virtual Observatory SED Analyser (VOSA; Bayo et al. 2008) to fit the spectral energy distributions (SEDs) of the 2001 candidates with single white dwarf models, single main-sequence star models and a combination of the two. Hence, only those with two-body fits and with VOSA parameter Vgf_b¹ < 15 were kept. During this exercise, Rebassa-Mansergas et al. (2021b) considered only objects with enough measurements to carry out the SED fitting and without wrong photometric values. This reduced the sample to 280 targets, of which 129 were further excluded due to bad fitting arising as a consequence of ultraviolet excess likely arising from magnetic activity of single main-sequence stars (3), bad SED fitting (90) or contamination from nearby sources (36). Closer inspection of the remaining 151 candidates revealed that 34 had unreliable Gaia astrometric solutions (RUWE > 2, or astrometric_excess_noise > 2 and astrometric_excess_noise_sig > 2) and were also excluded. Finally, 5 objects turned out to be cataclysmic variables thus leaving 112 good WDMS candidates.

The VOSA two-body fits yielded the following stellar parameters for both components of these 112 objects: effective temperature, radius and bolometric luminosity. White dwarf masses (hence surface gravities) were obtained interpolating the effective temperatures and radii in the La Plata cooling sequences (e.g., Camisassa et al. 2016; Camisassa et al. 2019). For 98 of the 112 objects, the measured parameters were considered as reliable by Rebassa-Mansergas et al. (2021b).

Fig. 1 Color-magnitude diagrams for the Gaia 100 pc sample. The sample construction process is exemplified with color-magnitude diagrams for the observational sample for: the different selection phases (left panels; see Table 1); the no-cuts samples (upper panels); the flux and parallax selection (upper-middle panels); the excess factor removal (lower-middle panels); and the final WDMS sample (bottom panels). The corresponding diagrams of the simulated sample are also shown (right panels). The WDMS region (Rebassa-Mansergas et al. 2021b) is marked with a solid blue line, and the dashed blue line depicts the magnitude limit of our simulations. Finally, the selection function (green line contours corresponding to 15%, 36%, and 49%) is derived from the observed sample from Rebassa-Mansergas et al. (2021 b, bottom-left panel, magenta dots; RM+21) and applied to a typical synthetic realization (bottom-right panel, magenta dots).

3 Population synthesis modeling

3.1 The Monte Carlo simulator

To generate a synthetic population of both single and binary stars in the Galaxy, we used the Monte Carlo simulator developed by our group (Torres et al. 1998; García-Berro et al. 1999). The code, named MRBIN, not only is capable of reproducing realistic single main-sequence and single white dwarf Galactic populations (e.g., Torres et al. 2001, 2019), but also allows us to reproduce main-sequence binary stars and their subsequent evolution (e.g., Torres et al. 2022). The binary stellar evolution is modeled following the binary stellar evolution (BSE) code developed by Hurley et al. (2002), with particular updates on white dwarf binaries (e.g., Camacho et al. 2014; Cojocaru et al. 2017; Canals et al. 2018).

The modeling of the entire single and binary population requires a series of input parameters that we adopted from Torres et al. (2022). In that work, we used the same code, MRBIN, to perform a population synthesis fit of the current populations of resolved WDMS and DWD systems within 100 pc that satisfactory reproduced the observed rates by Gaia, namely 6.31% for WDMS systems, and 1.18% for DWDs. This was achieved by considering a main-sequence binary fraction of 32%; an initial mass function from Kroupa (2001) for the range 0.08–50 M_⊙, and extending that distribution with a positive slope according to Sollima (2019), to reproduce the primary mass distribution the more massive star in the main-sequence binary M₁–; an initial thermal eccentricity distribution following Heggie 1975; an orbital separation of the form f(a) = a⁻¹, and a mass ratio distribution following n(q) ∼ q^−1.13 (with q = M₁/M₂, where M₂ is the secondary, less massive, star mass). The rest of parameters can be found in Torres et al. 2022.

Moreover, the simulator randomly assigns each single and binary star to a Galactic component according to the ratio 85:10:5 for the thin disk, thick disk and halo, respectively, following Torres et al. (2019). Objects from the different Galactic components are modeled following different criteria regarding their ages, metallicities, spatial distributions and kinematics. In what follows, we summarize the rest of Galactic inputs, while a detailed description and their corresponding references can be found in the aforementioned paper. For the thin disk stars, we adopted an exponential spatial distribution with a scale height of 250 pc and a scale length of 2.6 kpc, and a constant star formation rate with a disk age of 9.2 Gyr. For thick disk stars we considered a burst of star formation 10 Gyr ago extended up to 12 Gyr with a scale height and scale length of 1.5 kpc and 3.5 kpc, respectively. For halo stars we assumed an isothermal spatial distribution with a constant burst of 1 Gyr that took place 12.5 Gyr in the past. In the three cases the age-metallicity relation follows a constant model (Rebassa-Mansergas et al. 2021a). That is, we randomly generated metallicity values following a Gaussian distribution with a mean and a dispersion of Z = 0.02 and 0.005 for the thin disk, Z = 0.01 and 0.003 for the thick disk, and Z= 0.0015 and 0.001 for the halo, respectively.

Wide binary systems evolve in the same way as single stars, avoiding mass transfer episodes. However, for the evolution of binary systems with shorter separations, our code, based on the BSE code, uses a wide set of parameters to model the mass transfer and angular momentum losses. In particular, we adopted a common-envelope treatment characterized by the α_CE parameter (Tout et al. 1997, Webbink 2008) during a common-envelope phase, which quantifies the efficiency with which orbital energy is transferred from the binary system to the envelope of the donor star to eventually eject it. In our reference model, we used a value of α_CE = 0.3 (Camacho et al. 2014).

Once the systems (both single and binary stars) are evolved to present time, the MRBIN code allows the derivation of magnitudes in a given photometric system. For the purpose of this paper (which requires a direct comparison to the observed Gaia data; see Sect. 2), we derived the Gaia magnitudes following a three spline interpolation of luminosity, effective temperature, and metallicity in the PARSEC (Bressan et al. 2012) evolutionary tracks for the main-sequence stars. In particular, for luminosities higher than 0.1 L_⊙, we used the effective temperature, T_eff, directly provided by the BSE code. However, for luminosities below that value, we applied an improved correction to T_eff by performing a multi-variable interpolation using luminosity and metallicity values to obtain the corresponding T_eff according to the PARSEC evolutionary models. In any case, to avoid issues with photometric uncertainties in the models, our final simulated sample is restricted to absolute magnitudes M_G < 14.2 mag. For those systems that become white dwarfs, their cooling time is accurately determined using the most up-to-date evolutionary cooling sequences provided by the La Plata group (e.g., Camisassa et al. 2016; Camisassa et al. 2019), which take the different core and chemical compositions for the whole range of possible masses into account. The Gaia magnitudes are interpolated within the updated grid of white dwarf atmospheric models by Koester (2010), assuming all white dwarfs to be hydrogen-rich.

Based on the Gaia performance, we distinguished between resolved and unresolved synthetic binaries: we considered a binary to be unresolved if its angular separation was smaller than 2″ (Torres et al. 2022) and was thus observed by Gaia as just one stellar object instead of two separate well-characterized stars. It is worth noting that Gaia's resolution can nominally be even smaller (∼ 0.5″ Gaia Collaboration 2016), especially for objects of equal brightness, which is not the general case for WDMS systems. Therefore, we adopted a more conservative value to ensure completeness, as discussed in Torres et al. (2022). For the purpose of studying these unresolved systems, the MRBIN code generates the angular separation for each system from the corresponding parallax and orbital separation, taking into account the eccentricity, inclination, and phase of the orbit. As a result, some systems may be classified as unresolved even though they have large periods (see Sect. 4.2). For these synthetic unresolved systems it is required to consider new magnitude values that are calculated as an addition of the fluxes of both stars in the binary system. To do so, the absolute magnitudes previously obtained for each component are transformed into apparent magnitudes and subsequently converted into fluxes, which can then be added together and turned back into magnitudes.

In a final step, with the aim of closely mimicking the Gaia observational data, we introduced both photometric and astrometric errors following Riello et al. (2021)² to the synthetic data previously generated. We also introduced a synthetic excess factor (Riello et al. 2021) to each target to quantify the expected flux excess arising from contamination from other sources. This excess is quantified from the observational data as C^* = C– f(G_BP, G_RP), where C=(I_BP + I_RP)/I_G (Evans et al. 2018) and f(G_BP, G_RP) is a function that gives the expected excess at a given color. Objects affected by no flux excess have a C^* value close to zero. In our simulations, we used the previously obtained fluxes to determine the value of C and f(G_BP, G_RP) for each object, and if C^* < 1 were considered unaffected by excess.

3.2 Implementing the selection criteria and the selection function

To ensure that our simulations apply the same selection criteria as those used for the observed Gaia sample and to understand how these criteria impact our synthetic samples, we illustrate in Fig. 1 how the most important cuts affect the observed and synthetic populations (left and right panels, respectively). These cuts include: the initial selection based on fluxes and parallaxes with relative errors under 10% (top panels), the resulting sample after excluding objects with a high excess factor parameter (middle panels) and the final sample after excluding single stars (bottom panels). In these plots, it is easy to distinguish between the population of single white dwarfs and main-sequence stars. Between these two, there is a distinct region (following the selection criteria of Rebassa-Mansergas et al. 2021b), where unresolved WDMSs are the predominant objects. We recall that our synthetic models are restricted to main-sequence stars with M_G absolute magnitudes above 14.2 mag. Hence, for a proper comparison, only objects with M_G absolute magnitudes below this value were considered. This cut is illustrated in Fig. 1 as a horizontal blue line.

In Table 1, we provide the number of objects within the WDMS region that survive each selection cut for both the observed and simulated samples. When comparing these values some differences arise. To begin with, the number of initially selected objects is considerably larger in the observed sample. This can be easily understood since observed data are associated with uncertainties that are challenging to simulate, such as wrong or spurious astrometric solutions, large values of the excess factor parameter due to chance alignments, etc. Additionally, Gaia observations reveal a noticeable cluster of objects in the lower intermediate region between white dwarfs and mainsequence stars (see the top left panel of Fig. 1). These objects are primarily located in the Galactic plane, and therefore affected by crowding, and our current simulation code is unable to replicate them accurately. All these issues result in the inclusion of non WDMSs in our studied region, which are not represented in our simulations. Indeed, when the excess factor cut is applied, the number of objects that are excluded is significantly larger for the observed sample (∼ 83%) than for the simulated sample (∼ 16%). It should be noted that the Gaia observed sample contains many contaminated sources in the intermediate region, resulting in high color excesses or large astrometric errors. Our simulations, limited to Gaussian-distributed magnitude uncertainties, already remove most systems with large errors via the flux cut. Thus, applying the excess factor cut has a minimal effect, as few systems remain to be removed. However, after excluding single star contamination, the WDMS selection is very similar for the two samples (on the order of 240 objects). These samples are further refined by taking into account the following aspects that affect the observed sample: lack of sufficient photometry for performing the VOSA fits, detected UV excess due to magnetic activity, contamination from nearby sources and bad SED fitting (Sect. 2). The final observed sample is, thus, reduced to 112 candidates ( 100 if we consider only those with M_G < 14.2 mag). Implementing all these effects, which vary in nature and specific cases, in detail in our simulator would be a complex and overwhelming task. As a practical approach, to take them into account, we introduced a selection function (Rix et al. 2021), S_C, which assigns a probability to each member of our synthetic WDMS sample C of belonging to the final sample. The selection function is built as a function of the luminosity and color, that is, S_C = S_C(M_G, G_BP–G_RP), proportionally to the density of observed WDMS candidates (magenta dots in bottom-left panel of Fig. 1) in the Gaia color-magnitude diagram. In the bottomright panel of Fig. 1 we depict the contour-levels (green lines) of our selection function S_C. We observe a higher concentration of objects toward the lower edge of the main-sequence region, as theoretically expected (see, for instance, Fig. 6 of Rebassa-Mansergas et al. 2021b), compared to the corresponding dots in the observed sample. This indicates that the selection bias is even stronger than implemented. However, the significance of this observation is minor, as our analysis focuses on the total number of objects. That is, from the 240 simulated WDMSs, the selection function selects 113 objects, which is in very good agreement with the observed sample of 100 objects.

4 Analysis of the unresolved WDMS population

The synthetic single and binary populations, for both mainsequence and white dwarf stars, obtained from the model described in Sect. 3.1 show reasonably good agreement with the 100 pc observed sample from Gaia (bottom-right and left panels, respectively, of Fig. 1; see also Fig. B.1 for a detailed comparison of the cumulative distributions for the magnitude M_G and the color G_BP–G_RP under different selection criteria). Moreover, according to Rebassa-Mansergas et al. (2021b), 100 Gaia objects with M_G < 14.2 mag can be classified as unresolved WDMSs in the intermediate region between white dwarfs and main-sequence stars, in good agreement with our final synthetic sample of 113 objects. These facts support the validity of the assumptions implemented in our model. Consequently, we took advantage of our MRBIN code to analyze the general content of the white dwarf population, with a particular focus on the unresolved white dwarf binary population. An analysis of the different aspects is presented as follows.

Fig. 2 Percentages of white dwarfs and main-sequence stars (top chart) and those of the white dwarf subpopulations (bottom chart) within 100 pc. The complete dataset can be found in Table A.1.

4.1 General content of the white dwarf population

The total number of stars observed at the current time in the Gaia color-magnitude diagram within 100 pc in the white dwarf region is 12718 (Jiménez-Esteban et al. 2023). This value is used as a constraint for the normalization of our simulations. That is, our MRBIN code generates stars until the synthetic number of objects within the white dwarf region reaches this limit.

As previously stated, the input parameters used in our model reproduce the main features of the Gaia color-magnitude diagram for the entire population of stars reasonably well. Moreover, the simulations are able to reproduce the constraints on the percentages of resolved WDMSs and DWDs in the white dwarf region (Torres et al. 2022), as well as the number of unresolved WDMS candidates in the intermediate region (Sect. 3.2). Hence, at this point we were able to analyze the contribution of the different subpopulations. We emphasize that, in this work we considered only our reference model (that is, restricted to a specific set of parameters), and we leave the analysis of possible dependences with the model's input parameters for a future study.

In Fig. 2 we display a pie chart indicating the percentage of each individual population within the entire 100 pc Gaia color-magnitude diagram. In the top chart we show the contribution of the main populations: 67% are single stars ( 61% main-sequence and 6% single white dwarfs), while the remaining 33% are binary systems (30.5% formed by two main-sequence stars, and 2.5% containing a white dwarf). It is worth noting that, within the white dwarf population, 15.8% of white dwarfs are in binary systems. Recalling that our normalization criteria is fixed to the number of total simulated objects in the white dwarf region, the previous percentages imply a total of ∼ 100 000 main-sequence stars generated within 100 pc. This value is in agreement with Gaia Collaboration (2021).

In the bottom chart of Fig. 2 we analyze the resulting content of the different white dwarf sub-populations. As previously stated, most of them (84.2%) are isolated white dwarfs, with 70.2% formed as single white dwarfs and 14% resulting from binary systems that have merged. This means that, of the total sample of isolated white dwarfs, around 17% of them are products of merger episodes, in agreement with previous estimates (e.g., Temmink et al. 2020).

Regarding white dwarfs in unresolved binary systems, most are, as expected, WDMS systems, comprising 7.66% of the total white dwarf population, while only a residual 0.82% are DWD systems. This percentage of WDMSs implies that nearly 1000 unresolved systems are WDMSs. However, in the WDMS region, only 100 objects are ultimately identified. Thus, only a 10% of this population is accessible by Gaia under the conditions imposed by Rebassa-Mansergas et al. (2021b, see our Sect. 2). The remaining ∼ 90% are mainly located in the main-sequence region of the Gaia color-magnitude diagram. It is challenging to identify these objects as WDMSs since the main-sequence components outshine the white dwarfs in the optical. To identify such systems, it is common to use alternative methods involving ultraviolet photometry and/or spectroscopy (Parsons et al. 2016; Rebassa-Mansergas et al. 2017; Ren et al. 2020; Anguiano et al. 2022; Nayak et al. 2024; Sidharth et al. 2024).

The resolved WDMS and DWD populations will not be analyzed here since they are virtually identical to those of Torres et al. (2022) and were thus the dedicated subject of study in that work. Worth noting is that their percentages ( 6.06 for WDMSs and 0.93 for DWDs), as indicated in Fig. 2, are in very good agreement with the observed values in the white dwarf region ( 6.31 for WDMSs and 1.18 for DWDs; Torres et al. 2022). Finally, the percentage of systems formed by a white dwarf and a black hole or a neutron star is less than 0.3%, as expected for a 100 pc volume-limited sample.

4.1.1 The common-envelope treatment

According to Rebassa-Mansergas et al. (2021b), 100 Gaia objects with M_G < 14.2 mag can be classified as unresolved WDMSs in the intermediate region between white dwarfs and main-sequence stars. On the other hand, as we have shown (see Sect. 3.2), our reference model predicts a similar value of 113 objects. This value is actually the mean obtained from around 100 different Monte Carlo realizations, resulting a dispersion of ± 8. From the dispersion derived from our simulation and the Poissonian error associated with the observed value (i.e., 100 ± 10 observed objects), we can conclude that our reference simulation is in agreement with the observed sample.

We recall that for our reference model we used a value of the common-envelope efficiency parameter of α_CE = 0.3 (Torres et al. 2022; Camacho et al. 2014). However, depending on the adopted value of α_CE, we can derive a different number of synthetic unresolved WDMSs.

To constrain the possible α_CE values that yield synthetic unresolved WDMS samples in agreement with the observed value (taken as 100 ± 10), we compared the total number of synthetic WDMSs (after applying the selection function described in Sect. 3.2) to the observed one. In Fig. 3, we depict the total number of synthetic WDMSs across the full range of α_CE, with the average shown as a red line and the dispersion indicated by a red shaded area. Similarly, the observed value is shown in gray. As it can be seen, α_CE values larger than 0.37 and lower than 0.02 yield synthetic samples that start to deviate from the observed number and, in principle, can be ruled out. This confirms the preference of low values for the common-envelope efficiency during a common-envelope episode, a result that has been found in the past by several authors (Zorotovic et al. 2010; Toonen & Nelemans 2013; Camacho et al. 2014; Cojocaru et al. 2017; Iaconi & De Marco 2019; Zorotovic & Schreiber 2022).

Furthermore, our simulation also provides relevant information about how many unresolved WDMS systems that can be observed by Gaia within the WDMS region have gone through a common-envelope phase. In Fig. 3 we show synthetic post common-envelope systems in blue and synthetic systems that evolved avoiding mass transfer episodes in green. As expected, the number of systems that do not evolve through a common envelope is practically the same, around 80 (initially 160 prior to the application of the selection function). However, the number of WDMS systems that passed through a common-envelope phase increases with the adopted value of α_CE. This is a simple consequence of increasing the number of mergers during a common-envelope event for lower α_CE values, which imply a farther reduction of the orbit to eject the envelope.

We attempted to further constrain α_CE by considering how many of the synthetic WDMSs are eclipsing and comparing these numbers to the observed value. Of the 100 observed systems, Zwicky Transient Facility (ZTF; Dekany et al. 2020) observations (constrained to a Dec limit of −30°) indicate that 6 are totally eclipsing (Brown et al. 2023). By applying the same declination limit to our synthetic samples, we can determine how many of them are eclipsing by considering the minimum angle for total eclipse, i_t, as $$i_{\mathrm{t}}={90}^{\circ }-\arcsin \left(\frac{R_{\mathrm{M}\mathrm{S}}-R_{\mathrm{W}\mathrm{D}}}{a}\right),$$(1)

where R_MS and R_WD are the radius of the main-sequence and the white dwarf, respectively, and a the semimajor axis of the binary system. If the inclination of a given binary as obtained from our code is larger than i_t, then the simulated WDMS can be considered an eclipsing system. We also note that this geometrical criterion promotes systems with small orbital separations (hence short periods, typically a few to tens of hours), thereby increasing the likelihood of eclipsing, as observed in such systems. There is a seventh system that is partially eclipsing according to ZTF data (Van Roestel, in preparation), which has not been considered in this work. It is also important to mention that the selection effects for identifying these systems are not yet well understood, and there no published estimates of their completeness. Consequently, we used the six total eclipse systems as a lower limit for the 100 pc population.

The number of synthetic unresolved eclipsing WDMSs (blue line) as a function of α_CE is illustrated in Fig. 4. As in Fig. 3, the shaded regions indicate the corresponding dispersions and the gray area corresponds to the observed value within the Poissonian error. Given the low number of eclipsing systems, the only constraint we can derive from Fig. 4 is that α_CE should be higher than 0.08. From the analysis described in this section, we constrained the α_CE value to be in the range 0.1−0.4.

Fig. 3 Number of unresolved WDMS systems as a function of the αCE parameter for binary systems in the WDMS region (red), showing also how many have undergone a common-envelope episode (blue) and how many have not (green). The vertical dashed orange lines show the lower and upper limits of αCE compatible with the observed value (gray).

Fig. 4 Number of eclipsing binaries as a function of the αCE parameter for binary systems in the delimited WDMS region. As in Fig. 3, the vertical dashed orange line represents the lower limit of αCE compatible with the observed value (gray).

4.1.2 Stellar parameter analysis

In this section we compare the stellar parameters, namely masses and effective temperatures, of the synthetic WDMSs with those derived from the observed sample of Rebassa-Mansergas et al. (2021b). Based on the results from the previous section, only synthetic models in which the α_CE is approximately 0.1–0.4 should be used, and we fixed the value to α_CE = 0.3, which corresponds to our reference model (note that no major differences were found when using other possible values). This approach allows us to study the behavior of other parameters of our simulated sample.

Figure 5 represents the mass-effective temperature plane of the simulated (blue dots) and observed (red dots) samples for the main sequence (top figure) and the white dwarf (bottom figure) components. For each figure, we plot the corresponding distribution of mass (right panels) and effective temperature (top panels) for both the observed and synthetic samples. We note that the simulated sample combines results from a large number (100) of realizations, providing a better characterization of the statistical properties of the synthetic sample.

Inspection of the top diagram of Fig. 5 reveals a good agreement between the simulated and observed main-sequence parameters, with the exception of a couple of observed objects (bottom-right region) that are clearly outliers. These are likely WDMSs with unreliable two-body fits from VOSA. It is also worth noting that ∼ 10% of the observed companions have mass estimates below 0.1 M_⊙ and are probably brown dwarf candidates, objects that are not simulated by our code. The lack of such objects is also revealed when comparing the normalized distributions of mass and effective temperature (right and top panels). Additionally, we applied a Kolmogorov-Smirnov (KS) tests to the cumulative mass and effective temperature distributions (see the top and bottom panels, respectively, of Fig. B.2) yield KS probabilities of 0.04 and 0.03, respectively. When comparing these values with a typical statistical significance level of α = 0.05, we can argue that the null hypothesis can be rejected, that is, that the distributions are not likely the same. However, given that the obtained p-values are close to the threshold, this fact suggests that there is evidence of a difference (as previously stated), but it is only moderate. We excluded both observed and synthetic objects below 0.1 M_⊙, and the resulting KS probabilities are 0.40 for mass and 0.05 for effective temperature, indicating a strong agreement between the observed and simulated distributions.

Analogously, the observed and simulated white dwarf populations occupy similar regions in the mass-effective temperature plane (see the bottom diagram of Fig. 5). Indeed, the normalized effective temperature distributions are rather similar, with a tendency of increasing objects for decreasing temperatures. The KS probability obtained in this case is 0.17 (see the cumulative distributions in the bottom panel of Fig. B.3); clearly indicating compatibility between the two distributions. However, the observed sample displays a clear overabundance of low-mass white dwarfs with respect to the synthetic sample. Indeed, the KS test applied to the cumulative mass distributions (top panel of Fig. B.3) gives a null probability that the two samples are drawn from the same parent population. Even if we restrict our analysis to objects with masses above 0.35 M_⊙, the probability remains null. Along with the possibility that other nonstandard channels may exist for the formation of these objects, we speculate this is likely due to issues with the VOSA fits (Brown et al., in preparation, will provide a more detailed analysis), where discrepancies in white dwarf masses between different models could be as large as 20–30%.

Fig. 5 Distribution of mass and Teff for the main-sequence stars (top panel) and white dwarf component (bottom panel) of synthetic (blue) and Gaia-observed (red) WDMS systems. For a better statistical comparison, 100 different synthetic realizations are plotted.

Fig. 6 Initial (red) and final (blue) period distribution of simulated WDMS systems within the observed region.

4.2 Orbital period distribution and percentage of unresolved WDMSs evolving through the common envelope

Although some differences arise between the simulated unresolved WDMS sample and the one observed by Gaia and studied by Rebassa-Mansergas et al. (2021b, see our Sect. 4.1.2), we considered the synthetic sample as representative of the observed population. This allowed us to estimate how many of the observed unresolved systems have evolved through a common-envelope phase.

A clear way of illustrating this is by representing the orbital period distribution, which is shown in Fig. 6. The initial mainsequence binary orbital period distribution is illustrated in red and the resulting WDMS orbital period distribution in blue. It becomes clear that the final period distribution of WDMSs is bimodal, with the systems that evolved through the common envelope having periods of ≲ 100 days and the non-interacting binaries having much longer periods. We note that in the latter case the orbital periods increase due to the mass loss of the white dwarf progenitors. Furthermore, if we consider only those systems that have evolved through a common-envelope phase, we also observe a certain tendency to bimodality in their final distribution, consistent with recent analyses (Ashley et al. 2019; Lagos et al. 2022).

For a typical realization of our reference model, 35 WDMSs have evolved through common-envelope evolution. Hence, our simulations predict that ∼ 30% of the observed unresolved WDMS systems have gone through a common-envelope episode during their evolution. Thus, we expect a similar percentage of observed WDMSs to be post-common-envelope binaries. Moreover, in Fig. 7, we plot the final period versus distance for a typical realization of our reference model. Even at closer distances, we can observe the existence of wide binary systems that remain unresolved. Specifically, for systems with long periods, log [T (days)] > 6, resolution is particularly further limited by the inclination and eccentricity of the orbit. As a consequence, we can conclude that even for unresolved WDMS systems at close distances, it cannot be guaranteed that they are close binary systems.

Fig. 7 Distribution of the period versus the distance for the final simulated unresolved WDMS sample represented in Fig. 6. Also plotted are the nominal Gaia angular resolution limit (0.5″; dashed green line) and the more conservative value adopted in our analysis (2″; dashed red line).

4.3 The unbiased unresolved population of single and binary stars

The main purpose of this work has been to simulate the population of unresolved WDMSs accessible by Gaia within 100 pc. To that end, we restricted our simulated samples considering all possible observational biases affecting the population studied by Rebassa-Mansergas et al. (2021b), which decreased the original synthetic sample by ∼ 60%. Hence, the considered observed WDMS sample is expected to be ∼ 40% complete within the WDMS region. The results have shown that the synthetic population not only has nearly the same number of observed WDMS systems, but also that the two samples have similar stellar parameters. This implies that the theoretical model adopted in our code for reproducing the population of Galactic single and binary stars can be considered as reliable. Consequently, we can argue that the unfiltered synthetic sample should be reliable of the underlying stellar population.

Figure 8 shows the percentages of single and binary stars (top panel) and the percentages of the different type of binaries (bottom panel) for the unfiltered synthetic population, that is, assuming ideal conditions with no errors in flux, parallax, or excess factor, and without applying any final selection function. Interestingly, these percentages are rather similar to those derived from the restricted synthetic samples (see Fig. 2).

Regarding the unfiltered sample of unresolved WDMSs, the total number of synthetic systems is 384 ± 20 within the defined region of Rebassa-Mansergas et al. (2021b), which lowers the completeness of the observed sample farther from ∼ 40% to ∼ 26%. Furthermore, the overall WDMS population consists of around 1500 systems, representing only about ∼ 9% of those accessible to Gaia. As previously stated, most of these objects are dominated by the flux of the main-sequence companion in the optical. Consequently, they are situated in the region of the colormagnitude diagram representative of single main-sequence stars, making them challenging to identify. A comparison between the stellar parameters of the synthetic WDMSs with and without filtering the population is provided in Fig. 9 (blue and red histograms, respectively). Additionally, we estimate the completeness of each bin, i, as the ratio C_i = N_filtered/N_unfiltered. Although the mean completeness (∼ 26%) is relatively low, the distribution of parameters derived from the Gaia sample can be considered representative of the WDMS population within that selection region.

Fig. 8 Same as Fig. 2 but for ideal conditions, that is, without considering any observational limits. The complete dataset can be found in Table A.2.

5 Conclusions

During the last two decades, WDMSs have proven to be excellent tools for improving our understanding of a wide variety of open issues in astronomy, ranging from the evolution of close compact binaries to the chemical and dynamical evolution of the Galactic disk. However, due to major observational selection effects, only a small percentage of WDMSs are available for proper characterization and analysis. This is because most WDMS samples are magnitude-limited and, furthermore, the main-sequence components generally outshine the white dwarfs in the optical, which makes them difficult to identify unless ultraviolet data are used. This issue also affects recent volume-limited samples of unresolved WDMS systems from Gaia.

We used a detailed Monte Carlo code, MRBIN, developed by our team to simulate the population of both single and binary main-sequence stars, as well as white dwarfs, in our solar neighborhood. Special emphasis was placed on the analysis of the unresolved WDMS population, which has been filtered strictly following Gaia observational selection criteria. We compared the resulting synthetic samples with those obtained from the most homogeneous volume-limited catalog of WDMSs from Gaia (Rebassa-Mansergas et al. 2021b). This allowed us to test the models used in the simulations and evaluate the level of completeness of the observed sample.

Taking as a starting point the physical inputs obtained in Torres et al. (2022) in the fitting of the Gaia-resolved WDMS and DWD populations, our synthetic sample aligns well with the distribution of observed unresolved WDMSs in the Gaia color-magnitude diagram. This fact reinforces the validity of our modeling approach, which provides a comprehensive description of the different subpopulations of our solar neighborhood. In this sense, it is also important to emphasize that our simulator yields synthetic samples of single main-sequence stars and white dwarfs, as well as resolved and unresolved WDMSs and DWDs, whose overall numbers and fractions perfectly match the observed values. However, some discrepancies emerged when comparing the stellar parameters of the observed and synthetic unresolved WDMSs, notably regarding the mass distribution of white dwarfs. We find that the synthetic population under-represents the lower-mass end of the observed distribution (< 0.5 M_⊙). This is likely related to uncertainties in the mass determination of these specific systems in the observed sample.

By considering the synthetic WDMS sample without applying any observational cuts, we found that Gaia's observational constraints limit the detectable fraction of unresolved WDMSs to around 25% in the intermediate region between white dwarfs and main-sequence stars in the Gaia color-magnitude diagram. As mentioned before, this is not unexpected and makes it clear that the majority of unresolved WDMSs (∼ 75%) are located in the main-sequence locus of the Gaia color-magnitude diagram.

Our comparison between the observed and synthetic WDMS samples also places useful constraints on the common-envelope efficiency parameter (α_CE), with best-fit values between 0.1 and 0.4. We find an overall agreement with observations for a value of α_CE = 0.3, which confirms previous studies that suggested low values for α_CE. Additionally, our analysis indicates that ∼ 30% of the observed Gaia unresolved WDMS population has undergone a common-envelope episode.

This research strengthens the predictive power of population synthesis in stellar astrophysics and provides a foundation for interpreting future Gaia data releases and other white dwarf binary surveys.

Fig. 9 Mass (and effective temperature distribution (left and right panels, respectively) of the main-sequence (top panels) and white dwarf (bottom panels) components of the unresolved synthetic systems, with (blue) and without (red) applying the observational filtering. The blue curves (right Y-axis) represent the completeness of each bin, i, in the distribution, obtained as Ci = Nfiltered/Nunfiltered. Mean values from a total of 100 realizations of our reference model were used.

Acknowledgements

We acknowledge the positive and constructive suggestions and comments of our anonymous referee. This work was partially supported by the MINECO grant PID2020-117252GB-I00 and the PhD grant PRE2021 100503 funded by MICIU/AEI/10.13039/501100011033 and ESF+. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

Appendix A Stellar population fractions within 100 pc

In Tables A.1 and A.2 we present the number of objects and the percentages (shown in brackets) of the different subpopulation of our simulation of the 100 pc sample. Results correspond to a typical realization of our simulator, shown for cases with Gaia selection effects applied (Table A.1) and without them (Table A.2). Physical input parameters are those from Torres et al. (2022) and described in Sect. 3. The three last columns correspond to the total population (i.e., the complete color-magnitude diagram), the WDMS region (as defined in Rebassa-Mansergas et al. 2021a), and the white dwarf region (as delimited in Jiménez-Esteban et al. 2023).

Appendix B Kolmogorov-Smirnov tests

Fig. B.1 Cumulative distribution function (CDF) of the synthetic sample (blue line) and the observed Gaia sample from Fig. 1 (red line) for the magnitude MG (left panels) and the color GBP−GRP (right panels) under different selection criteria. For a more quantitative comparison, the distance (DIST; green line) between the CDFs is also plotted. As expected, the agreement between the observed and synthetic samples improves as the different selection criteria are applied.

Fig. B.2 KS tests for mass (top) and effective temperatures (bottom). In each panel, the top subpanels represents the cumulative distribution function (CDF) of the main-sequence component of the observed Gaia sample (Rebassa-Mansergas et al. 2021b, red line) and the synthetic sample (blue line), and the bottom subpanels display the distance (DIST; green line) between the CDFs.

Fig. B.3 As in Fig. B.2, but for the white dwarf component.

References

All Tables

All Figures

Fig. 1 Color-magnitude diagrams for the Gaia 100 pc sample. The sample construction process is exemplified with color-magnitude diagrams for the observational sample for: the different selection phases (left panels; see Table 1); the no-cuts samples (upper panels); the flux and parallax selection (upper-middle panels); the excess factor removal (lower-middle panels); and the final WDMS sample (bottom panels). The corresponding diagrams of the simulated sample are also shown (right panels). The WDMS region (Rebassa-Mansergas et al. 2021b) is marked with a solid blue line, and the dashed blue line depicts the magnitude limit of our simulations. Finally, the selection function (green line contours corresponding to 15%, 36%, and 49%) is derived from the observed sample from Rebassa-Mansergas et al. (2021 b, bottom-left panel, magenta dots; RM+21) and applied to a typical synthetic realization (bottom-right panel, magenta dots).

In the text

	Fig. 2 Percentages of white dwarfs and main-sequence stars (top chart) and those of the white dwarf subpopulations (bottom chart) within 100 pc. The complete dataset can be found in Table A.1.
In the text

	Fig. 3 Number of unresolved WDMS systems as a function of the αCE parameter for binary systems in the WDMS region (red), showing also how many have undergone a common-envelope episode (blue) and how many have not (green). The vertical dashed orange lines show the lower and upper limits of αCE compatible with the observed value (gray).
In the text

	Fig. 4 Number of eclipsing binaries as a function of the αCE parameter for binary systems in the delimited WDMS region. As in Fig. 3, the vertical dashed orange line represents the lower limit of αCE compatible with the observed value (gray).
In the text

	Fig. 5 Distribution of mass and Teff for the main-sequence stars (top panel) and white dwarf component (bottom panel) of synthetic (blue) and Gaia-observed (red) WDMS systems. For a better statistical comparison, 100 different synthetic realizations are plotted.
In the text

	Fig. 6 Initial (red) and final (blue) period distribution of simulated WDMS systems within the observed region.
In the text

	Fig. 7 Distribution of the period versus the distance for the final simulated unresolved WDMS sample represented in Fig. 6. Also plotted are the nominal Gaia angular resolution limit (0.5″; dashed green line) and the more conservative value adopted in our analysis (2″; dashed red line).
In the text

	Fig. 8 Same as Fig. 2 but for ideal conditions, that is, without considering any observational limits. The complete dataset can be found in Table A.2.
In the text

Fig. 9 Mass (and effective temperature distribution (left and right panels, respectively) of the main-sequence (top panels) and white dwarf (bottom panels) components of the unresolved synthetic systems, with (blue) and without (red) applying the observational filtering. The blue curves (right Y-axis) represent the completeness of each bin, i, in the distribution, obtained as Ci = Nfiltered/Nunfiltered. Mean values from a total of 100 realizations of our reference model were used.

In the text

Fig. B.1 Cumulative distribution function (CDF) of the synthetic sample (blue line) and the observed Gaia sample from Fig. 1 (red line) for the magnitude MG (left panels) and the color GBP−GRP (right panels) under different selection criteria. For a more quantitative comparison, the distance (DIST; green line) between the CDFs is also plotted. As expected, the agreement between the observed and synthetic samples improves as the different selection criteria are applied.

In the text

	Fig. B.2 KS tests for mass (top) and effective temperatures (bottom). In each panel, the top subpanels represents the cumulative distribution function (CDF) of the main-sequence component of the observed Gaia sample (Rebassa-Mansergas et al. 2021b, red line) and the synthetic sample (blue line), and the bottom subpanels display the distance (DIST; green line) between the CDFs.
In the text

	Fig. B.3 As in Fig. B.2, but for the white dwarf component.
In the text