Unexpected Short-Period Variability in Dwarf Carbon Stars from the Zwicky Transient Facility

Dwarf carbon (dC) stars, main sequence stars showing carbon molecular bands, are enriched by mass transfer from a previous asymptotic-giant-branch (AGB) companion, which has since evolved to a white dwarf. While previous studies have found radial-velocity variations for large samples of dCs, there are still relatively few dC orbital periods in the literature and no dC eclipsing binaries have yet been found. Here, we analyze photometric light curves from DR5 of the Zwicky Transient Facility for a sample of 944 dC stars. From these light curves, we identify 34 periodically variable dC stars. Remarkably, of the periodic dCs, 82\% have periods less than two days. We also provide spectroscopic follow-up for four of these periodic systems, measuring radial velocity variations in three of them. Short-period dCs are almost certainly post-common-envelope binary systems, since the periodicity is most likely related to the orbital period, with tidally locked rotation and photometric modulation on the dC either from spots or from ellipsoidal variations. We discuss evolutionary scenarios that these binaries may have taken to accrete sufficient C-rich material while avoiding truncation of the thermally pulsing AGB phase needed to provide such material in the first place. We compare these dCs to common-envelope models to show that dC stars probably cannot accrete enough C-rich material during the common-envelope phase, suggesting another mechanism like wind-Roche lobe overflow is necessary. The periodic dCs in this paper represent a prime sample for spectroscopic follow-up and for comparison to future models of wind-Roche lobe overflow mass transfer.


INTRODUCTION
Carbon (C) stars are those that show molecular absorption bands of C, such as C 2 , CN and CH, in their optical spectra (Secchi 1869). Intrinsic C stars have experienced C enrichment via dredge-up of fusion products from their cores. During the thermally pulsing (TP) phase of the asymptotic giant branch (AGB), shell He flashes cause strong convection in the shell regions bringing C into the atmosphere -the third dredge-up. mer AGB companion would since have become a white dwarf (WD) and cooled until no longer detectable alongside the dC. Indeed this seemed to have been confirmed when Dearborn et al. (1986) found G77-61 to be a radial velocity (RV) binary with an orbital period of 245.5 d.
Today, G77-61 is no longer the only known dC, with close to 1000 known in the literature. The majority of these dCs come from the Green (2013) and Si et al. (2014) samples which were found from all-sky spectroscopic surveys. This has included almost a dozen "smoking gun" systems in which the WD companion is sufficiently hot to be visible in the optical spectra as a spectroscopic composite (Heber et al. 1993;Liebert et al. 1994;Green 2013;Si et al. 2014). These samples have shown that the dC stars are actually the most common type of C star in the Galaxy.
Their carbon-enriched atmospheres make dCs the most likely progenitors of the carbon-enhanced metalpoor, CH, and possibly the Ba II stars, all showing carbon and s-process enhancements (Jorissen et al. 2016;De Marco & Izzard 2017). These stars, being more luminous than dCs, have been studied via RV campaigns, which have shown increased binarity compared to normal O-rich stars, indicating they have likely also experienced mass transfer from an unseen companion (Sperauskas et al. 2016). Barium dwarfs and CH subgiants show periods from RV analysis of 1-20 years (Escorza et al. 2019).
Blue straggler stars are another class similar to dCs in that they may have experienced mass transfer from a previous AGB companion. As discussed by Gosnell et al. (2019), blue straggler stars in a cluster color-magnitude diagram are more luminous and bluer than the main sequence turnoff. While some are likely formed in mergers or collisions, most blue straggler stars are thought, like dC stars, to be the result of mass transfer from a giant to a main sequence dwarf. Most blue straggler stars are found in wide binaries with periods of order 1000 days, consistent with expectations of mass transfer from an AGB star onto a main-sequence companion (Chen & Han 2008;Gosnell et al. 2014), which leaves a CO-core WD remnant (Paczyński 1971). Those blue straggler stars that form after mass transfer from an red giant branch star, yields a blue straggler in a shorter binary period (of order 100 days; Chen & Han 2008) leaving a He-core WD companion. The salient point relevant to dC stars is that significant mass is typically gained in these encounters.
While dC stars are now known to be numerous, details of their properties remains sparse. This is especially true of their orbital properties. Currently, there are only six orbital periods for dCs in the literature.
The first is of the dC prototype G77-61, found to be a single-line spectroscopic binary with an orbital period of 245.5 d (Dearborn et al. 1986). The central source of the Necklace Nebula was found to be a binary with a dC, which has a photometric period of 1.16 d (Corradi et al. 2011;Miszalski et al. 2013). The three longest period dCs in the literature are those from Harris et al. (2018) who found astrometric binaries with periods of 1.23 yr, 3.21 yr, and 11.35 yr. Margon et al. (2018) found a dC with a photometric period of 2.92 d and confirmed this as the orbital period with spectroscopic follow-up.
In their recent work, Whitehouse et al. (2021) conducted an RV survey of seven dCs with Hα emission, finding short orbital periods for all of them (six new periods). In addition, they found photometric periods with similar lengths as the orbital periods in the range of 0.2-5.2 d. Their light curve modeling suggests that the source of variability in their dCs is from stellar rotation and spots. As with the new photometrically periodic dCs in this paper, these dCs must have experienced a common-envelope phase with the former AGB companion.
There have also been large sample few-epoch spectroscopy campaigns of dCs. Whitehouse et al. (2018) conducted a few-epoch survey of 28 dCs finding RV variability in 21 of them, implying a high binary fraction. Roulston et al. (2019) conducted a larger survey of 240 dCs using few-epoch spectroscopy from the Sloan Digital Sky Survey (SDSS; Blanton et al. 2017). They found that dCs are consistent with 100% binarity, with separations of order 1 au and periods of order 1 yr. Both the Whitehouse et al. (2018) and Roulston et al. (2019) surveys lacked enough spectral epochs to fit individual orbits and instead relied on statistical analysis to describe the dC population as a whole.
de Kool & Green (1995) modeled the space density of dCs, and they predicted the dC period to be bimodal with peaks near 10 2 -10 3 d and 10 3 -10 5 d, consistent with the known periods listed above. de Kool & Green (1995) also found that the production of dCs is strongly dependent on metallicity, finding no dCs should be formed in systems with initial metallicity greater than half of solar (i.e. [Fe/H]> −0.3). This is in agreement with metallicity measurements of G77-61, where Plez & Cohen (2005) found [Fe/H] = −4.0, making G77-61 one of the lowest metallicity stars known. This also has been supported by Farihi et al. (2018) who found 30-60% of dCs to be halo objects, which are metal poor.
Until this year, the known dC periods spanned from ∼ 1 d to ∼ 4100 d, likely indicating different formation pathways. The longest dC periods that are of order 10s of years likely experienced only standard Roche-lobe overflow (RLOF) or wind-RLOF (WRLOF). These periods are consistent with other types of post-mass-transfer systems, such as the blue straggler stars.
The dCs with periods 1 d would have likely experienced common-envelope (CE) evolution (Paczynski 1976;Ivanova et al. 2013), since the TP-AGB envelope expands to several 100s of solar radii. Of interest is how CE evolution affects dC formation. Once the CE phase has started, the plunge-in of the lower-mass companion (in our case, the future dC) could truncate the evolution of the AGB by ejection of its envelope. If this happens before the TP-AGB phase and the third dredge-up, it is likely that the C enhancement needed for dC formation will not occur. However, if the CE begins after the AGB companion has already become a C-giant, then it may be possible for the main-sequence companion to accrete enough C-rich material from the CE to become a dC (depending on the accretion efficiency). If the accretion efficiency is not high enough, however, the main-sequence companion will not accrete enough material from the CE alone, requiring some combination of CE evolution with efficient mass transfer before the CE phase that is sufficient to transform an O-rich main-sequence star into a C-rich dC.
Significant accretion of mass and angular momentum from the AGB companion could result in significant spin-up and subsequent activity in some dCs . If there are dCs left in very tight orbits with the WD remnant, they may show tidally locked rotation periods (synchronous rotation), as well as tidal distortions causing ellipsoidal variations in photometric light curves. A search for periodicity in photometrically variable dCs could reveal some systems useful for constraining their evolution.
Another motivation to study variability in dCs is that no dC masses have yet been measured, because there are no known eclipsing dC systems. We can estimate dC masses from optical or infrared (IR) colors (see Section 6.2), but these estimates have uncharacterized systematics, due to differences between normal O-rich stars and C stars in the optical and IR regions.
This lack of eclipsing dC systems highlights the importance of photometric surveys to search for the first well-characterized eclipsing dC systems. These systems could, when combined with RV follow-up, provide the first reliable dC mass measurements, and help us understand more about the amount, and composition, of accreted mass needed to form a dC. Margon et al. (2018) conducted a search for periodic dCs using the Palomar Transient Factory (PTF; Law et al. 2009;Rau et al. 2009), finding just one periodic dC. However, they clearly highlighted the potential for large photometric surveys to find periodic dCs, particularly dCs with short periods that should have experienced the strongest phases of CE mass transfer.
In this paper we report on a unique sample of close binary dCs -implicating them as post-common envelope binaries (PCEBs) and likely pre-CVs -discovered from their periodic photometric variability in the Zwicky Transient Facility. In Section 2 we describe the sample of dCs we selected to search and in Section 3 we detail our process for cleaning and preparing the raw light curves. In Section 4 we describe our process for finding which dCs have detected periodic signals. In Section 5 we present spectroscopic follow-up for four of the periodic dCs in this paper. Finally, in Section 6 we present comparisons of these short period dCs to binary population synthesis models to understand how a commonenvelope phase relates to dC formation.

