Physical Properties of 29 sdB+dM Eclipsing Binaries in Zwicky Transient Facility

The hot subdwarf stars are located between the main sequence and the white dwarf (WD) sequence in the Hertzsprung–Russell diagram (Heber 2016). They are classified as B-type (sdB, $T_{\rm eff}<40000$ K) and O-type (sdO, $T_{\rm eff}>40000$ K) according to their different spectral types (Hirsch et al. 2008). sdB stars are located in the blue tail of the horizontal branch (HB), which is also known as the extreme horizontal branch (EHB; Heber 1986). sdB has a Helium-burning core and an extremely thin residual Hydrogen envelope (mass $<$ 0.01 $M_{\odot}$). Its progenitor is likely an evolved star that loses most of its Hydrogen envelope when it ascends the red giant branch. sdB will not enter the asymptotic giant branch phase, instead, it will become a white dwarf directly. The average lifetime of EHB is in the order of $10^{8}$ yr (Dorman et al. 1993). Maxted et al. (2001) suggested the binary ratio of sdB is $\sim$ 60$\%$. About half of these binary systems are close binaries with periods of a few days (Allard et al. 1994; Maxted et al. 2001), while the others are wider binaries with periods of several years (Stark & Wade 2003; Vos et al. 2018; Chen et al. 2013).

The evolution of sdB is still controversial (Geier & Heber 2012), although the binary evolution scenario is favored. Han et al. (2002, 2003) proposed three binary evolution channels for sdBs (common envelope evolution, stable Roche-lobe overflow, and the merger of two He white dwarfs). Subsequently, a number of works tried to testify these theories through observations (Luo et al. 2019, 2020; Kupfer et al. 2020; Kramer et al. 2020; Pelisoli et al. 2020). For single sdB stars, the formation has always puzzled us. About $40\%$ of the sdB stars have been found to be single. Although it has been argued that all sdB stars were evolved from binary stars (Pelisoli et al. 2020), but several facts suggest that the merger channel of two He WDs can not satisfactorily explain the single sdB stars (Burdge et al. 2020; Van Grootel et al. 2021). First, the number of low-mass white dwarf binaries is small (Ratzloff et al. 2019b). Second, the predicted broad mass distribution of single sdB stars from the merger channel seems to be inconsistent with observations (Fontaine et al. 2012). Third, observations show that most of the single sdB stars have a very slow rotation, which against the fast rotation predicted by the merger channel (Geier & Heber 2012; Charpinet et al. 2018; Geier & Heber 2012). Recently, a new channel for the formation of single sdB stars has been proposed (Meng et al. 2020).

When an sdB+dM binary system is eclipsed, it exhibits an HW Vir-type light curve (LC). The number of such stars is relatively small, and since they possesses a characteristic LC, they provide a direct way to measure the physical properties of components by modeling LC and the radial velocity curve, which helps to test the theory of formation and evolution of such systems. For a more detailed review on HW Vir-type stars, see Heber (2016). Such stars have been studied continuously and extensively over the last few decades, and their number is increasing (see e.g. Kilkenny et al. (1978); Menzies & Marang (1986); Drechsel et al. (2001); Østensen et al. (2007); For et al. (2010); Barlow et al. (2013); Schaffenroth et al. (2014); Kupfer et al. (2015); Almeida et al. (2017); Ratzloff et al. (2019a); Koen (2019); Ratzloff et al. (2020); Sahoo et al. (2020a); Baran et al. (2021)). The orbital periods of HW Vir-type stars are easily obtained from LCs, while the physical parameters such as mass and radius of sdB and the companion stars can be determined based on radial velocity measurements and LC analysis. Based on these parameters, it is possible to constrain the evolution of HW Vir-type stars and the mass–radius relations of two components. The mass–radius relation of brown dwarfs was studied by Baraffe et al. (2003) using the COND model from Chabrier et al. (2000). The mass–radius relation of low-mass stars was studied by Knigge et al. (2011) using the cataclysmic variables. The HW Vir-type stars are suitable objects for updating the physical properties of low-mass stars and brown dwarfs.

