The ELM Survey. IX. A Complete Sample of Low Mass White Dwarf Binaries in the SDSS Footprint

We present the discovery of 17 double white dwarf (WD) binaries from our on-going search for extremely low mass (ELM)<0.3 Msun WDs, objects that form from binary evolution. Gaia parallax provides a new means of target selection that we use to evaluate our original ELM Survey selection criteria. Cross-matching the Gaia and Sloan Digital Sky Survey (SDSS) catalogs, we identify an additional 36 ELM WD candidates with 17<g<19 mag and within the 3-sigma uncertainties of our original color selection. The resulting discoveries imply the ELM Survey sample was 90% complete in the color range -0.4<(g-r)_0<-0.1 mag (approximately 9,000 K<Teff<22,000 K). Our observations complete the sample in the SDSS footprint. Two newly discovered binaries, J123950.370-204142.28 and J232208.733+210352.81, have orbital periods of 22.5 min and 32 min, respectively, and are future LISA gravitational wave sources.


INTRODUCTION
This paper presents new low mass WD observations designed to assess the completeness of the ELM Survey. The ELM Survey is a spectroscopic survey targeting extremely low mass <0.3 M ⊙ , He-core WDs. To form such low mass WDs within the current age of the Universe requires binary evolution (Iben 1990;Marsh et al. 1995). Different evolutionary pathways can form low mass WDs, and both common-envelope and Roche lobe overflow evolution channels contribute to the observed low mass WD population (Istrate et al. 2016;Sun & Arras 2018;Li et al. 2019). The ELM Survey is responsible for the majority known low mass WDs and the discovery of over 100 detached double degenerate binaries (Brown et al. 2020a,b;Kilic et al. 2021).
Importantly, WD binaries with periods P 1 hr are potential multi-messenger sources detectable with both light and gravitational waves with the future Laser Interferometer Space Antenna (LISA, Amaro-Seoane et al. 2017). Theoretical binary population models predict that the majority of multi-messenger WD binaries will contain low mass, He-core WDs (Korol et al. 2017;Lamberts et al. 2019;Breivik et al. 2020;Thiele et al. 2021). Observational surveys provide an opportunity to anchor these population models with actual numbers. The difficulty is that observations are based on optical magnitude and color, not model parameters like WD mass or binary period. A complete and well-defined ob-servational sample provides the best basis for comparison.
Gaia astrometry now provides a powerful means to target WDs. Gentile Fusillo et al. (2021) use Gaia early Data Release 3 (eDR3) to identify hundreds of thousands of high-confidence WDs over the entire sky. Pelisoli & Vos (2019) used Gaia DR2 to identify candidate low mass WDs. The original ELM Survey target selection, for comparison, was made on the basis of magnitude and color using SDSS passbands. Here, we use Gaia astrometry to identify ELM WD candidates in the SDSS footprint within 3-σ of our original color selection criteria, and present the observational results.
We discuss the details of the target selection in Section 2. We use de-reddened magnitudes and colors indicated with a subscript 0 throughout. In Section 3, we present the observations and derive physical parameters. These measurements complete the color-selected ELM Survey sample in the SDSS footprint. In Section 4, we discuss the ELM Survey completeness and compare it with other Gaia-based WD catalogs. We close by highlighting J1239−2041 and J2322+2103, which are among the shortest period detached binaries known, and discuss the gravitational wave properties of the full sample of 124 WD binaries published in the ELM Survey.