SAMPLE SELECTION
To search for variability in as many dCs as possible, we compiled a list of all dCs from the current literature. The largest contributor (747 dCs, 79%) is the Green (2013) sample of carbon stars from the SDSS. We also selected a smaller number of dCs from Si et al. (2014), who found 96 new dCs using a label propagation algorithm from SDSS DR8, and from Li et al. (2018) who selected carbon stars from the Large Sky Area Multi-Object Fiber Spectroscopic Telescope survey (LAMOST; Cui et al. 2012) using a machine learning approach. Our resulting final sample consists of 944 dCs.
With our compiled sample, to ensure that any periodic candidate was indeed a dwarf carbon star, we used Gaia EDR3 parallaxes, proper motions (Gaia Collaboration et al. 2021) and distances (Bailer-Jones et al. 2021). We required that each periodic C star had M G > 5 mag from Gaia EDR3 based either on (1) significant parallax / err > 5 (27/34 periodic dCs) or (2) a significant proper motion (µ/σ µ > 5) which sets an upper limit on the dC distance by limiting its transverse velocity to be less than an assumed Galactic escape velocity (Smith et al. 2007) of about 600 km/s (7/34 periodic dCs).

LIGHT CURVE PROCESSING
Using our list of dCs, we cross-matched our sample to the Zwicky Transient Facility DR5 (ZTF; Bellm et al. 2019;Masci et al. 2019;Graham et al. 2019). We required a match to be within 2 of our target coordinates and each star having ≥ 10 epochs in the available ZTF filters.
From the resulting matches detected within each filter, we grouped all sources within the match distance to We checked for any epochs which appear to be discrepant by performing an outlier removal on all the light curves. We first select from the raw light curve the brightest and faintest 5% of epochs. Within these brightest and faintest 5%, we calculate the median magnitude of each (i.e. the median of the 5% brightest and 5% faintest) and the mean error of that same brightest and faintest 5%. We then removed any outliers that were 2σ brighter than the median of the brightest 5%, and removed those 2σ fainter than the median of the faintest 5%. If this selection dropped the number of epochs below 10, we removed that light curve from our analysis. This treatment rejects most artifacts without removing genuine astrophysical variability.
We checked the light curves for each dC, in each filter, to determine if each dC had detected variability by examining how the mean magnitude changed across the observed light curve time span. A small number of dCs which show no periodic variability in our analysis in Section 4 (and a few periodic dCs) show signs of nonperiodic variability, as well as secular, long-term trends. These non-periodic but variable dCs are of interest and may be signs of flaring, variable obscuration, or perhaps accretion onto the WD companion. They warrant further investigation, but we do not discuss them further in this paper.
The light curves that show long-term trends of brightening or dimming on 100s of days timescales cause the mean magnitude to vary over the entire time span of the light curve. This variable mean magnitude can cause issues with our period search. Therefore, we removed these long term trends by fitting out a third-order polynomial to the raw light curve.