There are currently about 200 HW Vir-type stars in the database of AAVSO International Variable Star Index¹¹1https://www.aavso.org/vsx/index.php (Watson et al. 2006, VSX), which are mainly collected from the periodic variable catalog of Zwicky Transient Facility (Chen et al. 2020, ZTF), the fourth phase of the Optical Gravitational Lensing Experiment (Udalski et al. 2015, OGLE), the catalog of the Asteroid Terrestrial-impact Last Alert System (Heinze et al. 2018, ATLAS). Schaffenroth et al. (2019) studied a sample of HW Vir-type stars in OGLE and ATLAS. ZTF is a northern-sky time-domain survey with limiting magnitudes of 20.8 mag in $g$ band and 20.6 mag in $r$ band (Masci et al. 2019). By monitoring the northern sky with several hundreds of detections for each object over a 3-year period, ZTF will discover a large number of HW Vir-type stars. The presence of more short-period ($P<0.2$ days) HW Vir-type stars will help constrain the evolutionary timescale and ending of such systems.

In this paper, we analyze the physical and orbital parameters of 29 sdB+dM eclipsing binaries using LCs from ZTF data release 5. By fixing sdB’s mass and temperature, we obtain the orbital inclination, the mass ratio, the radius of sdB and dM, and the dM temperature. We also determine the mass–radius relation of dM and investigate how sdB+dM evolves with the shortening of the orbital period. In Section 2, we describe the selection of sdB+dM eclipsing binaries and their LCs. We explain how to use PHOEBE to obtain orbital parameter solutions from LCs, and the reliability of our determined parameters in Section 3. In Section 4, we present our results and make a comparison with previous works. We also discuss the physical properties and the evolution stage of our sdB+dM eclipsing binaries in this section. Finally, we conclude this work in Section 5.

We cross-matched the HW Vir-type stars in VSX database with the ZTF catalog of periodic variable stars and obtained a sample of 31 HW Vir-type stars. Two other HW Vir-type stars not included in the VSX were identified by eye and added to the sample. About two-thirds of this sample have periods smaller than 0.2 days. We compared them to other $\sim 900$ short-period variables from ZTF periodic variable catalog through the $\sl Gaia$-band intrinsic color vs. absolute magnitude diagram (CMD). This helps determine the type of primary component of HW Vir-type stars. In Figure 1, the red dots represent the 33 HW Vir-type stars, and the gray dots are the ZTF short-period variables with $P<0.2$ days. We adopted the extinction values of Green’s 3D extinction map (Green et al. 2019) and converted them to $A_{G}$ and $E(G_{\rm Bp}-G_{\rm Rp})$ using the updated extinction law (Wang & Chen 2019). The absolute magnitude and intrinsic color are estimated by $M_{G}=G-A_{G}-5\log(1/\varpi)-10$ and $(G_{\rm Bp}-G_{\rm Rp})_{0}=(G_{\rm Bp}-G_{\rm Rp})-E(G_{\rm Bp}-G_{\rm Rp})$, respectively. $\varpi$ is the corrected Gaia early data release 3 (EDR3) parallax in $m$as (Lindegren et al. 2021). The corrected EDR3 parallax was found to reduce the systematic bias to $<10$ $\mu$as (Ren et al. 2021).