ELM WD TARGET SELECTION
The original ELM Survey used de-reddened (u − g) 0 color as a proxy for surface gravity to select candidate low mass WDs. Normal log g ∼ 8 DA-type WDs differ from normal log g ∼ 4 A-type stars by up to 1 magnitude in (u − g) 0 color in the temperature range −0.4 < (g − r) 0 < −0.1 (about 9,000 K to 22,000 K) due to the effect of surface gravity on the Balmer decrement in the spectrum of the stars (Brown et al. 2012b). This region of color space thus provides an efficient hunting ground for log g ∼ 6 ELM WDs.
Color selection does not identify all possible ELM WDs, however, because of inevitable confusion with background halo A-type stars and foreground DA WDs with similar colors. During the course of the ELM Survey, we learned that the color selection becomes inefficient at <9000 K temperatures, or (g − r) > −0.1 mag, where the hydrogen Balmer lines lose their sensitivity to temperature and gravity. Low mass WDs at <9000 K temperatures exist, but subdwarf A-type stars are 100 times more common at these temperatures (Kepler et al. 2015(Kepler et al. , 2016. Subdwarf A-type stars are mostly field blue stragglers and metal poor halo stars (Brown et al. 2017;Pelisoli et al. 2017Pelisoli et al. , 2018bYu et al. 2019).
Fortunately, ELM WDs also have distinctive ∼0.1 R ⊙ radii that fall in-between normal ∼0.01 R ⊙ DA WDs and ∼1 R ⊙ A-type stars. These radii mean ELM WDs are found in a distinct range of absolute magnitude at a fixed color (temperature). Gaia now provides a direct means of selecting ELM WD candidates on the basis of parallax.
To test the ELM Survey completeness for low mass WDs, we select ELM WDs using Gaia astrometry in the SDSS footprint as described below. Our spectroscopic observations began with Gaia DR2, and then transitioned to using Gaia eDR3 when it became available. We consider the results of both Gaia data releases.

Matching SDSS to Gaia DR2
Our initial approach was to cross-match a wide range of stars with possible low mass WD colors in SDSS against the Gaia DR2 catalog, so that we properly sample the SDSS footprint. SDSS provides ugriz photometry for about 10,000 deg 2 of the northern hemisphere sky, mostly located at high Galactic latitudes |b| 30 • .
We use SDSS DR16 photometry (Ahumada et al. 2020), which is identical to SDSS DR13 values used by the Gaia Collaboration. We correct all magnitudes for reddening using values provided by SDSS. Since our WDs have a median distance of 1 kpc, and background halo stars are even more distant, applying the full reddening correction is important for accurate color selection.
We then performed an exploratory 10 arcsec radius positional cross-match with Gaia DR2. The exploratory cross-match yields 577,253 matches, 16,680 of which are duplicates and 10,469 of which have no parallax or proper motion. Accounting for proper motion, we find that a 10.7 yr difference in epoch minimizes the position residuals for the blue stars, consistent with the mean epochs of SDSS and Gaia DR2.
Finally, we perform a clean cross-match between SDSS and Gaia DR2.
We require ruwe<1.4, astrometric n good obs al>66, and duplicated source=0 following the guidance of Lindegren et al. (2018). Our goal is to end up with a clean set of low mass WD candidates, even if the ruwe cut may select against wide binaries like HE 0430−2457 (Vos et al. 2018). We accept objects that match within 1 arcsec in position, corrected for 10.7 years of proper motion, and that match within 1 magnitude between Gaia G band and the average SDSS g + r band magnitude. We reject objects with no parallax or proper motion. This cross-match yields astrometry for 543,534 blue SDSS stars, 94% of the original sample.

Matching Gaia eDR3 to SDSS
Following the release of Gaia eDR3, which has improved completeness for faint blue stars and 2× better astrometric errors than Gaia DR2, we transitioned to using Gaia's best neighbour cross-match with SDSS. This approach bases the selection on the Gaia catalog, not on the SDSS catalog, and so imposes an unseen astrometric figure of merit to the target selection. The advantage is that the Gaia best neighbour cross-match makes use of the full 5 parameter astrometric covariance matrix.
Starting with the same list of 578,833 blue SDSS stars, we use the gaiaedr3.sdssdr13 best neighbour catalog to find Gaia eDR3 matches for 99.9% of the stars. We reject 3,063 objects with no parallax, with matched positions that differ by >1 arcsec, or with number of mates>0. We then create a clean catalog by applying the Gentile Fusillo et al. (2021) data quality cuts for high Galactic latitude white dwarfs 1 , excluding the parallax over error and pm over error criteria. We will do our own parallax and proper motion selection in the next sections. The final Gaia eDR3 cross-match catalog contains clean astrometry for 567,769 blue SDSS stars, 98.4% of the original sample.
With the Gaia cross-match in hand, we are ready to make various cuts of ELM WD candidates beyond the original ELM Survey. Gaia eDR3 astrometry supersedes Gaia DR2 in accuracy and precision, so we use Gaia eDR3 values throughout this paper.

ELM WD Candidates: Photometric Selection
We consider the apparent magnitude range 17 < g < 19 mag to maximize our volume for finding ELM WDs. We find that Gaia parallax is reliable to around 19 mag for our targets (see below). We adopt −0.4 < (g − r) 0 < −0.1 mag to fix the temperature range to approximately 9,000 K to 22,000 K and avoid subdwarf A-type stars.
Our first ELM WD selection is based on (u−g) 0 color. For reference, Figure 1 plots (u − g) 0 versus (g − r) 0 for every star in SDSS with 17 < g < 19 mag and E(V − B) < 0.1 mag. The lower band of dots are normal A-type stars: halo blue stragglers and blue horizontal branch stars. The upper band of dots are normal DA-type WDs, stars with the same temperatures but much higher gravities. ELM WDs fall in-between. Magenta lines plot Althaus et al. (2013) evolutionary tracks for 0.16, 0.19, and 0.32 M ⊙ WDs ( Figure 1). Blue lines mark the approximate color-selection of the HVS Survey (Brown et al. 2012a) and the ELM Survey (Brown et al. 2012b), which together formed the basis for our WD observations. Solid blue diamonds are the previously published clean ELM sample with 17 < g < 19 mag minus J0935+4411, a P = 20 min ELM WD binary that was not selected by color .
1 Gentile Fusillo et al. (2021) high Galactic latitude data quality cut used: (astrometric sigma5D max < 1.5 OR (ruwe ≤ 1.1 AND ipd gof harmonic amplitude < 1)) AND (phot bp n obs > 2 AND phot rp n obs > 2) AND (|phot br rp excess factor corrected| < 0.6) AND (astrometric excess noise sig < 2 OR (astrometric excess noise sig ≥ 2 AND astrometric excess noise < 1.5)) Figure 1. Color-color plot for every star in SDSS with 17 < g < 19 mag. For guidance, magenta lines plot theoretical 0.16, 0.19, and 0.32 M⊙ WD tracks (Althaus et al. 2013). Blue lines mark the approximate color selection of the HVS Survey and ELM Survey; blue diamonds mark previously published ELM WD binaries with 17 < g < 19 mag. Red dashed lines mark the new color selection; red crosses mark new ELM WD candidates, and red stars mark new ELM WD binaries.
The dashed red line in Figure 1 plots the color selection. There are 3094 candidates 17 < g < 19 mag in this color selection region. The large number of candidates reflects the overlap with the main locus of WDs in color space.

ELM WD Candidates: Parallax Selection
We can also select ELM candidates on the basis of Gaia parallax, independent of (u − g) 0 color. Our second approach is to use Gaia parallax to calculate absolute g-band magnitude M g = g 0 − 10 + 5 log(π), where π is in mas and corrected for the global Gaia parallax zero point of −0.017 mas (Lindegren et al. 2021). We do not bother with higher-order zero point corrections because the low mass WD candidates have order-of-magnitude larger ±0.17 mas median parallax errors. We consider the identical magnitude range 17 < g < 19 mag, temperature range −0.4 < (g − r) 0 < −0.1 mag, and require parallax over error > 3.
For reference, Figure 2 plots M g versus (g − r) 0 for every star in SDSS with 17 < g < 19 mag. The lower band of dots are normal WDs, intrinsically faint foreground objects with M g 10 mag. The upper cloud of dots are A spectral-type stars in the halo, intrinsically lu- Absolute g-band magnitude computed from Gaia eDR3 parallax versus de-reddened (g − r)0 color for every star in SDSS with 17 < g < 19 mag. Magenta lines plot theoretical 0.16, 0.19, and 0.32 M⊙ WD tracks (Althaus et al. 2013); discontinuities are where we clip large excursions due to shell flashes in the 0.19 M⊙ track for the sake of clarity. Symbols are the same as before. Red dashed lines mark the Gaia eDR3 selection. Red stars outside the dashed line were selected with Gaia DR2. minous background objects. ELM WDs fall in-between. Magenta lines plot the same Althaus et al. (2013) evolutionary tracks for 0.16, 0.19, and 0.32 M ⊙ WDs as before ( Figure 2). For the sake of clarity, we clip the shell flash loops in the 0.19 M ⊙ track that briefly jump to high temperature and luminosities. Blue stars mark previously published ELM Survey discoveries with 17 < g < 19 mag.
In Gaia eDR3, the median absolute magnitude error at g = 19 mag is 0.14 mag in our (u − g) 0 cut ( Figure 1). The median error is substantially worse (parallax over error = 0.3) outside of our (u − g) 0 cut because of the large number of distant halo stars at similar temperatures. Thus we consider the 0.32 M ⊙ track that sits 1 mag below the faintest ELM Survey discovery, the same track that we used for the color selection. A linear fit, valid over −0.4 < (g − r) 0 < −0.1 mag, yields M g = 11.58 + 5.83(g − r) 0 . We use this line as our M g cut, bounded by the same line 4 mag more luminous (smaller absolute magnitude) in M g , as seen by the dashed red line in Figure 2. There are 192 can- Figure 3. Reduced proper motion Hg computed from Gaia eDR3 proper motions versus de-reddened (g − r)0 color for every star in the SDSS with 17 < g < 19 mag. Symbols are the same as before. Red dashed lines mark a possible selection that we do not use because of its poor efficiency.
didates in this M g selection region, 57 of which are in common with the (u − g) 0 selection.

ELM Candidates: optional Reduced Proper
Motion Selection Parallax provides a powerful but not a perfect selection criteria. Given Gaia measurement errors, WDs more distant than about 1 kpc become confused with halo stars on the basis of parallax alone. Proper motion provides an independent constraint. Given that halo stars are an order of magnitude more distant than ELM WDs, we can use reduced proper motion as an optional third selection criteria.
We compute g-band reduced proper motion H g = g 0 + 5 log µ − 10, where µ is the Gaia proper motion in mas yr −1 . In Gaia eDR3, the median reduced proper motion error at g = 19 mag is 6% within our color cut. The precision is deceptive, however, as the intrinsic velocity dispersion of stars causes overlap in reduced proper motion space.
For reference, Figure 3 plots H g versus (g − r) 0 for every star in SDSS with 17 < g < 19 mag. The lower cloud of dots are normal WDs, nearby objects that exhibit larger proper motions than more distant stars of the same temperature. The upper clouds of dots are A spectral-type stars in the halo. ELM WDs fall approxi- mately in-between the foreground WDs and background halo stars, but with considerable spread. The problem is that known ELM WDs are a mix of disk and halo populations, and so exhibit a large spread in reduced proper motion space. We initially considered an empirical reduced proper motion selection 14 + 12(g − r) 0 < H g < 19 + 12(g − r) 0 that includes previously published discoveries and minimizes the contribution of background halo stars ( Figure  3). There are 5267 candidates in this H g selection. However, the H g selection has substantial overlap with the population of normal WDs, and legitimate ELM WDs are excluded. Given the inefficiency of the H g selection compared to the other selection criteria, we do not consider it further.

Final Selection
Parallax selection provides the most efficient means for identifying low mass WDs candidates with minimal contamination from foreground and background sources, followed by color selection. We combine the parallax and color selection to optimize our observing program.
Fifty seven objects jointly satisfy the Gaia eDR3 parallax and SDSS color selection criteria with 17 < g < 19 mag. These objects are the highest-probability ELM WD candidates in the SDSS footprint. Table 1 presents the list. We include 3 additional objects (J0130−0530, J0820+4543, and J2151+2730) in Table 1 that were selected and observed using Gaia DR2 but do not satisfy the Gaia eDR3 parallax cut (i.e., the 3 objects that fall just outside selection region in Figure 2). We include them here for completeness.
Twenty one objects are previously published (Brown et al. 2020a,b) and labeled Pub=1 in Table 1. We present spectroscopy for the remaining targets in this paper, plus orbital solutions for 17 new binary systems. One of the new binary systems, J1239−2041, was identified as an ELM white dwarf by Kosakowski et al. (2020), but its orbital parameters were unconstrained. Binaries are labeled Bin=1 in Table 1.

Spectroscopy
Our observational strategy is to obtain a spectrum for each target, and then follow-up likely ELM WDs with time-series spectroscopy. We began observations in late 2018 after the release of Gaia DR2. Follow-up spectroscopy started a year later, in late 2019. We lost the 2020 observing seasons to telescope closures. Observations resumed the second half of 2021.
We obtained spectra for all 39 of the previously unobserved objects at the 6.5m MMT telescope with the Blue Channel spectrograph (Schmidt et al. 1989). We configured the spectrograph with the 832 l mm −1 grating and the 1.25 ′′ slit to provide 3550 -4500Å coverage at 1.2Å spectral resolution. Spectra are flux-calibrated using a nightly standard star observation. All observations are paired with a comparison lamp spectrum for accurate wavelength calibration.
For 7 southern targets, we obtained additional spectra at the 6.5m Magellan telescope with the MagE spectrograph and at the 4.1m SOAR telescope with the Goodman High Throughput spectrograph (Clemens et al. 2004) to improve the time-baseline for radial velocity variables J0806−0716 and J1239−2041 and to validate the stellar atmosphere fits for the other five southern targets. At Magellan we used the 0.85 arcsec slit, providing a resolving power of R = 4800. For SOAR observations, we used the blue camera with the 1 arcsec slit and the 930 lines mm −1 grating resulting in a spectral resolution of 2.5Å. SOAR observations were obtained as part of the NOAO Programs 2019A-0134 and 2021A-0007.
We also obtained time-series spectroscopy of J0101+0401 at the Apache Point Observatory 3.5m Effective temperature versus surface gravity for the final sample of 57 + 3 low mass WD candidates. Cyan box marks our clean ELM WD sample, the region that maximizes the overlap between the data and the ELM WD tracks (magenta lines, Istrate et al. 2016). Symbols same as before.
telescope with the Dual Imaging Spectrograph (DIS). We used the B1200 and R1200 gratings with a 1.5 arcsec slit for a dispersion of 0.6Å per pixel.
We measure radial velocities by cross-correlating the observations using the RVSAO package (Kurtz & Mink 1998) with high signal-to-noise templates of similar stellar parameters and known velocities obtained in the identical spectrograph set-up. Any systematic observational errors (i.e., in the wavelength solution) are common to the targets and the templates, allowing us to measure precise relative velocities. The radial velocities typically have 10 to 15 km s −1 precision.

Stellar Atmosphere Fits
We measure the effective temperature and surface gravity of each target by fitting pure hydrogenatmosphere models to its summed, rest-frame spectrum. The model grid spans 4000 < T eff < 35,000 K and 4.5 < log g < 9.5  and includes the Stark broadening profiles from Tremblay & Bergeron (2009).
We apply the Tremblay et al. (2015) 3D stellar atmosphere model corrections, which are relevant to 10, 000 K objects, when needed. We present the corrected stellar atmosphere parameters for all targets in Table 1.  Istrate et al. (2016) models for the final ELM WD analysis because the tracks are computed for both disk and halo progenitors. We see both disk and halo objects in the observations, as discussed below. The cyan box is the region we refer to as the clean ELM WD sample: the region that maximizes the overlap between the ELM WD models and the data. Figure 4 shows that our choice of −0.4 < (g − r) 0 < −0.1 mag color is effective in selecting 9000 < T eff < 22,000 K objects and excluding stars at cooler temperatures. New objects (red symbols) are mostly found at both higher and lower log g compared to the original ELM Survey, consistent with the expanded color selection, but nine of the new objects are found in the clean ELM region. We explore these objects further below.
Interestingly, the overall distribution of points in Figure 4 does not uniformly fill T eff and log g space. The onset of hydrogen shell flashes in the WD evolutionary tracks around ∼0.2 M ⊙ creates a gap in the tracks in T eff and log g space. Observed WDs largely overlap the WD evolutionary tracks in T eff and log g space, and do not uniformly fill the plot.
The joint parallax and color criteria are very efficient at finding ELM WDs: half (29) of our targets fall in the clean ELM WD region that directly overlap the ELM WD tracks in the range 5.5 < log g < 7.2.

Disk and Halo Populations
Before selecting an evolutionary track for the analysis, we must first assess whether an observed target is a disk or halo object. We use 3 dimensional space velocities to make this assessment. In brief, we compute space velocities using our systemic radial velocities, corrected for the gravitational redshift of the observed WD, and Gaia proper motions. We then compare the velocity of the low mass WD with the velocity ellipsoid of the disk and halo as described in Brown et al. (2020b). The result is essentially the same as found by Brown et al. (2020b): 70% of the observed objects have motions consistent with the thick disk velocity ellipsoid and 30% have motions consistent with the halo.

White Dwarf Parameters
We then estimate WD mass and luminosity by interpolating the measured T eff and log g through WD evolutionary tracks. For 34 objects with log g < 7.2, we use Istrate et al. (2016) tracks with rotation and diffusion. We apply Z = 0.02 tracks to objects with disk kinematics, and Z = 0.001 tracks to objects with halo kinematics. We then estimate distances from the de-reddened SDSS apparent magnitude and the absolute magnitude derived from the tracks, d spec = 10 ((g0−Mg )/5−2) kpc. Figure 5 compares our d spec estimates derived from Istrate et al. (2016) tracks against the d geo = 1/π distances measured geometrically by Gaia. We plot the 29 objects with distance ratio uncertainties better than 30%. The mean ratio is d spec /d geo = 1.02 ± 0.04. Thus the Istrate et al. (2016) tracks provide a remarkably accurate measure of mass and radius for ELM WDs.
The most discrepant point in Figure 5 is one of the new ELM WD binaries, J2104+1712, with d spec /d geo = 0.69 ± 0.15. J2104+1712 is a P = 5.7 hr binary. Our spectra show no evidence for a double-lined spectroscopic system or emission in the cores of the Balmer lines, however it is possible that such effects could artificially inflate our log g measurement. If the actual log g were 0.5 dex lower (the WD had larger radius), the derived absolute magnitude would be ∼1 mag more luminous for the same temperature and bring the spectroscopic distance into agreement with the geometric distance.
For the remaining objects that do not overlap the Istrate et al. (2016) tracks, we estimate WD parameters using Althaus et al. (2013) helium-core WD tracks for objects up to log g = 7.5, and Tremblay et al. (2011) WD tracks for objects above log g = 7.5. We note that d spec estimates for the high-gravity WDs are systematically smaller than the Gaia distances, very similar to J2104+1712, but we defer a detailed investigation. We focus our attention on the well-measured ELM WDs.

Binary Orbital Solutions
Follow-up time-series spectroscopy reveals significant 200 km s −1 to 1100 km s −1 radial velocity variability in every 5 < log g < 7.2 target. Velocity variability in log g > 7.5 targets is not significant given our observational errors. For the new targets with 5 < log g < 7.2, we acquire a total of 10 to 20 spectra spaced over a few different nights. This sampling is sufficient to determine binary orbital elements for 14 targets, plus another 3 targets at slightly higher log g. Only one object in the clean ELM region, J1240−0958, lacks sufficient coverage to determine its orbital parameters. Table 2 presents the 17 binaries, and Figure 6 presents the radial velocities phased to the best fit orbital solution. Our computational approach is to minimize χ 2 for a circular orbit fit to each set of radial velocities. We estimate errors by re-sampling the radial velocities with their errors and re-fitting the orbital solution 10,000 times. We report the 15.9% and 84.1% percentiles of the parameter distributions in Table 2. Systemic velocities are corrected for gravitational redshift using the WD parameters in Table 2. One target, J0745+2104,   J1239−2041 is the shortest period (P = 22.5 min) binary in the new discoveries, and among the shortest period detached WD binaries known. Ultra-compact P < 1 hr WD binaries commonly exhibit ellipsoidal variation, reflection effects, and eclipses in their light curves, all of which constrain the inclination of the binary (e.g. Hermes et al. 2014;Burdge et al. 2020b).
We used Gemini South GMOS to acquire time-series imaging of J1239−2041 on UT 2021 June 28 as part of the program GS-2021A-DD-108. We obtained 161 back-to-back g-band exposures over a 66 min baseline to test for photometric variability. The short time-series is sufficient to rule out the presence of eclipses (Figure 7).
A frequency analysis finds marginally significant photometric variability at the radial velocity orbital period. We perform a simultaneous, non-linear least-squares fit that includes the amplitude of the Doppler beaming (sin φ), ellipsoidal variations (cos 2φ), and reflection (cos φ). The best-fit model, plotted as the red line in Figure 7, has relativistic beaming = 0.45 ± 0.12%, ellipsoidal variation = 0.53 ± 0.11%, and reflection effect = 0.00 +0.04 −0.0 %. We also use LCURVE to fit the mass ratio, scaled radii, effective temperatures, and velocity scale used to account for Doppler beaming and gravitational lensing.
We use fixed quadratic limb-darkening values and fixed gravity-darkening values (Claret et al. 2020) based on our spectroscopic atmospheric parameters for the ELM WD and a putative 9000 K, log g=8 companion. We place Gaussian priors on the ELM WD temperature, mass, and velocity semi-amplitude based on our spectroscopic and radial velocity measurements. Our MCMC algorithm uses 256 walkers, each performing 1000 fit iterations. We remove the first 100 of each as burn-in, leaving 230,400 iterations from which we created the final parameter distributions plotted in Figure 8. Note that LCURVE refers to the parameters of the unseen companion -the more massive star -with a subscript 1 and the ELM WD with a subscript 2.
We derive a 0.68 +0.11 −0.06 M ⊙ companion mass given the fitted mass ratio and our spectroscopic 0.291 ± 0.013 M ⊙ ELM WD mass. As seen in Figure 8, the companion's 12500 +10100 −7100 K temperature and 0.007 +0.005 −0.003 R ⊙ radius are poorly constrained by our short time-series observations, but the parameters are consistent with WD tracks. LCURVE returns a radius of 0.031 +0.001 −0.001 R ⊙ for the ELM WD in perfect agreement with input priors.
The lack of eclipses limits the inclination of J1239−2041 to i < 81 • . However, the 0.53% ellipsoidal amplitude implies an approximately edge-on orientation. The LCURVE fit yields i = 71 +8 −10 • . We also detect Doppler beaming, and the predicted 0.38% beaming amplitude is consistent within the errors of the observed 0.45% value. A longer baseline light curve is needed to more accurately constrain J1239−2041's parameters.

Observational Completeness for ELM WDs
A primary goal of this paper is to assess the observational completeness of the published ELM Survey sample. We search for probable low mass WD candidates by expanding the original color selection region by the 3-σ SDSS color uncertainties, and then applying Gaia parallax constraints to identify 57 candidates with 17 < g < 19 mag across the entire SDSS footprint in the color range −0.4 < (g − r) 0 < −0.1 mag.
Follow-up observations are complete. Twenty one candidates were previously published, and observations for the remaining 36 candidates are presented in the previous section. All but one of the 57 candidates are WDs. (We consider J2338+1235 a non-WD on the basis of its spectroscopic log g = 4.44 ± 0.13. J2338+1235 has one of the poorest parallax constraints in our sample, parallax over error = 3.8, and its −299 ± 11 km s −1 heliocentric radial velocity suggests it is a halo star.) In terms of target selection, the ELM Survey missed a dozen low mass WD candidates with −0.4 < (g − r) 0 < −0.1 mag in SDSS (Figure 1). Nearly all of the missing candidates are in blue half defined by the HVS survey (Brown et al. 2012b) and are missing because they failed the selection criteria at the time of observation. The HVS Survey was executed prior to SDSS DR8 when the SDSS photometric calibration and sky coverage changed every data release. Candidates near the hot edges of the color selection in Figure 1 were particularly sensitive to the changing calibrations. J1313+5828 had (g − r) 0 < −0.40 in DR7 and earlier data releases, for example, and so was not included in the original sample. Five candidates located at southern Galactic latitudes or low declination did not appear in SDSS sky coverage until DR8. Four more candidates were actually observed by the HVS Survey, but were not included in ELM Survey follow-up because the first spectra indicated log g ≤ 5 surface gravity and the objects were not considered probable WDs.
In terms of sample completeness for ELM WDs, we estimate that the published ELM Survey was 90% complete. Our observations identify a total of 9 new ELM WDs, 10 if you count J2306+0224 just outside the clean ELM WD region in Figure 4. Only three of the new ELM WDs fall in the original color selection criteria: a 9% increase in the clean sample of ELM WDs with 17 < g < 19 mag and −0.4 < (g − r) 0 < −0.1 mag. This color-selected sample should now be complete over the entire SDSS footprint.
Of course ELM WDs exist outside our original color selection criteria, as demonstrated by the broader se-lection used in this paper. Yet most of the candidates from the expanded color selection in Figure 1 are normal WDs. Torres et al. (2022) estimate that only 1% of all WDs are He-core objects based on the Gaia 100 pc sample. Applying a joint parallax constraint is thus important for identifying real ELM WDs in a crowded region of parameter space. It is interesting to explore how other Gaia-based catalogs of WDs and ELM WD candidates compare.  (2019) Gaia DR2 ELM WD candidate catalog. We restrict the figure to objects within SDSS with 17 < g < 19 mag to make a one-to-one comparison.

Comparison with Gaia WD catalogs
Gentile Fusillo et al. (2021) select 1.3 million WD candidates from the Gaia eDR3 catalog from which they identify 390,000 high-confidence P WD > 0.75 sources. The 21,500 high-confidence sources with SDSS 17 < g < 19 mag are drawn as black dots in Figure 9. The majority of sources fill the locus of normal WDs.
(2021) photometric hydrogen atmosphere parameter estimates are not very reliable for our class of objects, however, on average 7 ± 16% lower in T eff and 0.4 ± 0.5 dex lower in log g than our spectroscopic values. Furthermore, the 20 objects in our sample classified as low- While there is clear overlap between our ELM WD selection and existing Gaia-based catalogs, it is less clear how to compare the underlying selection criteria. The Gentile Fusillo et al. (2021) catalog is focused on normal log g ∼ 8 WDs, and so may provide a better means to exclude normal WDs (and possibly background halo stars, which should not exist in the WD catalog) from a search for ELM WDs than the other way around. Pelisoli & Vos (2019) intentionally target ELM WDs but with low efficiency. If we restrict the Pelisoli & Vos (2019) catalog to candidates with SDSS 17 < g < 19 mag and −0.4 < (g − r) 0 < −0.1 mag -our ELM WD sample should be essentially complete in this magnitude and color range -then 4.5% (16 of 357) of the Pelisoli & Vos (2019) candidates are ELM WDs with 5.5 < log g < 7.2. Figure 9 suggests that the remaining candidates are likely WDs with higher >0.3 M ⊙ masses. We conclude that the combination of multi-passband photometry and parallax makes for a more efficient selection of ELM WDs.

New LISA Verification Binaries
One reason ELM WDs are interesting is that they reside in ultra-compact binaries. Nearly two-thirds (35 of 57) of the candidates in the current sample are detached, double-degenerate binaries. Two of the newly discovered binaries have orbital periods less than 32 minutes and are potential LISA gravitational wave sources.
J1239−2041 was first identified as a 0.29 M ⊙ ELM WD by Kosakowski et al. (2020) in ELM Survey South observations. Follow-up motivated by the object's SDSS colors and Gaia parallax reveal it to be a P = 22.514 ± 0.186 min binary. The K = 557.2 ± 10.4 km s −1 semiamplitude requires a M 2 > 0.61 M ⊙ minimum mass companion at ≃0.25 R ⊙ separation. This is a double degenerate binary, one that will merge in < 4 Myr due to loss of energy and angular momentum via gravitational wave radiation. Located at a heliocentric distance of d = 0.89 ± 0.11 kpc, the NASA LISA Detectability Calculator 2 estimates LISA will detect J1239−2041 with a 4-yr SNR=12.
J2322+2103 is a 0.25 M ⊙ ELM WD in a P = 31.971 ± 0.361 min binary. Its K = 248.1 ± 4.3 km s −1 semi-amplitude requires a M 2 > 0.19 M ⊙ minimum mass companion with a < 30 Myr gravitational wave merger time. Located at a heliocentric distance of d = 1.04 ± 0.16 kpc, J2322+2103 has an estimated LISA 4-yr SNR=3.5. Figure 10 plots the distribution of characteristic strain versus gravitational wave frequency, f GW = 2/P , for the current sample (solid symbols) and all other ELM Survey binaries (open symbols) compared to the LISA 4-yr sensitivity curve (Robson et al. 2019). New binaries are red, previously published binaries are blue. The overall distribution reflects the fact that these detached WD binaries look like constant frequency sources to LISA, so their SNR is proportional to √ N orbits in 4 yr. Burdge et al. (2020b) target similar sources with P < 1 hr variability in the Zwicky Transient Factory and have discovered even shorter period systems (Burdge et al. 2019(Burdge et al. , 2020a. Four of the 35 binaries in the current sample are likely LISA sources: J1239−2041 and J2322+2103, plus J1235+1543 (Kilic et al. 2017) and J2322+0509 (Brown et al. 2020a). All four are labeled in Figure 10. Interestingly, their LISA SNR may be underestimated. Recent work shows that if the observed metallicitydependence of the stellar binary fraction is applied to WD models (Thiele et al. 2021), or if the separation distribution of observed low mass WD binaries is applied to all classes of WDs (Korol et al. 2022), then the predicted gravitational wave foreground may be reduced by a factor of 2 for resolved LISA sources. Regardless, ELM WD binaries are a significant source of LISA "verification binaries," known binaries observable with both light and gravitational waves.

CONCLUSIONS
We use Gaia parallax and SDSS multi-passband photometry to identify ELM WD candidates. A joint parallax and color selection yields 57 ELM WD candidates with 17 < g < 19 mag over the entire SDSS footprint in the color range −0.4 < (g − r) 0 < −0.1 (approximately 9,000 K < T eff < 22,000 K). Observations presented here complete the sample.
The joint target selection is highly efficient for finding ELM WDs: 50% of the candidates fall in the clean ELM WD region 5.5 < log g < 7.2 that overlaps the ELM WD evolutionary tracks. All of the clean sample of ELM WDs are detached, double-degenerate binaries with P < 1 day orbits.
The 17 binaries published here bring the ELM Survey sample to 124 WD binaries. We find that 10% of low mass WD binaries have P < 1 hr orbits relevant to LISA. This implies that doubling the sample of low mass WDs, for example by surveying the southern sky with SkyMapper as we have done with SDSS in the north, should yield a dozen new multi-messenger WD binaries observable with LISA.
We thank M. Alegria, B. Kunk, E. Martin, and A. Milone for their assistance with observations obtained at the MMT telescope, and Y. Beletsky for their assistance with observations obtained at Magellan telescope. We thank the referee for constructive feedback that improved this paper. This research makes use the SAO/NASA Astrophysics Data System Bibliographic Service. This work was supported in part by the Smithsonian Institution and in part by the NSF under grant AST-1906379. The Apache Point Observatory 3.5-meter telescope is owned and operated by the Astrophysical Research Consortium.
Based on observations obtained at the Southern Astrophysical Research (SOAR) telescope, which is a joint project of the Ministério da Ciência, Tecnologia e Inovações (MCTI/LNA) do Brasil, the US National Science Foundation's NOIRLab, the University of North Carolina at Chapel Hill (UNC), and Michigan State University (MSU).
Based on observations obtained at the international Gemini Observatory, a program of NSF's NOIRLab, which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation. on behalf of the Gemini Observatory partnership: the National Science Foundation (United States), National Research Council (Canada), Agencia Nacional de Investigación y Desarrollo (Chile), Ministerio de Ciencia, Tecnología e Innovación (Argentina), Ministério da Ciência, Tecnologia, Inovações e Comunicações (Brazil), and Korea Astronomy and Space Science Institute (Republic of Korea).