PERIODIC VARIABILITY
For each light curve, we searched for periodic signals down to a minimum period of 0.1 d using the the Lomb-Scargle periodogram (LS;Lomb 1976;Scargle 1982). We used the Astropy (Astropy Collaboration et al. 2018) implementation of the LS algorithm (VanderPlas et al. 2012;VanderPlas & Ivezić 2015). We selected the highest peak, and if this peak corresponds to an observational alias (1d, 29.5d, 1yr, etc.) or a harmonic of one of these aliases (1/2, 1/3, 1/4, 1/5, 2, 3, 4, 5), we removed that signal from the light curve and recalculated the periodogram until the highest-power frequency was not an alias (we counted a frequency not as an alias if it was more than 150 frequency bins away from the pure alias frequency, i.e. more than 0.005 d −1 away from an alias).
For the highest remaining peak, we calculated the false-alarm probability (VanderPlas 2018). We required that log (FAP) ≤ −5 in at least one filter for us to select a specific dC as a periodic candidate, more conservative than e.g., the log (FAP) ≤ −3 used in the recent ZTF periodic variable catalog of Chen et al. (2020).
For the dCs which have light curves selected as periodic candidates, we checked for any possible harmonic confusion in the found period. For each dC, in each filter, we plot a power spectrum from the LS analysis. This is used to determine how strong the highest-power frequency is compared to the log (FAP) limit and the background peaks. Figure 1 shows an example power spectrum for an object with a very strong periodic signal and shows clear peaks (with 1-d aliasing) above the background, and the resulting phased light curve. The complete figure set (90 figures) is available in the online journal. In some cases, the strongest peaks were aliases, typically harmonics of 1 month, that overwhelmed the power spectrum. For these dCs, we inspected each power spectrum in conjunction with phased light curves. If another non-alias peak (i.e., with a frequency more than 0.005 d −1 away from an alias) was found in the power spectrum meeting our FAP limit, that new peak was selected as the period for that dC. If no non-alias peaks could be found, the dC candidate was removed from our sample.  .9. This dC shows a clear and strong signal at 3.3 d −1 (with 1-d aliasing) that stands out above the low noise background in the power spectrum. The red horizontal line represents the peak height needed for a signal to meet our log (FAP) ≤ −5 criterion. Grey vertical dashed lines mark the 1-d aliases caused by gaps in data collection, and the best period from our analysis is marked by an arrow. The highest significance peak is used to fold the observed light curve, yielding the phase-folded light curve in the top panel. Each light curve is plotted twice in phase to clearly show the periodic variability. The data are shown as the black scatter points with their respective error bars, and the best fitting model (see Section 4) is marked by the red solid line. The residuals are shown below the light curve. The complete figure set (35 figures) is available in the online journal.
Some dCs show strong periodic signals in one filter, but do not reach our FAP limit in the other available filters. For these dCs, if one filter has a period that meets our FAP limit and that period is visible in the other filter, we include that second filter even if its FAP does not meet our limit. This makes it possible for some dCs to have a log (FAP) ≥ −5 in a filter if they have log (FAP) ≤ −5 in another filter.
For all periodic dC candidates selected after inspection of their power spectra, we plotted phased light curves folded on the highest selected peak period. In addition, we plotted 8 different harmonics of that period (1/2, 1/3, 1/4, 1/5, 2, 3, 4, 5) to check for aliases caused by gaps in the observational coverage. Using this plot, we calculated model-fit statistics (χ 2 ) and selected which period harmonic has the best model fit. We used the best period to phase the light curve, to which we fit the final periodic model.
Our best-fit models were computed using the automatic Fourier decomposition (AFD) method, as detailed in Torrealba et al. (2015). We set an upper limit on the number of Fourier series terms of n max = 6 to reduce over-fitting. No significant nonharmonic terms were included; though one dC, LAM-OST J062558.33+023019.4, showed different peaks in its power spectrum between the g and r filters with the second highest peak in each filter being the highest peak in the other. The best AFD model was used to calculate the amplitude and epoch of brightest time (t 0 ) for each light curve. We removed any dC for which the folded light curve shows no clear periodic signal or for which the amplitude of the variability was less than 0.005 mag. For dCs with multiple filters, we use the period from the filter with the strongest detection as the selected period and force the fits in the other filters to this fixed period. Some filters may not have a clear detection from the signal found in another filter, resulting in model fits with large errors for that filter. Figure Set 1 contains the folded light curves, models, and power spectra for all the periodic dC candidates. Table 2 contains the properties for this final periodic dC sample. We estimated errors for the best found period using a Markov-Chain Monte Carlo (MCMC) method. For each dC, in each filter with a detected period, we started 50 MCMC walkers in a Gaussian around the detected period. We sampled the walkers for 10,000 steps each, at each step using the phase dispersion minimization technique (Stellingwerf 1978) to calculate the likelihood at each walker position. We used the 1σ of the marginalized period distribution as the photometric period error for that dC. However this is only a statistical error, and does not account for the possibility that we have selected an alias rather than the true period.
Our final dC sample contains 34 individual dCs that are periodic in at least one ZTF filter. Given the wide initial orbits necessary for progenitor dC systems to avoid truncation of the TP-AGB phase before enough C-rich material can be transferred, it is remarkable that 19 (56%) of these dCs have periods < 1 d (and 28 (82%) of these dCs have periods < 2 d), indicating they should have experienced a common-envelope (CE) phase. The likely origins of the variability in these dCs include spot rotation on the dC or tidal distortion of the dC atmosphere from being in a close orbit with a WD. Since many of these dCs have short periods, we assume that these systems would have experienced a CE phase and have circularized and synchronized (Hurley et al. 2002). However, if the light curve variability is from the dC being tidally distorted, our detected period would be half the orbital period (even with 2× longer true orbital periods, these systems should still have experienced a CE phase).
In addition, the 1 d aliasing caused by the observing window function causes peaks at frequencies of ±1 d −1 . These alias peaks can also meet our log (FAP) limit (see Figure 1), and while we take the highest significance peak from the filter which produces the best fitting model (via a χ 2 fit) there is a possibility this is the wrong period. This can only be solved by either low cadence photometry or confirming the photometric period with RV follow-up. For example, the dC SBSS 1310+561 has ±1 d −1 aliases the meet our log (FAP) limit. In our initial search using ZTF DR3 data, the g and r filters had different highest peaks, with the best period in the g filter being ∼5.18 d and the best period in the r filter being ∼0.838 d. These are separated exactly by the 1 d −1 aliasing of the window function, with the r filter providing a better model fit. However, using the newest (and larger) ZTF DR5 data set results in both the g and r filters having the same highest peak, at 5.1878 ± 0.0012 d. Whitehouse et al. (2021) confirmed this as the period with their RV observations of this dC.
The recent work by Whitehouse et al. (2021) included modeling of light curves for their sample of periodic dCs. They examined whether spot rotation, tidal distortion or irradiation by a hot WD companion could be the source of the photometric variability in their dC sample. They found that for periods near and longer than 1 d, both tidal distortion and irradiation are reduced to levels that would not be detectable in the light curves. While tidal distortion could be detectable for our shortest period dCs in this paper (0.1 < P < 0.2 d), the predicted amplitudes at these periods from Whitehouse et al. (2021) are larger than those we find in our sample (as was the case for their sample). The irradiation modeling of Whitehouse et al. (2021) (which assumes a WD temperature from 30,000 K to 20,000 K) predicts amplitudes large enough to be detected. However, the majority of our dC amplitudes are smaller than predicted by the models, suggesting irradiation is not the source of variability for most of our dCs. Additionally, the majority of the dCs in our sample are mid-and latetype dCs (see Roulston et al. 2019) which do not have a visible WD in their optical spectra. This sets a limit that for these types, that WD must be cooler than about 10,000 K, reducing the irradiation effects below our detection limits. Only the composite dC+WD systems which have a hot WD (like SDSS J151905.96+500702.9) may have detectable irradiation effects. This leaves spot rotation as the most likely source of variability in our periodic dC sample. However, the origin of the variability in these dCs is not truly confirmed without comparison to spectroscopic RV follow-up.
Finally, as the ZTF survey continues to accure more data we expect to find more photometrically variable dCs from the current sample of known dCs. However, even given a favorable inclination (say i > 85 • ) the ZTF errors are too large for the detection of an eclipse of a cool WD in these systems. Using our estimated dC radii and luminosities and assuming a WD companion with a standard mass of 0.6M and temperature of 7000K (we expect to see the WD component in the optical spectra if it is any hotter than this) we would only expect an average primary eclipse depth of 0.005 mag. This is below our detection threshold with ZTF, with our dCs having median errors of 0.019 mag, compared to their median amplitude of 0.059 mag. The Vera Rubin Observatory's LSST survey (Ivezić et al. 2019) is expected to have errors of approximately 0.005 mag for a point source with r = 19.0, which may allow for the detection of dC eclipses. However, the majority of known dCs reside outside the LSST footprint, with 17% below the declination cut of of δ = +2 • (Ivezić et al. 2019). Detection and characterization of the first eclipsing dC system will likely require dedicated observations with high cadence and low photometric noise.   Rows marked with an * have a light curve with no detectable variability at the given period in this filter. The model fit for this filter is unreliable. dCs marked with †, while meeting all our criteria to be included, have suspect periods due to having more than 1% of the periodogram above the power needed to have log (FAP) ≤ −5.
Note-Periodic properties of dwarf carbon stars found in this paper. For each dC, we list light curve properties of the observed ZTF filter, the mean magnitude and mean magnitude error. We include the selected best period, the logarithm of the false-alarm-probability for that period, the amplitude of variability from the best fit model at that period, and the time of light curve maximum brightness. Finally, we include a few diagnostics including the number of terms in our model fit and the resulting reduced χ 2 of the model.