From Figure 1, we found that HW Vir-type stars concentrate as a clump with an absolute magnitude around $M_{G}=5$ mag. The clump is bluer than the main-sequence variables and brighter than the cataclysmic variables. This means our HW Vir-type stars are sdB+dM eclipsing binaries rather than WD+dM eclipsing binaries. To further confirm that primary components our sample are all sdB stars, we adopt the method in (Gentile Fusillo et al. 2015), where the authors calculated the reduced proper motion defined as $H$ = $G$ + 5 log$\mu$ + 5. $G$ and $\mu$ are $G$-band magnitude and proper motion. We found that $H$ of our sample are all less than 14.5, which proves that they are sdBs rather than WDs (Schaffenroth et al. 2019). The distribution of HW Vir-type stars is consistent with constraints of $-1.0<M_{G}<7.0$ and $(G_{BP}-G_{RP})_{0}<0.7$ (Geier et al. 2019; Geier 2020), but more concentrated in $M_{G}$. Four sdB+dM are brighter than the others, and these deviations are due to the large uncertainties propagated from Gaia parallaxes. We showed their error bars in Figure 1. The $1\sigma$ uncertainties of the absolute magnitude for other sdB+dM are less than 0.35 mag. We also noted that in CMD, there are dozens of gray dots distributed in the location of our sdB+dM eclipsing binaries. By examining their LCs, we found that they are non-eclipsing sdB+dM binaries, and their brightness variations are due to reflection. Reflection occurs when the hot component illuminates part of the spherical surface of the cooler component when the temperature difference between the two components is large. The size of the illuminated surface varies with the projection angle, which leads to an eventual sinusoidal-like light curve without eclipse. This suggests that many sdB binaries, which orbit at such a small inclination that we cannot see the eclipse. The HW Vir-type star is a subtype of detached eclipsing binaries (EA-type) that have similar LC shapes. The specific distribution in the intrinsic color vs. absolute magnitude diagram helps to separate HW Vir-type stars from other main-sequence eclipsing binaries. With these properties, we expect to find hundreds of HW Vir-type stars in future ZTF-based variable star searches. Besides, CMD is useful to separate sdB+dM eclipsing binaries and WD+dM eclipsing binaries when Gaia parallaxes are accurate enough.

We downloaded the $g,r$-band LCs of 33 sdB+dM eclipsing binaries from the ZTF DR5 database. To ensure the quality of LCs, we selected only good photometric data with ‘catflag $<$ 10’. ZTF DR5 contains photometric data obtained during a two-and-a-half year survey, with more than 200 detections per band for most objects. Given that typical photometric uncertainties are around 0.015 mag, we believed that a reliable analysis can be performed on the ZTF LCs. LCs were folded with periods from Chen et al. (2020), and the ephemeris of the primary minimum $T_{0}$ was estimated by lc_geometry package in PHOEBE 2.3 (Conroy et al. 2020). We also performed a careful visual inspection of all light curves to make sure that the period and $T_{0}$ were not significantly biased. We excluded photometric outliers (less than 10) for each LC according to the fit line of the Gaussian processes. Magnitude was converted into normalized flux by assuming the maximum flux in $g,r$ bands is equal to 1. The primary eclipse of four sdB+dM is fainter than the detection limit of ZTF ($r>20.5$ mag). We excluded them and performed LC analysis on the remaining 29 sdB+dM.

sdB+dM eclipsing binaries are the best objects to study the physical properties of sdB and dM. However, only the radial velocity curve of the primary component is available, so the mass ratio cannot be directly determined from radial velocity curves. One way to solve this problem is to assume a typical mass of sdB, e.g., $M_{1}=0.47M_{\odot}$. Based on this assumption, parameters of two components can be well constrained by the analysis of the single-line radial velocity curve and light curves. The shortcoming of this method is that $M_{1}$ cannot be studied. The other way is to constrain the mass ratio $q$ based on $\log g$ and fix it in the subsequent analysis (For et al. 2010; Schaffenroth et al. 2021; Drechsel et al. 2001; Heber et al. 2004; Vučković et al. 2007; Schaffenroth et al. 2014; Almeida et al. 2017; Barlow et al. 2013). In this method, the accuracy of $q$ depends on the extent to which it is constrained by LCs. When $q$ does not converge in the iterations, the uncertainty is relatively larger. In this work, we used the PHOEBE 2.3 code to study the probability distribution of parameters for 29 sdB+dM and to evaluate how well $q$ is constrained by LCs.

With some exceptions, most sdB stars have a mass around 0.47 $M_{\odot}$, which is also known as the typical mass (Han et al. 2002, 2003; Sahoo et al. 2020b). Most sdB stars have a temperature in range of $27,000-31,000$ K (Geier et al. 2012). 3 of 29 sdB+dM were studied previously with spectra and sdB’s temperature is in this range. The three sdB+dM are ZTF J162256.66+473051.1 with a temperature $T_{1,\rm eff}=29000\pm 600$ K (Schaffenroth et al. 2014), ZTF J153349.44+375927.8 with a temperature $T_{1,\rm eff}=29230\pm 125$ K (For et al. 2010), and ZTF J223421.49+245657.1 with a temperature $T_{1,\rm eff}=28500\pm 500$ K (Almeida et al. 2017) or $T_{1,\rm eff}=28370\pm 80$ K (Østensen et al. 2007). The companion star of sdB+dM is an M dwarf or a brown dwarf, and its temperature is around 3000 K (Menzies & Marang 1986; Schaffenroth et al. 2014; Vučković et al. 2014).

Our idea is to assume that the primary components of sdB+dM are all canonical sdB star with a mass of 0.47 $\pm\ 0.03$ $M_{\odot}$, a radius of 0.175 $\pm\ 0.025$ $R_{\odot}$, and temperature about 30,000 K (Van Grootel et al. 2021). We fixed the sdB mass $M_{1}$ to 0.47 $M_{\odot}$, the sdB temperature $T_{1}$ to 30,000 K in our LC analysis. Small deviations in $T_{1}$ will directly affect the poorly constrained $T_{2}$ (typically with $30\%$ uncertainty), but have little effect on $R_{2}$ and other parameters. With these assumptions, we can still study the properties of the companion star and the evolution stage of the system in statistics. We set the orbital eccentricity $e$ to 0 since all of our sdB+dM eclipsing binaries have a very short period (less than 0.5 days). We set the linear limb-darkening coefficient of sdB to 0.25 and 0.2 in the $g$ band and $r$ band, respectively, while the companion is 1.0 in both the $g$ band and $r$ band (Claret & Bloemen 2011). we set the gravity-darkening factor to 1.0 and 0.32 for sdB and the companion, respectively (Claret & Bloemen 2011). The reflection factor of sdB and the companion was set to 1.0. The range of the initial orbital inclination $i$ is $60-90$ degrees, and the other parameters are unrestricted. All preset parameters are list in Table 1.

We used the Markov chain Monte Carlo (MCMC) sampler based on EMCEE (Foreman-Mackey et al. 2013) to determine parameters of 29 sdB+dM eclipsing binaries. For sampling, we used 50 walkers and 2000 iterations, and each run costs 30 hours on a 4-core CPU. We also tested 5000 iterations for 5 sdB+dM and 20,000 iterations for one sdB+dM and found the results were almost the same as those based on 2000 iterations. In the first run, the mass ratios $q_{i}$ are not well constrained for most of the objects, as it is difficult to obtain the true $q$ from many local optimal mass ratios. In contrast, the radii $R_{1,i},R_{2,i}$ and inclination $i_{i}$ are well constrained in two-band LC analysis. From the $R_{1,i}-q$ and $R_{2,i}-q$ distributions, $R_{1,i},R_{2,i}$ only increase by $20\%$ when $q$ increases from 0 to 1. So we decided to use $R_{2,i}$ to establish a prior distribution of $q$. According to the mass–radius relationship of dM, the ratio of mass to radius is slightly less than 1 in the solar unit. Given that $q$ and $R_{2}$ are usually overestimated in the first run ($q$ of the sdB+dM is usually less than 0.3), we assumed a prior distribution of $q$ as [0, $R_{2,i}$/0.47]. We again performed MCMC estimation and obtained the reliable parameter distributions for most of sdB+dM. 11 sdB+dM eclipsing binaries have temperatures $T_{2}>>3000$ K are not consistent with their masses and radii. For these objects, we added another prior $T_{2}<4000$ K and performed the MCMC estimation again to obtain the final distributions of parameters.