Spectroscopic Set-Up
The dCs J1519 and J1230 were observed with the Binospec spectrograph on the MMT telescope . For all observations, we used the 0.85 slit with the 600 l mm −1 grating centered on 7250Å, giving coverage from 6000Å to 8000Å covering Hα and the CN bands. The reduced spectra have a dispersion of 0.61Å pix −1 with R≈3590. All Binospec data were reduced using the standard Binospec reduction pipeline 1 (Kansky et al. 2019).
The dC J0625 was observed with the Magellan Echellette (MagE; Marshall et al. 2008) spectrograph on the Magellan Baade Telescope. All observations used the 0.85 slit and were reduced using the MagE reduction pipeline 2 (Chilingarian 2020). The reduced spectra cover from about 3200Å to 10000Å with R≈4500.
Observations for SBSS 1310+561 were acquired at the 1.5m Fred Lawrence Whipple Observatory (FLWO) telescope with the FAST spectrograph (Fabricant et al. 1998) using the 600 l mm −1 grating and the 1.5 slit, which provides wavelength coverage from 6000Å to 8000Å at 1.5Å spectral resolution.

SDSS J151905.96+500702.9
One of the more interesting dCs in our periodic sample, with photometric periodicity detected with highest significance, is SDSS J151905.96+500702.9 (also known as CBS 311; we use J1519 in the rest of this paper), a dC+DA spectroscopic composite binary. J1519 was discovered by Liebert et al. (1994) and has been studied on numerous occasions (Farihi et al. 2010;Green 2013;Whitehouse et al. 2018;Ashley et al. 2019;Roulston et al. 2019;Green et al. 2019). However, this is the first reporting of its periodic variability. J1519 (r = 17.3 mag) has four epochs of optical spectra in the SDSS, with the most recent spectrum shown in the top panel of Figure 2. The spectrum of J1519 shows a dC with a hot DA WD companion, as well as Hα emission. Whitehouse et al. (2018) and Roulston et al. (2019) found RV variability using few-epoch spectroscopy with ∆RV max of 46.8 ± 15.8 km s −1 and 44 ± 20 km s −1 , respectively. Farihi et al. (2010) conducted a study of WD-red dwarf systems, including J1519, using the Hubble Space Telescope. They found J1519 to be unresolved, placing the constraint on its separation of < 10 au.

J1519 WD Model Fits
Since J1519 is a spectroscopic dC+DA composite, we can fit WD model atmospheres to the WD component to fit T eff and log(g) using the SDSS spectra. Bédard et al. (2020) fit WD models and found fit values of 31230±210 K and 7.97±0.05 respectively. Farihi et al. (2010) found that spectroscopically fit WD parameters are often biased due to a cool companion. To update the fits of Bédard et al. (2020), we performed our own model atmosphere fits to the DA component of J1519 using the synthetic WD model atmospheres of Levenhagen et al. (2017). We first fit the late type dC (dCM) template of Roulston et al. (2020) to the SDSS spectrum of J1519 by finding the best-fit velocity, shifting the template, and then scaling it to the flux near Hα. We then removed the dC spectrum from the total spectrum, leaving just the WD component. We then fit the visible Balmer lines from Hβ and blue-ward to the entire grid of WD model spectra. We interpolated the grid of WD model spectra to include half-steps in the model space. Our best-fitting model parameters for T eff and log(g) were 31000±500 K and 7.85±0.05, respectively, and can be seen in Figure 2. The black line is the single SDSS spectrum with the highest S/N shifted to the rest-frame, and the blue line is the best fit WD model spectrum. We did not use Hα for the WD fit as the dC component contributes most to the spectrum in emission. In addition, we did not use the H9 line, as only half of the line is visible in the SDSS spectrum.
The WD temperature of our fit is in good agreement with that of Bédard et al. (2020). However, our log(g) is 0.12 dex lower, resulting in both our WD mass and cooling age being lower than those in Bédard et al. (2020). For the purposes of this paper, we adopt our fit values of log(g) and T eff . The WD properties we use can be found in Table 3, with the mass, radius, and cooling age coming from the models of Fontaine et al. (2001).