Figure 3 shows the parameter distributions of ZTF J204046.56+340702.8 as an example of one class. We can see that $i$, $q$, $R_{1}$ and $R_{2}$ obey Gaussian-like distributions and the correlations between these parameters are small. Parameters of the other three sdB+dM eclipsing binaries are similarly distributed (See Figures in Appendix A). Figure 3 shows the parameter distributions of ZTF J224547.36+490824.7. It differs from Figure 3 in that $q$ does not satisfy the Gaussian distribution. For this class, $q$ is uniformly distributed, or has a wide tail, or is biased towards one side in the prior range (See Figures in Appendix B). Considering that the prior range is narrow, i.e., $q_{i}$ is roughly distributed between 0 and $0.2-0.5$, the deviation between the maximum likelihood and the true $q$ is not large. This deviation is already included in the error of $q$ (mean uncertainty is 0.08). Besides, the error in $R_{2}$ hardly increases with the error in $q$ due to the weak correlation between $R_{2}$ and $q$. Based on multi-band LC analysis, the determined masses and radii of the companions in sdB+dM eclipsing binaries are reliable for studying their statistical properties. Nevertheless, in most cases, uncertainty in $q$ or $M_{2}$ is larger.

In both Figure 3 and Figure 3, $T_{2}$ does not show a Gaussian distribution, which means that $T_{2}$ is difficult to be constrained in the light curve analysis. The typical uncertainty of $1500$ K also indicates that $T_{2}$ is rather uncertain. The upper panels of Figures 5 and 5 show a comparison of the modeled and observed light curves of ZTF J204046.56+340702.8 and ZTF J224547.36+490824.7, while the lower panels show the residual diagrams. Diagrams for other sdB+dM eclipsing binaries are supplemented in Appendix A and B.

In this section, we present a table including parameters of 29 sdB+dM eclipsing binaries and make a comparison with previous works. We also discuss how to use our result to constrain the evolutionary stage of sdB+dM eclipsing binaries and the physical properties of companion dM stars.

In this work, we have selected a sample of 33 HW Vir-type stars in ZTF DR5. Based on Gaia EDR3 parallax and extinction correction, we found that these HW Vir-type stars are concentrated in a clump in the intrinsic color vs. absolute magnitude diagram. Their locations in the CMD imply that they are sdB+dM eclipsing binaries. By fixing the mass and temperature of sdB to $M_{1}=0.47M_{\odot}$ and $T_{1}=30,000$ K and setting the prior of $q$ and $T_{2}$, we analyzed LCs using PHOEBE 2.3 code to obtain the probability distributions of parameters $q$, $R_{1}$, $R_{2}$, $i$, and $T_{2}$. We obtained LC solutions for 29 sdB+dM eclipsing binaries with full primary eclipse detection. $R_{1}$, $R_{2}$, and $i$ obey a Gaussian-like distribution and have little correlation with other parameters, which means that they are well constrained by the LC analysis. $q$ does not show a Gaussian distribution in most cases, and the mean uncertainty is 0.08.

Our parameters for three sdB+dM are comparable with previous works that used both LCs and the single-line radial velocity curve if the mass of sdB is $M_{1}=0.47M_{\odot}$. This means that parameters of sdB+dM determined from LCs are suitable for statistical analysis. Based on 29 sdB+dM, we found that both $q$ and $R_{2}$ decrease with the decreasing of the orbital period. sdB+dM with larger mass companions are likely to merge early during the shortening of the orbit. It is worth mentioning that companions of all three sdB+dM are brown dwarfs when the orbital period is less than 0.075 days. The masses and radii of the companion stars are consistent with the mass–radius relation for low-mass stars and brown dwarfs. We found a gap between $R=0.12R_{\odot}$ and $R=0.15R_{\odot}$ which can be explained as the boundary between low-mass stars and brown dwarfs.

sdB+dM eclipsing binaries are important objects to study the nature of sdB and dM, and their evolutionary endings are very interesting. With deeper photometry and more detections, ZTF will detect and classify hundreds of short-period sdB+dM eclipsing binaries. Before the availability of large-scale and deep spectroscopic surveys, the statistical properties of sdB and dM can be obtained from the LC analysis of large samples.