J1519 Radial Velocities
Although RV variability has been detected in J1519, there are no published RV orbital fits for this system. Based on our photometric analysis, we found a period of 0.302356 ± 0.000021 d (∼ 7 hr) for J1519. Therefore, we conducted a spectroscopic monitoring of J1519 using Note-Best fit model parameters for the DA component of SDSS J151905.96+500702.9. Each parameter lists the source used: (1) this paper (2) from evolutionary models of Fontaine et al. (2001) the MMT spectroscopic setup as was described in Section 5.1. On the nights of 2020 August 19 and 20, we observed a sequence of 21×200 s exposures, on the night of 2020 August 22 we observed 27×200 s exposures, and on the nights 2021 April 21 and 23 we observed 24×230 s exposures. The exposures on each night were then coadded in threes, resulting in seven final epochs on the first two nights, nine epochs on the third night, and eight on each of the last two nights for a combined total of 39 epochs (with about 600 s total exposure each), with an average S/N ≈ 5 for all epochs in the continuum region near Hα.
Since the full spectrum includes both stellar components, we measured the RV from the Hα emission line, presumed to come from the dC atmosphere alone. First, for each epoch, we re-scaled the late-type (dCM) dC template of Roulston et al. (2020) to the flux in our MMT spectrum in the region of 6300-6500Å. We then used this as the model for the dC continuum level of that epoch, which was used to calculate the Hα emission line center, equivalent width and associated errors. The RV measurements have an average error of approximately 5 km s −1 , and the equivalent width measurements have an average error of approximately 0.18Å. Figure 2 shows the measured RV (middle) and Hα equivalent widths (bottom) for J1519. To fit the RV curve, we used the rvfit program which uses a simulated adaptive annealing procedure, the details of which can be found in Iglesias-Marzoa et al. (2015). We left all parameters free to be fit, with the solution quickly converging to a circular orbit. We therefore refit the RV curve leaving all parameters free again except for the eccentricity, which we fix to e = 0.0. The resulting bestfit model can be seen in Figure 2 (red curve) and the fit parameters can be found in Table 4.
The best fitting orbital period from the RVs (0.327526 ± 0.000012 d) is longer than the best photo-metric period by 0.025170 ± 0.000024 d (about 36 minutes). Fixing the period in the RV fitting procedure to that of the photometric period results in a poorer model fit, with the longer period model being a better fit at the 3.2σ level. The best-fit semi-amplitude of K 2 = 33.3±1.4 km s −1 (∆RV max = 2K = 66.6 km s −1 ) is in agreement with the RV variations found by Whitehouse et al. (2018) and Roulston et al. (2019), as their random epochs likely did not catch the true RV amplitude. However, the low measured semi-amplitude suggests an extremely low inclination of this system, with i ≈ 10 • if we take our estimated dC mass of 0.41 M from Section 6.2.
One possible explanation for a longer orbital period than photometric period is that J1519 was spun up by the accretion that it experienced, and has not yet synchronized the rotation and orbital periods in the approximate 8 Myr since mass transfer stopped (assuming the mass transfer ceased at the same time the WD formed). Green et al. (2019) analyzed Chandra observations of J1519 (as well as five other Hα emission dCs) and found it to show X-ray emission consistent with having a short rotation period, which would lend support to the accretion spin-up scenario. Deeper photometric imaging and RV follow-up, particularly of the WD component, could even better characterize this system. It is clear, however, that this dC has both photometric and RV variability on a <0.33-d timescale, indicating it most likely has a short orbital period and formed through a CE event.

SDSS J123045.53+410943.8
Another interesting dC is SDSS J123045.53+410943.8 (J1230), whose SDSS spectrum is shown in Figure 3. J1230 shows the C 2 and CN lines typical of late type dC stars, but also shows strong absorption lines of K and a strong CaH band near 6800Å. Additionally, this dC shows strong emission lines of Hα, Hβ, Hγ, Hδ, and Ca H and K. Unlike J1519, there is no visible WD component in the spectrum of J1230.
We observed J1230 (r = 16.6 mag) with the Binospec spectrograph on the 6.5m MMT using the setup described in Section 5.1. We took 6×266 s spectra on the nights of 2021 February 3, 5, 7, 8, 9, 11 and 15. This resulted in a total of 42 spectra, with an average S/N of 13 in the continuum region near Hα. We measured the line center and equivalent widths of the Hα line, as well as for the two K lines visible in our MMT spectra of J1230. The emission and absorption lines have the same velocities, indicating the are coming from the same region. We average these velocities to measure the RV for each of the 42 epochs.  In the same method as J1519, we used the rvfit program to fit the RV curve of J1230. For this dC, we left all parameters free for fitting, with the resulting best period fit matching that of the photometric light curve (P = 0.882519 ± 0.000020 d). We therefore fix the period to the photometric period and the eccentricity to 0.0, and refit the RVs. The resulting best fit can be found in Table 4 and the phased RV curve in Figure 3. As with J1519, we find the parameter errors using a MCMC method centered around the best fit parameters.
The best fit gives a circular orbit with a semiamplitude of K 2 = 123.0 ± 0.7 km s −1 . If we use our estimated dC mass from Section 6.2 (0.25 M ) and an assumed WD mass of 0.6 M the implied inclination of this system is around i = 56 • .
The presence of multiple emission lines of H and Ca suggest that the photometric variability of J1230 is coming from re-processing of the WD flux on the surface of the dC. Even though the WD companion to J1230 is not hot enough to be seen in the optical spectrum, it may be warm enough to still heat the surface of the dC. If this is true, we could expect the dC to be at maximum brightness when the WD-facing side is pointed toward us maximally, i.e. when the dC is moving transversely on the sky, between the ascending and descending nodes. Comparing the light curve of J1230 to the RVs however shows this is not the case, as if it were, we would expect the RV to be moving to the descending node after the photometric maximum, which the RV curve for J1230 is 0.33 out of phase with. This suggests that the photometric variability is not coming from re-processing, but rather from spot rotation on an active dC (with the emission lines indicating chromospheric activity). We do note that the uncertainties on our epoch of maximum brightness and period may cause us to incorrectly predict the phase when our spectroscopy was collected in 2021 February by up to 0.03 cycles.

SBSS 1310+561
We observed SBSS 1310+561 (r = 14.1 mag) using the FAST spectrograph on the 1.5m telescope at FLWO using the setup described in Section 5.1. We took 3×300 s spectra during the nights of 2021 February 10 and 11,  The blue curve is the best fit single sinusoid model to the data, while the grey dotted line is the average equivalent width value. The y-axis has been inverted so that smaller equivalent width values (more emission) are up. We do not detect any significant variation in phase for the equivalent width, to a limit of <1.50Å. and 6×300 s spectra on the night of 2021 February 12 for a total of 12 spectra with an average S/N of 32 in the continuum region near Hα.
Unlike with our MMT and Magellan observations, because of observing time constraints on our awarded FAST time, we chose to obtain these spectra close to the quadrature phases based on the ZTF photometry (P = 5.1878 ± 0.0012 d). We assumed that the photometric period corresponds to the orbital period, and used t 0 from the light curve to calculate the expected times that SBSS 1310+561 should be at the quadrature phases (φ = 0.25 and φ = 0.75). Our actual observations were taken at phases φ = 0.27±0.02 and φ = 0.47±0.01.
From these spectra, we measured the RV at φ = 0.27 to be V r = −79.7 ± 9.5 km s −1 and the RV at φ = 0.47 to be V r = −19.3 ± 5.5 km s −1 . Taking the difference in these two velocities for this system (∆RV = 60 ± 11 km s −1 ) can place a lower limit on the semi-amplitude of K 2 > 30 ± 6 km s −1 . Using our estimated mass from Section 6.2 (0.46M ) and an assumed WD mass of 0.6 M , this constrains the inclination to i ≥ 25 • (if i = 60 • , then we would expect K 2 = 109 km s −1 ). Since the phase difference between our two epochs is quite small, this ∆RV suggests that SBSS 1310+561 is in a tight orbit and is very likely a PCEB.

LAMOST J062558.33+023019.4
The dC LAMOST J062558.33+023019.4 (hereafter J0625, r = 13.9 mag) was observed on the nights of 2021 January 11 and 12 using the Magellan MagE instrument setup described in Section 5.1. Each night, we observed 15×300 s exposures. The final reduced spectra consists of 30 epochs with an average S/N of 22 each in the continuum region near Hα.
Using the Hα emission line, we measured the RV of J0625 for each epoch. We found no evidence for RV variability, nor any variability in Hα equivalent width. We found the RV to vary with only a standard deviation of 3.9 km s −1 , and with a maximum ∆RV = 12.1 ± 3.2 km s −1 . In addition, cross-correlation of the spectra across epochs resulted in no significant measured RV variations. Using our estimated mass from Section 6.2 (0.84 M ) and an assumed WD mass of 0.6 M , this places a constraint on the inclination of i ≤ 7 • (if i = 60 • , then we would expect K 2 = 44 km s −1 ). This may suggest that the photometric variability is not related to the orbital period in this system since such a low Note-Fit parameters from the radial velocity follow-up. The value for each parameter is given as the median of the marginalized distribution of the MCMC samples. The errors for each parameter are the 1σ values from the marginalized distribution of the MCMC samples. Additionally, derived values for the orbital separation (a sin i) and mass function (f (m1, m2)) are given.
inclination (and low semi-amplitude) is unlikely if the photometric period of 7.6080 ± 0.0014 d represents the orbital period. Hence, this system adds weight to the evidence that the photometric variability in dCs may often be due to spot rotation.

COMMON ENVELOPE CONNECTION
For the progenitor of the dC companion to become a C giant, it must enter the third-dredge up phase (Iben 1974). AGB stars have a degenerate CO core with a double-shell (moving outward from the core) of helium and hydrogen. As the hydrogen shell (which produces most of the energy) continues to fuse H into He, the helium shell surrounding the core continues to grow. Eventually, the helium shell experiences runaway fusion, driving expansion of the envelope material above. This He-shell "flash" and expansion means the star is now in the thermally pulsing AGB (TP-AGB) phase. Helium shell fusion causes the inter-shell region to become strongly convective, dredging helium fusion products to the surface, i.e., the third dredge-up. As the expansion continues, the pressure in the helium shell will drop, eventually stopping its energy production. The layers contract again with hydrogen shell fusion resuming, and the cycle repeats.
Each successive thermal pulse becomes stronger, reaching deeper into the intershell zone, and the stellar radius increases (Iben & Renzini 1983). As helium shell fusion products are brought to the surface, it is possible that the envelope carbon abundance increases until C/O > 1. Since C preferentially binds with O, C 2 and CN bands only appear when C/O> 1, forming a C giant star.
AGB stars going through the TP-AGB phase can reach radii of 800 R (3.7 au) as they experience successively stronger thermal pulses (Marigo et al. 2017). Assuming an AGB mass of 2.5 M , AGB radius of 800 R , and a dC mass of 0.4 M , this system would experience the beginning of a common-envelope (CE) with an initial period of ≈ 4.2 yr (if the dC mass is 1.0 M instead, then P≈ 3.8 yr). Therefore, dCs with initial periods ≈ 4 yr (1500 d) or less will very likely have experienced a CE phase, corresponding to the shorter-period peak modeled by de Kool & Green (1995). The dCs in this paper with P< 1 d are most certainly the result of a CE spiral-in. Of the six dC periods in the current literature, two of them have P< 3 d, so have likely experienced a CE. It seems then that many dCs may have experienced a CE phase.
Dell'Agli et al. (2021) recently studied the extreme AGB stars (those AGB stars which have extremely red mid-IR colors, e.g. Gruendl et al. 2008) and showed that the excess dust and outflow densities of these stars may be explained by envelope stripping in a CE event. Their models suggest that these extreme AGB stars are actually post-common-envelope binaries (PCEBs) with orbital periods of order 1 d, matching the periods for dCs in our sample. Dell' Agli et al. (2021) also found that the CE in their models starts after the rapid growth of the AGB radius, once the C/O ratio increases past unity, which corresponds well with the requirements for producing the short-period dCs we find. This makes these extreme AGB stars potential progenitors systems of the dCs that are in the CE phase currently.
However, is mass accretion during a CE phase the most likely mass transfer mechanism to form dCs? We can address this question by looking at our periodic dC sample in the context of models that simulate expected binary populations.

Binary Population Synthesis Models
We used the binary population synthesis (BPS) models of Toonen & Nelemans (2013) to see if the observed population of dCs can be reproduced by theory. The full details of the BPS models can be found in Toonen & Nelemans (2013) and are briefly described here.
These BPS models were created using the SeBa (Portegies Zwart & Verbunt 1996;Nelemans et al. 2001;Toonen et al. 2012;Toonen & Nelemans 2013) population synthesis code. This code generates an initial population of binaries and simulates their evolution, taking into account processes such as stellar winds, magnetic breaking, mass transfer, common-envelope, and angular momentum loss. The initial stellar population is generated from the classical BPS distributions found in Toonen & Nelemans (2013) via a Monte Carlo method. The resulting binaries are then convolved with a Galactic model including a star formation history that depends on time and location in the Milky Way based on Boissier & Prantzos (1999) so that the simulated binaries can be compared to our observed sample.
For the synthetic populations used here, the commonenvelope phase is modeled on the basis of the energy budget i.e. the classical α-formalism of Tutukov & Yungelson (1979). We discuss the results of two different models here that account for two different CE efficiencies: model αα and αα2 which have αλ of 2 and 0.25, respectively. The parameter λ is the structure parameter of the envelope to calculate the envelope binding energy (Paczynski 1976;Webbink 1984;de Kool et al. 1987;Livio & Soker 1988;de Kool 1990;Xu & Li 2010). The α parameter describes the efficiency with which orbital energy is consumed to unbind the CE. A smaller value of α implies less efficient usage of orbital energy, and therefore a stronger shrinkage of the orbital period during the CE-phase. We do not consider the orbital angular momentum method of Nelemans et al. (2000), as this model does not reproduce the observed characteristics of the general PCEB (WD/main sequence) population (Toonen & Nelemans 2013).
Furthermore, the BPS models here allow for accretion during the common-envelope phase. The accretion rate is limited by the thermal timescale of the accretor times a factor that is dependent on the stellar radius and the corresponding Roche lobe (Portegies Zwart & Verbunt 1996;Toonen et al. 2012) following Kippenhahn & Meyer-Hofmeister (1977); Neo et al. (1977); Packet & De Greve (1979); Pols & Marinus (1994). The total accreted mass is then given by the integral of the accretion rate times the timescale of the CE event, which here is taken to be 100 yr. This timescale is consistent with hydrodynamical simulations (Ricker & Taam 2008;Ivanova et al. 2013) and observations of hot subdwarf binaries (Igoshev et al. 2020), although cataclysmic variables may suggest a longer CE timescale, up to 10 4 yr (Michaely & Perets 2019;Igoshev et al. 2020).

BPS Comparison to Observed dC Sample
We use the resulting model population for a direct comparison to our observed sample of short-period dCs, assuming the photometric period is the current orbital period. To do this, we estimated dC masses based on their infrared absolute magnitude M K in the K band.
Comparisons of M dwarf spectra (Ivanov et al. 2004) to C star spectra (Tanaka et al. 2007) reveal them to be much more similar in the infrared than in the optical region. We used K s band (2.159 µm) magnitudes from the Two Micron All-Sky Survey (2MASS; Skrutskie et al. 2006). Six of our periodic dCs do not have K s band magnitudes. For these, we first fit Gaia absolute G band (M G ) to the dCs that do have K s band magnitudes. This fit was then used to convert the Gaia M G into M Ks for the dCs lacking K s band magnitudes. We then fit M Ks for our dCs to stellar masses using data from Kraus & Hillenbrand (2007). This fitting also provides us with bolometric luminosities for the dCs in our sample. Comparing our bolometric luminosities to those provided in Green et al. (2019) (who used a spectral energy distribution method fitting 0.35 − 12.5 µm) for the four dCs that overlap, we find our luminosities agree within 3%, indicating our dC mass estimates should be reliable.
Our mass estimates can be found in Table 5. We find that none of our dCs are fit with masses > 1 M or < 0.2 M , in agreement with the range for which detectable C 2 , CN, and CH bands are expected. We note that some of the lowest mass dCs may have been brown dwarfs or even planets before they accreted significant C-enriched material from their former AGB companion.
Using the mass-radius relationship for main-sequence stars of Eker et al. (2018), we estimate the radius for these periodic dCs as well, which are included in Table 5. Using these estimated radii we calculate the Roche-lobe filling factor (RLFF), using the equation of Eggleton (1983) to find the Roche radii. Six out of 34 of our periodic dCs may be experiencing RLOF back onto the WD (all have a RLFF > 1 in Table 5). However, we caution that physical parameters are derived from Orich main-sequence models, which may not accurately represent all dCs. For example, (1) we do not know the mass of the unseen WD companion and assume it is 0.6 M (2) we assume these mass-radius and M Kmass relations hold for dCs, as they do for normal Orich stars (3) dCs are thought to be of a lower metallicity population and studies have found that low metallicity M dwarfs may have smaller radii (Kesseli et al. 2019) and (4) since dCs may have increased activity and magnetic fields due to their mass accretion, their radii may be inflated (Kesseli et al. 2018). We see no obvious evidence of flickering or accretion outbursts in any of our ZTF light curves that might indicate current RLOF back onto the WD. Note-Distances and magnitudes for the periodic dCs in our sample. We use the Gaia distances, colors, and magnitudes, as well as the 2MASS absolute K magnitudes to estimate masses and bolometric luminosities for our dCs. For the solar bolometric luminosity, we adopt the value log 10 L = 33.58. We also calculate the Roche-lobe filling factor (RLFF) under the assumption of a 0.6 M WD companion. We calculate the mass errors to be of order 0.05 M , the log 10 (L bol /L ) errors to be of order 0.1, and the radius errors to be of order 0.05 R . However, we caution that physical parameters are derived from O-rich main-sequence models, which may not accurately represent all dCs.
To compare the BPS models directly to our observed dC sample, we applied a series of cuts and selection effects to the models, as follows: (1) P < 100 d (2) r mag < 19.5 (3) M dC ≤ 1 M (4) 1.0 <M ZAMS < 4 M (5) the initial primary must be a TP-AGB star at the onset of the CE phase. Here, M dC is the current mass of the main sequence companion in the BPS models, and M ZAMS is the initial mass of the primary at the beginning of the models (which will become the AGB donor).
We show the resulting BPS models in Figure 4 and Figure 5 -models αα2 and αα, respectively. In both figures, the BPS models are shown as the colored 2D histogram in mass and period (note that the histogram color scale is logarithmic and its range is different for each plot), and the periodic dC stars from this paper represented as the red scatter points (with KDE contours). The dashed black line represents the RLOF boundary, with systems occupying the region to the left (shaded in grey) filling their Roche lobes, under the as-  Figure 4. Binary population syntheses results for model αα2, which is the model with lower common-envelope efficiency (αλ = 0.25) in which less orbital energy goes to unbinding the CE. The background heat map shows the number of systems (colored on a log scale) in model αα2 going through a CE phase while the primary is in the TP-AGB phase. The plotted masses are for the secondary (i.e. the dC) and periods are for the currently observed system. Each panel uses a different magnetic braking formalism: (a) magnetic braking from Rappaport et al. (1983) (b) magnetic braking from Ivanova & Taam (2003) (c) magnetic braking from Knigge et al. (2011). The red plus scatter markers are our observed short period dC sample without Hα emission, while the purple triangles are dCs which show Hα emission. The red contours are a KDE contour for the entire periodic dC sample. The dashed black curve represents the boundary for current dC systems to be experiencing RLOF, where systems to the left in the grey shaded region may be experiencing RLOF. The bottom panel shows a histogram of the estimated mass accreted by the secondary during the CE phase. While this model (αα2) reproduces the period distribution better than model αα (see Figure 5), the estimated accreted mass is even lower than that of αα because of the lower CE efficiency. Since >0.03 M of mass transfer is likely required to create a dC, a CE is not likely to be the primary mechanism for accretion to form a dC.
sumption of a 0.6 M WD companion. Both figures also show a histogram of the estimated mass accreted during the CE phase (assumed to last 100 yr). Figure 4 shows model αα2 (αλ = 0.25) and includes three different magnetic braking prescriptions. Panel (a) uses the magnetic braking of Rappaport et al. (1983), panel (b) that of Ivanova & Taam (2003) and panel (c) that of Knigge et al. (2011). Again, the color scale is logarithmic and its range is different for each sub-figure.
Model αα2, however, does not reproduce the mass distribution of our dCs very well, generating low mass systems than observed (still under the assumption that our physical parameters derived from O-rich main-sequence models apply to dCs). While it may be that model αα2 does not produce low mass dCs, we have not considered our sample selection effects in this comparison. Our observed sample is likely biased toward lower mass dCs as (1) they have stronger C 2 and CN bands, and (2) their variability is fractionally larger and so easier to detect.
Model αα ( Figure 5) uses a higher CE efficiency (αλ = 2) similar to classical BPS studies and includes the stan-dard magnetic braking of Rappaport et al. (1983). From Figure 5, it is seen that this model is unable to reproduce the short periods of our observed dC sample. This is in agreement with the conclusions based on the SDSS PCEBs (WD+MS systems; Zorotovic et al. 2010;Toonen & Nelemans 2013;Camacho et al. 2014), where Toonen & Nelemans (2013) found that standard efficiency (αλ = 2) CE was also unable to reproduce the observed periods, as it generated too many long-period PCEBs.
A crucial shortfall is that the estimated mass accreted for all models is too small to convert a main-sequence star into a dC (see Section 7 for a discussion). Miszalski et al. (2013) suggest that to shift the secondary envelope from approximately solar (C/O) i ∼ 1/3 to the observed (C/O) f > 1 requires ∆M 2 = 0.03-0.35 M for M 2 = 1.0-0.4 M . The predicted mass accretion in our BPS models is lower than this by 2−3 orders of magnitude. Together with the strong mismatch between the modeled and observed dC period-mass distributions, it seems clear that there must be further mass accretion outside the brief CE phase. The background heat map shows the number of systems (colored on a log scale) in model αα going through a commonenvelope phase while the primary is on the TP-AGB phase. The plotted masses are for the secondary (i.e. the dC) and periods are for the currently observed system, including standard magnetic braking from Rappaport et al. (1983 Qualitatively it is also possible to argue that accretion during CE evolution is rarely significant for nondegenerate companions. The common envelope itself typically possesses much higher specific entropy than the surface of the accretor, with the consequence that matter accreted by the companion star reaches pressure equilibrium at the surface of that star with much higher temperature, and vastly lower density, than the accre-tors initial surface layer. A temperature inversion or roughly isothermal layer is expected to bridge this entropy jump with the result that, over the duration of the CE evolution (which is much shorter than the thermal time scale of the accretor), the accretor is thermally isolated from the common envelope, while the common envelope itself becomes increasingly tenuous.
A solution to the under-prediction of accreted mass may be that more mass accretion may take place before the CE phase, but during the TP-AGB phase, by wind accretion during wind-RLOF (WRLOF; Mohamed & Podsiadlowski 2007). WRLOF is a mass transfer state that lies between standard wind mass transfer and standard RLOF, where the wind of the primary star is focused toward the secondary star. This results in increased mass transfer efficiency as compared to the standard Bondi-Hoyle-Lyttleton case (Hoyle & Lyttleton 1939;Bondi & Hoyle 1944).
In the WRLOF regime, the primary is technically not filling its Roche lobe. However, low-velocity wind matter is funneled through the Roche lobe to the companion, allowing for mass transfer to take place in binaries with wider orbits than the classical RLOF case. WRLOF would boost dC formation since, if the initial orbital separation is too small, the expanding AGB atmosphere can cause a CE before the third dredge-up can turn the AGB into a C star.
A variety of simulations (Abate et al. 2013;Saladino et al. 2018) have shown that WRLOF in binaries with AGB primaries can have mass-transfer efficiencies of 40-50%. For an average AGB wind mass loss rate of 10 −7 -10 −4 M yr −1 (Höfner 2015), a main sequence companion could accrete enough material (∼ 0.35M ) in only 10 3 -10 6 yr, within the time an AGB star is expected to stay a C star (10 6 yr; Marigo et al. 2017). WRLOF has also been shown to efficiently tighten the orbit (Saladino et al. 2018;Chen et al. 2018) so that more systems could be driven towards orbits with the periods we find.
Indeed, it appears the WRLOF may be the dominant mass transfer mechanism for many chemically peculiar stars. Abate et al. (2013) showed that simulations of AGB binaries with WRLOF were better able to reproduce the observed formation rates of the CEMP-s stars.  and  also performed simulations finding AGB binaries with WRLOF were consistent with the observed properties of the CEMP-s, CH, and Ba stars. This further strengthens the connection between dCs, WRLOF, and the other chemically peculiar stars.
However, detailed modeling of the WRLOF in the specific case of the TP-AGB phase with a C star is needed to understand how the larger stellar radii of AGB-C stars and the increased dust formation often found in their winds affect the WRLOF mass transfer efficiencies. The periodic dCs in this paper represent a prime sample that is ready for spectroscopic follow-up and for comparison to future models of WRLOF mass transfer.

SUMMARY
We searched ZTF light curves of a sample of 944 dCs for periodic signals. We found 34 dCs with signs of significant photometric variability, with 82% having P< 2 d. The most likely origins of this periodicity is either from spot rotation or surface heating of the dC from the close WD companion. Even if the detected ZTF periods arise from ellipsoidal variations and represent half the orbital period, such short periods are surprising for dC stars, which require significant accretion from a TP-AGB C giant to turn them into the C-enriched dwarfs we see today.
Spectroscopic follow-up is needed to determine the source of the variability in each of these periodic dCs, especially to confirm for the case of spot rotation that the system is circularized and tidally locked so that the rotation period can be assumed to equal the orbital period (i.e. that our reported photometric periods correspond to both rotational and orbital periods). The periodic dCs in this paper provide a rich new sample to target for spectroscopic follow-up, as well as to study dC formation and properties. We have confirmed the photometric period as the orbital period for one dC for which we have obtained spectroscopic follow-up. In two other dCs we have confirmed that they must have short orbital periods from their RVs (not confirming the photometric period as the orbital period however). In all three cases, these short (P< 1 d) orbits indicate these dCs have indeed experienced a CE phase.
These short periods indicate that at least some dCs will experience a CE at some point in their formation, with P< 1 d dCs having experienced experiencing substantial plunge-in. We used binary population synthesis models to show that the observed sample of dCs is not well-reconstructed by mass transfer during the commonenvelope phase alone, since the dCs in our sample require at least 0.03 M of mass accretion but our models predict 2-3 orders of magnitude less transfer during the CE phase, suggesting mass accretion before the CE phase. However, some systems such as cataclysmic variables indicate CE timescales an order or two longer (Michaely & Perets 2019;Igoshev et al. 2020) than our assumed 100 yr (based on Ricker & Taam 2008;Igoshev et al. 2020), which may substantially increase the amount of accreted material to the point that the CE alone could provide enough mass to form a dC.
Hydrodynamical simulations of CE evolution typically find that accretion onto a non-degenerate companion is not common, because of the entropy barrier between the companion and the surrounding material (e.g. Ivanova et al. 2013, for a review). In fact, even in the case of neutron star companions, which can accrete more efficiently due to neutrino cooling, accretion is limited to 0.1M (MacLeod & Ramirez-Ruiz 2015). Further modeling of the CE phase involving C-AGBs may provide further insight.
dC systems that begin as very wide binaries would experience stable mass transfer and a widening of the orbit. Systems that initially are close would begin a CE phase either during the red giant branch or during the AGB before the TP-AGB and, without experiencing the third dredge-up during the TP-AGB, would not produce dCs. Therefore, it seems that the most likely mass transfer mechanism to form dCs is WRLOF.
Further modeling of WRLOF in binaries with a TP-AGB star are needed to fully test this formation pathway of dCs. Additionally, further work is needed to understand the relationship between initial dC metallicity and mass to constrain the amount (and composition) of material that needs to be accreted to form a dC. This would be an important step in constraining the mass-transfer efficiency in the WRLOF case.  (Hunter 2007), Numpy (Harris et al. 2020), Scipy (Virtanen et al. 2020), Scikit-learn (Pedregosa et al. 2011), TOPCAT (Taylor 2005)