NewReferences

# A dense $\mathbf{0.1M_{\rm\odot}}$ star in a 51-minute orbital period
eclipsing binary

Kevin B. Burdge1,2 ∗ Kareem El-Badry3,4,5 Thomas R. Marsh6 Saul Rappaport1,2
Warren R. Brown3 Ilaria Caiazzo7 Deepto Chakrabarty1,2 V. S. Dhillon8,9 Jim
Fuller7 Boris T. Gänsicke6 Matthew J. Graham7 Erin Kara1,2 S. R. Kulkarni7 S.
P. Littlefair8 Przemek Mróz10 Pablo Rodríguez-Gil9,11 Jan van Roestel7 Robert
A. Simcoe1,2 Eric C. Bellm12 Andrew J. Drake7 Richard G. Dekany13 Steven L.
Groom14 Russ R. Laher14 Frank J. Masci14 Reed Riddle13 Roger M. Smith13 &
Thomas A. Prince7

###### Abstract

In over a thousand known cataclysmic variables (CVs), where a white dwarf is
accreting from a hydrogen-rich star, only a dozen have orbital periods below
75 minutes[1, 2, 3, 4, 5, 6, 7, 8, 9]. One way to achieve these short periods
requires the donor star to have undergone substantial nuclear evolution prior
to interacting with the white dwarf[10, 11, 12, 28, 14], and it is expected
that these objects will transition to helium accretion. These transitional CVs
have been proposed as progenitors of helium CVs[28, 14, 15, 16, 17, 27].
However, no known transitional CV is expected to reach an orbital period short
enough to account for most of the helium CV population, leaving the role of
this evolutionary pathway unclear. Here we report observations of ZTF
J1813+4251, a 51-minute orbital period, fully eclipsing binary system
consisting of a star with a temperature comparable to that of the Sun but a
density 100 times greater due to its helium-rich composition, accreting onto a
white dwarf. Phase-resolved spectra, multi-band light curves and the broadband
spectral energy distribution allow us to obtain precise and robust constraints
on the masses, radii and temperatures of both components. Evolutionary
modeling shows that ZTF J1813+4251 is destined to become a helium CV binary,
reaching an orbital period under 20 minutes, rendering ZTF J1813+4251 a
previously missing link between helium CV binaries and hydrogen-rich CVs.

Department of Physics, Massachusetts Institute of Technology, Cambridge, MA
02139, USA

Kavli Institute for Astrophysics and Space Research, Massachusetts Institute
of Technology, Cambridge, MA 02139, USA

Center for Astrophysics, Harvard & Smithsonian, 60 Garden Street, Cambridge,
MA 02138

Harvard Society of Fellows, 78 Mount Auburn Street, Cambridge, MA 02138

Max-Planck Institute for Astronomy, Könighstuhl 17, D-69117 Heidelberg,
Germany

Department of Physics, University of Warwick, Coventry CV4 7AL, UK

Division of Physics, Mathematics and Astronomy, California Institute of
Technology, Pasadena, CA, USA

Department of Physics & Astronomy, University of Sheffield, Sheffield S3 7RH,
UK

Instituto de Astrofísica de Canarias, Vía Láctea s/n, La Laguna, E-38205
Tenerife, Spain

Astronomical Observatory, University of Warsaw, Al. Ujazdowskie 4, 00-478
Warszawa, Poland

Departamento de Astrofísica, Universidad de La Laguna, E-38206 La Laguna,
Tenerife, Spain

DIRAC Institute, Department of Astronomy, University of Washington, 3910 15th
Avenue NE, Seattle, WA 98195, USA

Caltech Optical Observatories, California Institute of Technology, Pasadena,
CA 91125, USA

IPAC, California Institute of Technology, Pasadena, CA, USA

Using data from the Zwicky Transient Facility (ZTF)[19], we recently conducted
a systematic period search of 1,220,038,476 unique sources using a graphics
processing unit (GPU) based implementation of the generalised Lomb-Scargle
algorithm[1], as part of an ongoing campaign to identify and characterize
short period astronomical variables[3] (Methods). This effort led to the
identification of ZTF J1813+4251, which exhibited strong ellipsoidal
variations due to the tidal deformation of a star[22] and an eclipse, as seen
in the lightcurves shown in Figure 1. The orbital period of this source is
just 51 minutes, despite a spectral energy distribution indicating that it
hosts a $~{}6000\rm\,K$ star of spectral type F. The _Gaia_ eDR3 astrometric
solution suggests that the source is nearly a kiloparsec away, making it
significantly less luminous than an F type main-sequence star given its
apparent brightness. This ultracompact binary system had been excluded from
previous searches for ultracompact binaries[3] because of its red color
relative to other ultracompact systems.

On June 4 and July 5 2021, we obtained phase-resolved spectroscopic
observations of ZTF J1813+4251 using the Low Resolution Imaging Spectrometer
(LRIS)[23] on the 10-m W. M. Keck I Telescope on Mauna Kea (Methods). Figure 2
illustrates these observations, revealing a late F type stellar spectrum
blanketed with metallic absorption lines characteristic of main sequence
stars. These absorption lines Doppler shift with a large velocity semi-
amplitude of $461.3\pm 3.4\rm\,km\,s^{-1}$ due to the short orbital period of
the binary. Additionally, double-peaked emission lines of hydrogen, helium,
and calcium reveal an accretion disk, most prominently visible in the red half
of the spectrum, as seen in Figure 3.

On June 12, 2021, we obtained high speed images of ZTF J1813+4251 with the
quintuple-beam high speed photometer HiPERCAM[5] on the 10.4-m Gran Telescopio
Canarias (GTC) on La Palma. This multi-color high signal-to-noise (SNR) light
curve, illustrated in Figure 1, revealed that the depth of the eclipse of the
white dwarf varies as a function of wavelength, with a depth of $\sim 50$
percent in the HiPERCAM-$u_{\mathrm{s}}$ filter centered at
$3526\,\mathrm{\AA}$, but only a $\sim 10$ percent depth in the
HiPERCAM-$z_{\mathrm{s}}$ filter centered at $9156\,\mathrm{\AA}$. This
dependence of the eclipse depth on wavelength indicates that the accreting
white dwarf contributes a larger share of the luminosity at short wavelengths,
due to it having a substantially higher surface temperature
($T_{\mathrm{WD}}=12600\pm 500\,\rm K$) than the donor star
($T_{\mathrm{Donor}}=6000\pm 90\,\rm K$). There are two photometric maxima per
orbit due to ellipsoidal variations, and the high signal-to-noise (SNR)
HiPERCAM light curve also revealed a pronounced “O’Connell” effect[25], in
which the peak flux of these alternating maxima is substantially different due
to a modulation component at the orbital frequency. This effect may be due to
a spot on the donor star, which is tidally locked, and thus completes a
rotation every orbital period.

The accretion disk seen in panels (d) and (e) of Figure 3 indicates that the
donor star in ZTF J1813+4251 is transferring matter to the white dwarf, and
thus is filling its Roche-lobe. The Roche-lobe has a scale-invariant geometry
dependent only on the mass ratio of the system, and thus the geometry of the
donor depends only on this mass ratio[21]. In a binary system undergoing a
total eclipse, the time between mid-ingress and mid-egress depends only on the
geometry of the donor and the inclination of the system. This means that the
eclipses in ZTF J1813+4251’s light curve constrain the system to a unique mass
ratio vs inclination relation[15]. Additionally, the depths of the primary and
secondary eclipse constrain a unique radius ratio between the white dwarf and
that of the donor, which, when combined with the overall geometry of the
primary eclipse, yields a unique solution for the radii divided by the semi-
major axis (the scaled radii), orbital inclination, and thus mass ratio,
because the system is Roche lobe filling. These constraints, together with the
precise donor radial velocity semi-amplitude measured from the spectra (see
the coadded spectrum in Figure 2), allow for a robust determination of the
system parameters on the basis of only Roche lobe geometry and Kepler’s laws
(Methods). We report these system parameters in Table 1. Because the donors in
transitional CVs are still transferring significant hydrogen and have not
switched to primarily helium accretion yet, as seen in a helium CV (also known
as AM CVns), these objects often exhibit a mixture of hydrogen disk lines
along with unusually strong helium lines in their spectra as a result of their
slowly transitioning to transferring mainly helium, and as result, we see in
panel (d) of Figure 3 that the disk also exhibits significant helium emission.

A Niels Gehrels Swift observatory (Swift) X-ray telescope (XRT) observation of
the system revealed no detectable X-ray flux, with the X-ray luminosity of the
system, $L_{\mathrm{x}}$, constrained to a $3\sigma$ upper limit of
$L_{\mathrm{x}}<1.22\times 10^{31}\,\rm erg\,s^{-1}$ given an assumed distance
of $891^{+174}_{-135}\,\mathrm{pc}$ inferred from _Gaia_ EDR3[28], an upper
limit consistent with the typical x-ray luminosity of a non-magnetic
cataclysmic variable[29] at these energies. The source was detected with the
Swift Ultraviolet Optical Telescope in the ultraviolet, constraining the
temperature of the accreting white dwarf (Methods).

As seen in Figure 3, archival ZTF data captured an outburst of the object on
September 21, 2019, brightening by over a factor of two relative to its
quiescent brightness in the ZTF-$g$ and ZTF-$r$ filters. The final observation
of the night during which the outburst was detected took place during the
totality of the primary eclipse, and as a result the flux was significantly
attenuated, indicating that this brightening must have originated from the
accretor. Because the white dwarf was completely eclipsed during this phase,
but the flux was still partially elevated relative to quiescence, we can infer
that the luminosity must be originating from the accretion disk, which was not
completely occulted at this orbital phase.

We used Modules for Experiments in Stellar Astrophysics[25] (MESA) to model
the past and future evolution of systems like ZTF J1813+4251, and Figure 4
illustrates these evolutionary tracks. Currently, ZTF J1813+4251 exhibits an
elevated temperature because the white dwarf is rapidly removing matter from
the donor as the orbital period shrinks, causing its radius to decrease at
roughly constant luminosity. We expect to be able to detect this orbital
period evolution over the next decade, as using the HiPERCAM light curve we
are able to measure the mid-eclipse time, $T_{0}$, to a precision of $\approx
0.16\,\rm s$, which will allow us to test whether the angular momentum loss in
the binary deviates from that expected purely due to general relativity.
Assuming that the system is evolving purely according to general relativistic
orbital decay, the eclipse time will deviate from that expected if the orbital
period were constant by about $7\,\rm s$ over the next ten years[31]. However,
evolutionary models (Methods) indicate that magnetic braking may still be
responsible for as much as two thirds of the orbital angular momentum loss in
the system, and there are other effects which may influence the orbital
evolution as well such as mass-transfer. As seen in panels (a) and (d) of
Figure 4, when the system evolves to shorter orbital periods, the last
remnants of hydrogen will be removed from the donor, eliminating its fuel
source, and its helium rich remnant will begin to transition to being
supported by electron degeneracy pressure. At this stage, which our models
predict will occur in approximately 75 million years, ZTF J1813+4251 will have
completed its transition to a helium CV system, reaching a period minimum
around 18 minutes. Due to it being supported by electron degeneracy pressure,
ZTF J1813+4251’s structure will undergo an important transition, in which it
will begin to expand in response to mass loss, rather than shrink. This will
cause further accretion to decrease the donor’s density as the system evolves
to longer orbital periods, and the donor will rapidly cool as a result of
adiabatic expansion in response to mass loss. We estimate that after the
period minimum, the system will spend three hundred million years evolving
back out to an orbital period of about 30 minutes as a helium CV, thus
spending a much larger fraction of its remaining life in this state than in
its current short-lived orbital decay phase as a transitional CV.

The discovery of ZTF J1813+4251 demonstrates that evolved CV donors can reach
orbital periods below 20 minutes, and likely form helium CVs. Other known
examples of candidate transitional CVs[2, 3, 4, 5, 7, 8, 9] below the period
minimum are either dominated entirely by the flux of an accretion disk, making
it difficult to characterize the properties of the donor star, or exhibit a
detectable donor which is not evolved enough to reach orbital periods short
enough to account for the full range of periods seen in the helium CV
population. Additionally, none of the known transitional CVs exhibit a clean
radial velocity semi-amplitude paired with a detectable full primary and
secondary eclipse like that seen in ZTF J1813+4251, making ZTF J1813+4251 the
only transitional CV with robustly characterized system parameters. As seen in
panel A of Figure 4, the most important property in predicting how evolved the
donor is – and thus what period minimum it can reach – is its effective
temperature. As seen in panel A of figure 4, we also modelled EI Psc[7], which
is one of the few transitional CVs with enough parameter estimates to
construct such models. Like ZTF J1813+4251, we expect this system to evolve
into a helium CV, however it reaches a period minimum of around 30 minutes
rather than below 20 minutes. ZTF J1813+4251, which has a shorter orbital
period and significantly hotter donor than any currently known transitional CV
candidate below the period minimum, demonstrates that the donors of these
systems can reach much more extreme states than seen in previously identified
systems, and solidifies the role of this evolutionary channel in the formation
of helium CVs, including those with orbital periods short enough to emit
significantly in the LISA gravitational wave band. ZTF J1813+4251, with its
well-characterized components, will serve as an anchor point for models of
binary evolution connecting the population of hydrogen accreting CVs with
evolved donors at longer orbital periods and the population of ultracompact
helium CVs.

Figure 1: a) The five color HiPERCAM light curve of ZTF J1813+4251, with filters $u_{\mathrm{s}}$, $g_{\mathrm{s}}$, $r_{\mathrm{s}}$, $i_{\mathrm{s}}$ and $z_{\mathrm{s}}$ arranged from top to bottom. The object exhibits a strong sinusoidal component at twice the orbital frequency due to the tidal deformation of the donor star. The system undergoes a full eclipse, in which the donor fully occults the accreting white dwarf. The depth of the eclipse varies dramatically with wavelength, with the white dwarf contributing half the luminosity in the $u_{\mathrm{s}}$ filter, while only about ten percent in the $z_{\mathrm{s}}$ filter, primarily because the white dwarf preferentially emits its radiation at shorter wavelengths than the cooler donor. b) Best fit light curve models of the primary and secondary eclipses of ZTF J1813+4251. These models, combined with the spectroscopically derived radial velocity semi-amplitude, allow us to robustly constrain the properties of both components in the system. Figure 2: a) A sinusoidal fit to the measured radial velocities of the donor star in ZTF J1813+4251, with a best-fit velocity semi-amplitude of $K_{2}=461.3\pm 3.4\rm\,km\,s^{-1}$ and systemic velocity of $\gamma=-36.8\pm 2.4\rm\,km\,s^{-1}$. b) The blue side of the Keck LRIS spectrum, coadded in the rest frame of the donor star. The spectrum exhibits a large number of narrow metal lines, including those from calcium, iron and magnesium, characteristic of main sequence stars with spectral types similar to the Sun. c) The trailed blue Keck LRIS spectra, illustrating the significant wavelength shifts in the spectrum caused by the large Doppler shifts of the donor. The trailed spectra reveal these Doppler shifts occur in the large number of narrow lines associated with the donor spectrum. Figure 3: a) The ZTF light curve of ZTF J1813+4251. The system has been in quiescence in all ZTF observations, except for those taken on September 21, 2019, in which the system brightened by over a factor of two. b) A zoom in of the ZTF data taken on September 21, 2019, during the outburst. The final of the four measurements during the outburst is significantly fainter because it occurred during the primary eclipse of the system. c) The phase-folded ZTF data, illustrating that the faint point which occurred during the outburst coincided with the primary eclipse in the system. d) The red side of the Keck LRIS spectrum, coadded in the rest frame of the accreting white dwarf. Double-peaked emission lines associated with the blueshifted and redshifted components of an accretion disk are clearly visible. This disk exhibits emission associated with several elements, including hydrogen, helium and calcium. e) The profile of the hydrogen H$\alpha$ emission feature, the strongest disk feature in the optical spectrum. Figure 4:  Evolutionary tracks of ZTF J1813+4251 generated using Modules for Experiments in Stellar Astrophysics (MESA). a) The donor’s temperature evolution as a function of orbital period in transitional CVs like ZTF J1813+4251. The wide dashed orange line, dotted red line and alternating dashed blue line indicate MESA evolutionary sequences corresponding to a transitional CV where mass loss from the donor commenced at 97, 95 and 94 percent of its main sequence lifetime, and the black star indicates the position of ZTF J1813+4252 on the red dotted track. These sequences host warmer donors than normal hydrogen CVs, whose location are indicated on the right of the plot, and evolve into the region occupied by helium CVs in the lower left corner. For comparison, we also illustrate a track for El Psc[7], another transitional CV candidate which will likely reach a period minimum of just under 30 minutes. b) The future orbital period evolution of the MESA tracks. The track corresponding to ZTF J1813+4251 will reach a period minimum of about 18 minutes in approximately 75 million years, and will spend the next 300 million years evolving back out to a period of about half an hour as a helium CV. c) The evolution of the donor mass as a function of orbital period, reaching just a few hundredths of a solar mass as the tracks evolve out to longer orbital periods as helium CVs. d) The hydrogen remaining in the donor over the course of its future evolution. At around the period minimum, the donor star loses all remaining hydrogen. (See also Kalomeni et al. (2016)[28] for a comprehensive discussion of the evolution of these systems.) $\rm Right\,\,ascension$ | $18^{\mathrm{h}}\,\,13^{\mathrm{m}}\,\,11.13^{\mathrm{s}}$
---|---
$\rm Declination$ | $+42^{\circ}\,\,51^{\prime}\,\,50.4^{\prime\prime}$
$\rm Proper\,\,motion\,\,in\,\,right\,\,ascension$ | $-12.32\pm 0.19\,\mathrm{mas\,yr^{-1}}$
$\rm Proper\,\,motion\,\,in\,\,declination$ | $-2.66\pm 0.19\,\mathrm{mas\,yr^{-1}}$
$\mathrm{Parallax}\,\,(\sigma)$ | $1.20\pm 0.16\,\mathrm{mas}$
$\mathrm{Distance}\,\,(D)$ | $891^{+174}_{-135}\,\mathrm{pc}$
$\mathrm{Systemic\,\,velocity}\,\,(\gamma)$ | $-36.8\pm 2.4\,\mathrm{km\,s^{-1}}$
$\mathrm{Orbital\,\,period}\,\,(P_{\mathrm{b}})$ | $3069.64398\pm-0.00015\,\mathrm{s}$
$\mathrm{Orbital\,\,period\,\,derivative}\,\,(\dot{P_{\mathrm{b}}})$ | $-0.29\pm 1.15\times 10^{-12}\,\mathrm{s}\,\mathrm{s}^{-1}$
$\mathrm{Time\,\,of\,\,superior\,\,conjunction}\,\,(T_{0})$ | $59377.1538184\pm 0.0000020\,\mathrm{BMJD_{TDB}}$
$\mathrm{Radial\,\,velocity\,\,of\,\,donor}\,\,(K_{\mathrm{Donor}})$ | $461.3\pm 3.4\,\mathrm{km\,s^{-1}}$
$\mathrm{Orbital\,\,inclination}\,\,(i)$ | $78.80\pm 0.18\,\mathrm{deg}$
$\mathrm{Semi\,\,major\,\,axis}\,\,(a)$ | $0.4000\pm 0.0041\,R_{\rm\odot}$
$\mathrm{White\,\,dwarf\,\,mass}\,\,(M_{\mathrm{W}D})$ | $0.562\pm 0.015\,M_{\rm\odot}$
$\mathrm{Donor\,\,mass}\,\,(M_{\mathrm{Donor}})$ | $0.1185\pm 0.0067\,M_{\rm\odot}$
$\mathrm{White\,\,dwarf\,\,radius}\,\,(R_{\mathrm{W}D})$ | $0.01374\pm 0.00023\,R_{\rm\odot}$
$\mathrm{Donor\,\,radius}\,\,(R_{\mathrm{Donor}})$ | $0.1017\pm 0.0019\,R_{\rm\odot}$
$\mathrm{White\,\,dwarf\,\,temperature}\,\,(T_{\mathrm{W}D})$ | $12600\pm 500\,\rm K$
$\mathrm{Donor\,\,temperature}\,\,(T_{\mathrm{Donor}})$ | $6000\pm 80\,\rm K$
Table 1: Table of parameters (The first five parameters are the astrometric
solution reported by _Gaia_ eDR3[32], at Epoch J2016.0 and Equinox J2000.0).

## References

  * [1] Green, M. J. _et al._ Spectroscopic and photometric periods of six ultracompact accreting binaries.  _MNRAS_ 496, 1243–1261 (2020).
  * [2] Kato, T. _et al._ Survey of Period Variations of Superhumps in SU UMa-Type Dwarf Novae.  _PASJ_ 61, S395–S616 (2009).
  * [3] Littlefield, C. _et al._ A New Sub-period-minimum Cataclysmic Variable with Partial Hydrogen Depletion and Evidence of Spiral Disk Structure.  _AJ_ 145, 145 (2013).
  * [4] Ramsay, G. _et al._ BOKS 45906: a CV with an orbital period of 56.6 min in the Kepler field?  _MNRAS_ 438, 789–795 (2014).
  * [5] Augusteijn, T., van der Hooft, F., de Jong, J. A. & van Paradijs, J.  V485 Centauri: a dwarf nova with a 59min orbital period.  _A &A_ 311, 889–900 (1996).
  * [6] Kato, T. _et al._ Survey of period variations of superhumps in SU UMa-type dwarf novae. V. The fifth year (2012-2013).  _PASJ_ 66, 30 (2014).
  * [7] Thorstensen, J. R. _et al._ 1RXS J232953.9+062814: A Dwarf Nova with a 64 Minute Orbital Period and a Conspicuous Secondary Star.  _ApJ_ 567, L49–L52 (2002).
  * [8] Breedt, E. _et al._ CSS100603:112253-111037: a helium-rich dwarf nova with a 65 min orbital period.  _MNRAS_ 425, 2548–2556 (2012).
  * [9] Kato, T. _et al._ Survey of period variations of superhumps in SU UMa-type dwarf novae. VII. The seventh year (2014-2015).  _PASJ_ 67, 105 (2015).
  * [10] Sienkiewicz, R.  The minimum period of hydrogen-deficient cataclysmic binaries.  _Acta Astron._ 34, 325–330 (1984).
  * [11] Rappaport, S. & Joss, P. C.  The lower main sequence and the nature of secondary stars in ultracompact binaries.  _ApJ_ 283, 232–240 (1984).
  * [12] Nelson, L. A., Rappaport, S. A. & Joss, P. C.  The Evolution of Ultrashort Period Binary Systems.  _ApJ_ 304, 231 (1986).
  * [13] Kalomeni, B. _et al._ Evolution of Cataclysmic Variables and Related Binaries Containing a White Dwarf.  _ApJ_ 833, 83 (2016).
  * [14] Podsiadlowski, P., Han, Z. & Rappaport, S.  Cataclysmic variables with evolved secondaries and the progenitors of AM CVn stars.  _MNRAS_ 340, 1214–1228 (2003).
  * [15] Nelemans, G., Yungelson, L. R., van der Sluys, M. V. & Tout, C. A.  The chemical composition of donors in AM CVn stars and ultracompact X-ray binaries: observational tests of their formation.  _MNRAS_ 401, 1347–1359 (2010).  0909.3376.
  * [16] Goliasch, J. & Nelson, L.  Population Synthesis of Cataclysmic Variables. I. Inclusion of Detailed Nuclear Evolution.  _ApJ_ 809, 80 (2015).  1607.06217.
  * [17] Ramsay, G. _et al._ Physical properties of AM CVn stars: New insights from Gaia DR2.  _A &A_ 620, A141 (2018).
  * [18] El-Badry, K., Rix, H.-W., Quataert, E., Kupfer, T. & Shen, K. J.  Birth of the ELMs: a ZTF survey for evolved cataclysmic variables turning into extremely low-mass white dwarfs.  _MNRAS_ (2021).
  * [19] Bellm, E. C. _et al._ The Zwicky Transient Facility: System Overview, Performance, and First Results.  _PASP_ 131, 018002 (2019).
  * [20] Zechmeister, M. & Kürster, M.  The generalised Lomb-Scargle periodogram. A new formalism for the floating-mean and Keplerian periodograms.  _A &A_ 496, 577–584 (2009).
  * [21] Burdge, K. B. _et al._ A Systematic Search of Zwicky Transient Facility Data for Ultracompact Binary LISA-detectable Gravitational-wave Sources.  _ApJ_ 905, 32 (2020).
  * [22] Burdge, K. B. _et al._ Orbital Decay in a 20 Minute Orbital Period Detached Binary with a Hydrogen-poor Low-mass White Dwarf.  _ApJ_ 886, L12 (2019).
  * [23] Oke, J. B. _et al._ The Keck Low-Resolution Imaging Spectrometer.  _PASP_ 107, 375 (1995).
  * [24] Dhillon, V. S. _et al._ HiPERCAM: A quintuple-beam, high-speed optical imager on the 10.4-m Gran Telescopio Canarias.  _MNRAS_ (2021).
  * [25] O’Connell, D. J. K.  The so-called periastron effect in close eclipsing binaries ; New variable stars (fifth list).  _Publications of the Riverview College Observatory_ 2, 85–100 (1951).
  * [26] Paczyński, B.  Evolutionary Processes in Close Binary Systems.  _ARA &A_ 9, 183 (1971).
  * [27] Chanan, G. A., Middleditch, J. & Nelson, J. E.  The Geometry of the Eclipse of a Pointlike Star by a Roche-Lobe Companion.  _ApJ_ 208, 512–517 (1976).
  * [28] Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., Demleitner, M. & Andrae, R.  Estimating Distances from Parallaxes. V. Geometric and Photogeometric Distances to 1.47 Billion Stars in Gaia Early Data Release 3\.  _AJ_ 161, 147 (2021).
  * [29] Mukai, K.  X-Ray Emissions from Accreting White Dwarfs: A Review.  _PASP_ 129, 062001 (2017).
  * [30] Paxton, B. _et al._ Modules for Experiments in Stellar Astrophysics (MESA).  _ApJS_ 192, 3 (2011).
  * [31] Burdge, K. B. _et al._ General relativistic orbital decay in a seven-minute-orbital-period eclipsing binary system.  _Nature_ 571, 528–531 (2019).
  * [32] Gaia Collaboration _et al._ Gaia Early Data Release 3. Summary of the contents and survey properties.  _A &A_ 649, A1 (2021).

K.B.B. is a Pappalardo Postdoctoral Fellow in Physics at MIT and thanks the
Pappalardo fellowship program for supporting his research. The design and
construction of HiPERCAM was funded by the European Research Council under the
European Union’s Seventh Framework Programme (FP/2007-2013) under ERC-2013-ADG
Grant Agreement no. 340040 (HiPERCAM). VSD and HiPERCAM operations are
supported by STFC grant ST/V000853/1. TRM and BTG acknowledge support from the
UK’s Science and Technology Facilities Council (STFC), grant ST/T000406/1.
This project has received funding from the European Research Council (ERC)
under the European Union’s Horizon 2020 research and innovation programme
(Grant agreement No. 101020057).

Based on observations made with the Gran Telescopio Canarias (GTC), installed
at the Spanish Observatorio del Roque de los Muchachos of the Instituto de
Astrofísica de Canarias, on the island of La Palma.

Based on observations obtained with the Samuel Oschin Telescope 48-inch
Telescope at the Palomar Observatory as part of the Zwicky Transient Facility
project. ZTF is supported by the National Science Foundation under Grant No.
AST-1440341 and a collaboration including Caltech, IPAC, the Weizmann
Institute for Science, the Oskar Klein Center at Stockholm University, the
University of Maryland, the University of Washington, Deutsches Elektronen-
Synchrotron and Humboldt University, Los Alamos National Laboratories, the
TANGO Consortium of Taiwan, the University of Wisconsin at Milwaukee, and
Lawrence Berkeley National Laboratories. Operations are conducted by COO,
IPAC, and UW.

Some of the data presented herein were obtained at the W.M. Keck Observatory,
which is operated as a scientific partnership among the California Institute
of Technology, the University of California and the National Aeronautics and
Space Administration. The Observatory was made possible by the generous
financial support of the W.M. Keck Foundation. The authors wish to recognize
and acknowledge the very significant cultural role and reverence that the
summit of Mauna Kea has always had within the indigenous Hawaiian community.
We are most fortunate to have the opportunity to conduct observations from
this mountain.

The authors declare that they have no competing financial interests.

Correspondence and requests for materials should be addressed to K.B.B.
(email: kburdge@mit.edu).

###  0.1 Period Search

We performed a period search using a GPU implementation of the generalised
Lomb-Scargle algorithmZechmeister2009, an analog of the Fourier transform
optimized for data with non-equispaced sampling. This search was based on
newly generated forced point spread function (PSF) photometry performed at the
coordinates of all Pan-STARRS1 sources in ZTF images over a period range of
$250\rm\,days$ down to just $2\rm\,minutes$. By systematically searching all
1,220,038,476 unique sources in the Pan-STARRS1 source catalog[2] with more
than 50 ZTF epochs, we eliminated biases of previous searches relying on
colour or astrometric based selections[3], which represents a major
advancement in probing for short period astrophysical variables across the
northern sky. We searched a total of 1,461,592 trial frequencies per source,
corresponding to an oversampling factor of 2. In total, given that we searched
1,220,038,476 sources, this means a total of $1.78\times 10^{15}$ trial
frequencies were searched, and each source had on average approximately $\sim
10^{3}$ epochs. This search was completed in six weeks on a single desktop
containing four Nvidia 2080 Ti GPUs, and this same desktop was used to perform
the forced PSF photometry of all Pan-STARRs sources in ZTF images.

###  0.2 Keck LRIS observation

We obtained the Keck LRIS observations using the 600/4000 grism as the
dispersive element on the blue arm and the 600/7500 grating on the red arm. We
used $2\times 2$ binning on both channels to reduce both readout noise and
readout time, and used an exposure time of $120\,\rm s$ to mitigate the impact
of orbital Doppler smearing on individual exposures. We reduced the Keck LRIS
observations using a publicly available LPIPE automated data reduction
pipeline[4].

###  0.3 Keck ESI observation

We obtained sixteen $120\,\rm s$ medium-resolution spectra using the
Echellette Spectrograph and Imager (ESI) on the 10-m W. M. Keck II Telescope
on Mauna Kea. These exposures were obtained in the cross-dispersed echelle
mode and were binned $2\times 2$ in order to reduce readout overheads and
readout noise. The data were reduced using the MAuna Kea Echelle Extraction
(MAKEE) pipeline.

###  0.4 HiPERCAM observation

We reduced the HiPERCAM data using the publicly available pipeline[5], and
selected several comparison stars in each filter to ensure accurate absolute
flux calibrations. The pipeline performed aperture photometry with a dynamic
full-width-half-maximum. We used $9,3,3,3$ and $3\,\rm s$ exposures in
$u_{\mathrm{s}}$, $g_{\mathrm{s}}$, $r_{\mathrm{s}}$, $i_{\mathrm{s}}$, and
$z_{\mathrm{s}}$, respectively. The charge coupled devices (CCDs) were
operated with conventional amplifiers using frame-transfer to effectively
eliminate readout overheads. The total observation lasted 150 minutes, or
about three orbital cycles.

###  0.5 Swift XRT observation

We obtained a 2,777 s Swift observation of ZTF J1813+4251 (observation ID
00014355002) in order to constrain the X-ray flux associated with the mass
transfer in the system. We do not detect the source with significance, with
just two counts within an 18 arcsecond aperture, corresponding to a $3\sigma$
upper limit of $0.0039\,\rm cnts\,s^{-1}$ at $0.2-10\,\rm keV$. We used
WebPIMMs to estimate the source flux. Assuming a $3\,\rm keV$ thermal
Bremsstrahlung model characteristic of non-magnetic CVs, and a column density
of $N_{H}=5.58\times 10^{19}\,\rm{cm}^{-2}$ based on the estimated reddening
of $E(g-r)=0.01$, we place a $3\sigma$ upper limit on the $0.2-10\,\rm keV$
X-ray luminosity of $L_{\mathrm{x}}<1.23\times 10^{31}\,\rm erg\,s^{-1}$. We
did not use the NASA HEASARC WebPIMMs column density estimator tool, as this
estimates the column density at infinity, and thus is a significant
overestimate for a nearby source like ZTF J1813+5242, and instead we based the
column density estimate on more detailed dust maps which take into account the
distance to the source[6].

###  0.6 Swift UVOT observation

The 2,777 s Swift observation (observation ID 00014355002) also included data
obtained with the UVOT instrument in the UVW2 filter. We computed an average
apparent magnitude of the source using the uvotsource tool and a 5 arcsecond
aperture centered on the source. Due to a nearby bright star, the background
near the source was significantly elevated, and thus we used an annulus
background aperture centered on the source, with an inner radius of 7.5
arcsecond and an outer radius of 15 arcsecond. We find an apparent magnitude
of $21.20\pm 0.1\mathrm{\,(statistical)}\pm
0.03\mathrm{\,(systematic)}\,m_{\rm AB}$.

###  0.7 Spectral energy distribution (SED) analysis

We estimated the spectroscopic properties of the donor in ZTF J1813+4251 by
fitting atmospheric models to the in-eclipse spectral energy distribution, as
our spectra are composed of contributions from both the white dwarf and the
accretion disk with exposures too long to have obtained a completely in-
eclipse spectrum. We use the in-eclipse apparent magnitudes of the system to
obtain an uncontaminated SED of the donor star, and fit the BT-NextGen (GNS93)
atmospheric models for low mass stars[7] to the SED, allowing the metalicity,
flux scale and temperature to vary, and fixing the surface gravity to
$\log(g)_{\mathrm{D}onor}=5.43\,\rm(cgs)$ based on an estimate derived from an
initial light curve model which used the _Gaia_ eDR3 estimated temperature of
$6088\,\rm K$ for the donor. Using the UltraNest kernel density estimator[8,
9, 10], we sampled over these parameters and obtained a donor temperature
estimate of $T_{\mathrm{Donor}}=6000\pm 80\,\rm K$. The metallicity preferred
the lowest possible value covered by the atmospheric model grids,
$Z_{\mathrm{Donor}}=-1.5\,\rm[Fe/H]$. Because of the evolved nature of the
donor, we consider this metallicity estimate based on main-sequence models to
be dubious, and thus do not report it in the table; we suspect the
uncertainties in the $T_{\mathrm{Donor}}$ are also likely underestimated,
though we are unable to appropriately quantify the model error associated with
using model atmospheres for main-sequence stars on such an extreme evolved
star. We overplot the best fit atmospheric models with the donor and accretor
SEDs in Extended Data Figure 1.

We are also able to precisely measure the contribution of the white dwarf to
the SED because of the eclipses in the light curve, and thus we performed an
SED fit to it as well in order to estimate $T_{\mathrm{W}D}$. A fit to the
optical SED using a white dwarf model[11, 12] with a fixed
$\log(g)_{\mathrm{WD}}=7.90\,\rm(cgs)$ based on an initial light curve model
yields an estimated effective temperature of $T_{\mathrm{WD}}=12600\pm
500\,\rm K$.

Extended Data Figure 1: a) A fit to the donor star’s SED. The blue vertical
lines represent the measured donor star SED as inferred from the apparent
magnitude during the primary eclipse. Because of the high SNR of the HiPERCAM
data, these uncertainties are very small, and we have added in a systematic
five percent uncertainty associated with contribution from the accretion disk.
The red diamonds represent the filtered averaged apparent magnitudes of the
best fit synthetic spectrum, and this synthetic spectrum is plotted in black.
b) A fit to the accreting white dwarf SED. The vertical blue lines represent
the SED of the white dwarf, as measured by modelling the amount of flux lost
during the primary eclipse, when it is occulted by the donor. We added a five
percent model error to this SED in the optical to account for contribution
from accretion features which might cause the white dwarf to deviate from a
standard hydrogen rich DA model spectrum[11, 12]. The surface gravity of the
white dwarf in the model spectrum was fixed using the values derived from the
lightcurve analysis. The red diamonds illustrate the filter averaged apparent
magnitudes from the synthetic spectrum, which is plotted in black. In both
panels, we have omitted the HiPERCAM $z_{s}$-band measurement, which exhibits
an excess likely associated with the accretion disk. The model deviates from
the observed spectrum around the hydrogen absorption lines due to emission
lines from the accretion disk around the white dwarf.

###  0.8 Atmospheric analysis

We summed the phase-resolved LRIS spectra in the rest frame of the donor and
performed fits with atmospheric models in order to estimate a temperature
independently of the SED based analysis. We broaden the spectra by the
predicted rotational velocity in the donor and convolve the synthetic spectra
with the instrument resolution. We fix the relative contributions of the WD
and donor by fixing the radius of each object in our spectroscopic models, as
well as the overall distance to the system. Using hydrogen rich white dwarf DA
models in combination with the BT-NextGen (GNS93)[13] low mass main sequence
models, we construct a grid of composite spectra to carry out the fitting
procedure. Extended Data Figure 2 shows our best fit, as well as a fit based
on the temperatures derived from the SED. We consider both fits acceptable,
though the estimated temperatures differ by $\approx 1200\,\rm K$. We consider
the temperature derived from the SED to be more reliable, as it uses the Swift
UVOT measurement, which strongly constrains the white dwarf temperature (the
estimated spectroscopic temperature of $13780\pm 380\,\rm K$ significantly
overestimates the Swift UVOT luminosity).

Extended Data Figure 2: a) A best fit spectral model consisting of a low mass
main-sequence model and a hydrogen rich white dwarf model to the observed
spectrum of ZTF J1813+4251. The model, shown in red, has been convolved with
the resolution of the spectrograph and rotationally broadened, and fit to the
spectrum coadded in the rest frame of the donor, shown in black. b) A model
spectrum fixed to the parameters derived from the SED analysis. We consider
this model more reliable because it takes into account the Swift ultraviolet
flux measurement, which strongly constrains the temperature of the white
dwarf. The model, shown in blue, slightly underfits the blue part of the
spectrum, plotted in black, which may be due to the presence of an accretion
disk, or simply a systematic error introduced in the data reduction.

###  0.9 Radial velocity analysis

We measured radial velocities from the Keck LRIS spectra by fitting nine
narrow absorption lines with Voigt profiles. We fit the Fe I lines at
$4271.7602,4325.7616,4383.5447,\\\
4404.7501,4891.4921,4920.5028,4957.5965\,\rm\AA$, the Mg I lines at
$4702.9909\,\rm\AA$ and the Ca I line at $4226.73\,\rm\AA$. We forced all
lines to share the same Doppler shift and full-width at half-maximum for both
the Gaussian and Lorentzian components of the Voigt profile as common free
parameters, but allowed the amplitude of each line to vary individually. We
used a least squares minimizer to fit each individual spectrum for a Doppler
shift relative to the rest frame, and then fit a sinusoidal function to these
velocities to derive a semi-amplitude and systemic velocity $\gamma$, as seen
in panel a) of Figure 2, with uncertainties derived from the covariance matrix
of the fit.

We validate the semi-amplitude and systemic velocity with an independent
cross-correlation analysis. We try two cross-correlation spectral templates: a
template matched to the spectral type, and a template created from summing the
ZTF J1813+4251 spectra shifted to the rest frame of the donor. We cross-
correlate the individual spectra from 4900 to 5600 Å, a wavelength range
filled with narrow absorption lines and dominated by the F-type star. We find
statistically identical results, $K=463.2$ km s-1, consistent within the
1$\sigma$ uncertainty of the published value (Table 1).

###  0.10 Rotational broadening analysis

We coadded the phase-resolved ESI spectra in the rest frame of the donor using
the radial velocity semi-amplitude derived from the LRIS spectroscopy. Due to
the highly non-linear dispersion in the instrument, only the red end of each
echelle order exhibited signal above the read-noise floor, and the SNR of the
coadded spectra was too low to obtain a useful constraint on the rotational
broadening of the donor. As illustrated in Extended Data Figure 3, we fit four
lines, the Fe I lines at $4920.5028,4957.5965\,\rm\AA$ and the Ca I lines at
$6122.22,6162.17\,\rm\AA$. These fits yielded measured rotationally broadened
velocities of $159,178,129~{}{\rm and}~{}184\,\rm km\,s^{-1}$, with a large
scatter due to the low SNR of the spectra, but with values roughly in
agreement of the predicted value of $145\,\rm km\,s^{-1}$ (they are on average
sligthly higher, but were not corrected for Dopper smearing, which should
inflate values by several tens of $\rm km\,s^{-1}$. Future higher SNR moderate
to high resolution spectroscopy could yield a more precise measurement, which
would serve as a test of the parameters reported in Table 1.

Extended Data Figure 3: Rotationally broadened atmospheric model fits to the
moderate resolution ESI spectra of ZTF J1813+4251. Due to the low SNR of the
spectra, we were unable to constrain the rotational broadening of the lines
with a precision better than $\sim 50\,\rm km\,s^{-1}$, but the measured
values are consistent with the predicted value of $145\,\rm km\,s^{-1}$.

###  0.11 Light curve analysis

The core of our analysis was based on the modelling of ZTF J1813+4251’s
HiPERCAM light curves using the LCURVE light curve modelling code[14]. LCURVE
simulates two stars in a Roche potential, determining the flux at each surface
element of the star, taking into account effects such as limb darkening and
gravity darkening, and computes the flux as seen from the perspective of an
observer at a given inclination and orbital phase. While this code is
numerical, below we outline a simplified analytic argument demonstrating why
the modelling code is able to robustly constrain the parameters in the system.
Extended Data Figure 4 illustrates this idealized example and provides all
relevant analytic expressions.

Extended Data Figure 4: A panel illustrating an idealized version of the basic
constraints we obtain by modelling the primary and secondary eclipses of ZTF
J1813+5251, which depend only on Roche geometry and Kepler’s laws. Using the
light curve, we are able to measure the orbital period, $P_{\mathrm{b}}$, the
in eclipse flux levels as a fraction of the out of eclipse flux, $F_{1}$ and
$F_{2}$, and the third and fourth contact phases of the eclipse, $\phi_{3}$
and $\phi_{4}$. By combining these five quantities we can determine from the
light curve, with the donor radial velocity semi-amplitude measured from the
spectra, $K_{2}$, and our knowledge that the donor is Roche-lobe-filling, we
are able to obtain a robust solution for the two component masses, $M_{1}$ and
$M_{2}$, the radii of the two components, $R_{1}$ and $R_{2}$, the surface
brightness ratio, $J$, the semi-major axis of the binary, $a$, and the orbital
inclination, $i$. The idealized expressions tying these constraints to the
observable quantities are listed in the right half of the figure. We would
like to emphasize that these are idealized expressions, and caution readers
that there are further important subtleties not discussed here (e.g. in a
Roche filling system, the relevant radii in the above equations $R_{1}$ and
$R_{2}$ are complicated to define because of the ellipsoidal deformation of
the components, however, lightcurve modelling codes such as LCURVE take these
effects into account–this means for example that the $R_{2}$ in expressions 3
and 4 is not exactly the same as the $R_{2}^{*}$ in expression 5, as in the
former case, $R_{2}$ is measured perpendicular to the line between the two
stars, whereas in the latter case, $R_{2}^{*}$ is a volume averaged radius).

Because ZTF J1813+4251 undergoes a total eclipse, in which the white dwarf
passes completely behind the donor, we are able to measure three parameters
from it: the fractional flux depth in eclipse with respect to the out of
eclipse flux, $F_{1}$, the third contact phase corresponding to the start of
egress relative to the mid-eclipse time, $\phi_{3}$, which measures the
duration of the flat bottomed portion of the eclipse when the white dwarf is
completely obscured, and the fourth contact phase, $\phi_{4}$, which measures
the end of egress with respect to the mid-eclipse time, and in combination
with the third contact phase, describes the duration of egress (one could
alternatively use the first and second contact phases, which measure the start
and end of ingress, but these are degenerate with the third and fourth contact
phases). In an idealized circular orbit involving two spherical objects, the
sum and differences of the scaled radii of the two components are constrained
by the following two expressions:

$\sin^{2}(i)\sin^{2}(2\pi\phi_{3})+\cos^{2}(i)=\left(\frac{R_{\mathrm{Donor}}}{a}-\frac{R_{\mathrm{WD}}}{a}\right)^{2}$
(1)
$\sin^{2}(i)\sin^{2}(2\pi\phi_{4})+\cos^{2}(i)=\left(\frac{R_{\mathrm{Donor}}}{a}+\frac{R_{\mathrm{WD}}}{a}\right)^{2}~{},$
(2)

where $\frac{R_{\mathrm{WD}}}{a}$, $\frac{R_{\mathrm{Donor}}}{a}$ are the
radii of the components, scaled to the semi-major axis. While the primary
eclipse gives three parameters to measure, its geometry depends on four
physical parameters: the orbital inclination, $i$, the scaled radii of the two
components and the surface brightness ratio of the two components, $J$. This
means that in order to constrain the parameters of the system from the primary
eclipse alone, it is necessary to invoke an additional constraint such as a
white dwarf mass-radius relation, or by constraining the surface brightness
ratio using temperature estimates of the components, which depend on model
atmospheres. However, we detect the secondary eclipse in the system, when the
white dwarf transits the donor. This provides one additional geometric
constraint, as we can measure the flux depth during the secondary eclipse as a
fraction of the out of eclipse flux, $F_{2}$ (the duration of ingress/egress
and the base of the secondary eclipse should be identical to that of the
primary, so the only additional parameter is the depth). Using the flux level
in the two eclipses, we can constrain the ratio of the radii of the two
components, using the following relation:

$\frac{R_{\mathrm{WD}}}{R_{\mathrm{Donor}}}=\sqrt{\frac{1-F_{2}}{F_{1}}}~{}.$
(3)

Additionally, the surface brightness ratio is also constrained by the relative
depths of the eclipses, and is given by:

$J=\frac{1-F_{1}}{1-F_{2}}~{}.$ (4)

The above expressions are only applicable for a highly idealized case of two
spherical objects transiting each other, and thus only approximate in the case
of ZTF J1813+4251, which hosts a highly distorted donor. The LCURVE modeling
simulates the fully distorted donor, and thus accounts for the complexities
accompanied by this.

The fact that the donor is filling its Roche lobe (which we know because the
system is mass transferring) allows us to place an additional constraint, as a
Roche-lobe-filling donor in a binary which undergoes a full eclipse must obey
a pair of non-linear parametric equations which map any given eclipse geometry
to a unique mass ratio, $q$, vs inclination relation[15]. Because we have
already solved for a unique inclination given the primary and secondary
eclipse geometry using the above expressions, this means we also obtain a
unique mass ratio as a result of the Roche-lobe-filling assumption. We have
also measured a velocity semi-amplitude of the donor, $K_{2}$, which is given
by the binary mass function:

$\frac{M_{\mathrm{WD}}\sin^{3}(i)}{(1+q)^{2}}=\frac{P_{\mathrm{b}}\,K_{\mathrm{Donor}}^{3}}{2\pi\mathrm{G}}~{},$
(5)

where G is the gravitational constant. Because we have already solved for a
unique mass ratio, inclination, and have determined the orbital period, this
expression yields a solution for $M_{\mathrm{WD}}$, which in combination with
the mass ratio, also yields a solution for $M_{\mathrm{Donor}}$. Finally,
since we have determined both component masses in the system, and know the
orbital period, we can use Kepler’s law to solve for the semi-major axis of
the system:

$\frac{a^{3}}{P_{\mathrm{b}}^{2}}=\frac{G(M_{R_{\mathrm{WD}}}+M_{R_{\mathrm{Donor}}})}{4\pi^{2}}~{},$
(6)

and thus, we use the semi-major axis in combination with the previously
determined scaled radii $\frac{R_{\mathrm{WD}}}{a}$,
$\frac{R_{\mathrm{Donor}}}{a}$ to compute the physical radii of the
components, $R_{\mathrm{WD}}$ and $R_{\mathrm{Donor}}$. Thus, using the
orbital period, four parameters measured from eclipse geometry, the Roche-
lobe-filling nature of the donor and its radial velocity semi-amplitude, we
are able to determine a unique solution for both component masses, radii, the
semi-major axis of the system, the surface brightness ratio of the components
and the orbital inclination.

When modelling the light curve, we used Gaussian priors constraining the
radial velocity semi-amplitude based on the spectroscopic analysis. We sampled
over the orbital inclination, $i$, the time of superior conjunction of the
white dwarf, $T_{0}$, the mass of the white dwarf, $M_{\mathrm{WD}}$, the mass
of the donor star, $M_{\mathrm{Donor}}$, the white dwarf to donor surface
brightness ratio in each filter, $J$, and the radius of the white dwarf,
$R_{\mathrm{WD}}$. We did not sample over the radius of the donor star, as we
instead fixed it to fill the Roche lobe. We fixed the limb darkening and
gravity darkening coefficients based on the temperature estimate of the white
dwarf and donor star derived from atmospheric fits[16, 17, 18, 19].

Extended Data Figure 5: A corner plot of some of the quantities derived from
our final model, illustrating clean convergence in the distributions for all
quantities. We would like to note that the radii are volume-averaged radii,
and that $R_{2}$ was not directly sampled (instead, the two masses+inclination
are sampled, and because the system is fully eclipsing $R_{2}$ is determined
uniquely by the mass ratio and inclination). We also sampled over the surface
brightness ratio J, but did not include it in these corner plots because we
have five different Js (one for each filter–we avoided using a common J
because doing this correctly requires atmospheric corrections in each
passband, and the solutions for the other free params are largely independent
of J). To ensure that using a different J for each filter was not influencing
the other free params, we modeled all 5 filters independently, and found that
they all converged to parameters in agreement with the combined fit.

The strong O’Connell effect and contribution of the disk to the lightcurve
makes it challenging to model the ellipsoidal variations in the system. We
subtracted a term at the orbital frequency out of the light curve to account
for this effect. However, we discovered there was a slight residual phase
offset between the time of superior conjunction of the white dwarf inferred
from these corrected ellipsoidal variations, and that indicated by the primary
eclipse.

Given that the system undergoes a full eclipse of the white dwarf and has a
Roche-lobe-filling donor, rather than attempt to model the full light curve
and account for the complexities in the ellipsoidal variations, we elected to
model only points around the primary and secondary eclipses (within $\pm
160\,\rm s$ of the mid-eclipse times), as the geometry of these contain all
the information needed to fully constrain the system, and provide a more
robust set of constraints than relying on ellipsoidal variations, which depend
on quantities like the gravity darkening coefficients of the exotic donor star
in each passband. The results of this modelling yielded the values reported in
Table 1, and Figure 6 illustrates the best-fit model overplotted with the
primary and secondary eclipses, as well as the full lightcurve, while Figure 1
illustrates the model fit to all the binned primary and secondary eclipses.
All of our models included a first order linear polynomial with both a slope
degree of freedom and arbitrary flux offset to minimize impacts of modulation
due to the O’Connell effect near the eclipses on our models. We wish to
emphasize that we did include ellipsoidal variations in our model, but also
ran a test model with a second order polynomial fit around each eclipse to
remove the sensitivity of our model to the ELVs (which can influence the
inferred mass ratio), and other effects such as contributions from the
accretion disk, effectively making the models only sensitive to the geometry
of the eclipses of the donor and white dwarf themselves (since those contain
all the information needed to constrain the system in a highly model-
independent manner). We found that the parameter estimates derived from the
models with the second order polynomial included were consistent with those
derived by simply using a first order polynomial, indicating that effects such
as the ELVs were not influencing the parameter estimates, and that the
eclipses indeed are the dominate features influencing the constraints in the
system.

Extended Data Figure 6: a) Our best fit model of the primary eclipse, with the
model shown as a solid black line, the binned data as red points, and the
model without an eclipse as the dashed black line. In this figure, a linear
polynomial component (of the functional form $y=a\times t$, where $t$ is the
time from mid-eclipse (we applied separate corrections for the primary and
secondary), has been subtracted out of both the model and the data for better
visualization. We constructed the model by simultaneously fitting data around
the primarily eclipse as well as data from the secondary eclipse (panel c). b)
A best fit of the eclipse-derived model constructed from the eclipse data
shown in panels a and c to the full HiPERCAM lightcurve. While the model
roughly reproduces the correct amplitude of ellipsoidal variations, it does
not to fully capture the structure seen in the full lightcurve. c) Our best
fit model to the secondary eclipse, with a linear correction subtracted out of
both the data and models. d) The residuals of the best fit of the eclipse-
derived model to the full dataset shown in panel b. As is readily apparent,
the strongest residuals occur out of eclipse, and likely arise due to a
combination of effects from the accretion disk, and an O’Connell effect
associated with the donor.

We used the Ultranest[8, 9, 10] package to perform the final parameter
estimation by combining the light curve modelling of the eclipses and the
donor’s radial velocity semi-amplitude constraint. The derived parameters are
reported in Table 1, and Extended Data Figure 5 illustrates the posterior
distributions of these parameters.

As a sanity check regarding our assumptions, we conducted an independent
modelling exercise, in which we did not assume a Roche-filling donor, but
instead invoked a white dwarf mass radius relation[20] for the accreting white
dwarf, and we found that the models converged to values consistent with those
derived from the models which assumed a Roche-filling donor. We expected this
model to converge to similar values, because the values of the white dwarf
mass and radius we found with our Roche-filling model agree remarkably well
with a white dwarf mass-radius relation.

###  0.12 Timing analysis

Because of the substantial contribution of magnetic braking to the predicted
orbital decay rate of the system, we undertook an effort to constrain the
orbital decay rate, $\dot{P}_{\mathrm{b}}$, as we predict this value should be
a factor of a few larger than what it would be due to purely general
relativistic orbital decay, $\dot{P}_{\mathrm{b}}^{\mathrm{GW}}=(-4.27\pm
0.24)\times 10^{-13}\rm\,s\,s^{-1}$. We use the HiPERCAM light curve in
combination with archival data from ZTF, the Palomar Transient Factory (PTF)
and the Catalina Real Time Transient Survey (CRTS) in order to constrain
$\dot{P}_{\mathrm{b}}$. We find a value consistent with zero, with an
uncertainty of $10^{-12}\rm\,s\,s^{-1}$, larger than the predicted value due
to pure general relativity, indicating that we need to continue monitoring the
system for the next several years to reach a degree of precision to test
whether the orbital evolution is consistent with pure general relativity.

###  0.13 Accretion disk analysis

In order to estimate the white dwarf radial velocity semi-amplitude, we tried
to obtain radial velocities from the centroid of the H$\alpha$ emission
component of the accretion disk by fitting a pair of Lorentzians with positive
amplitude and a symmetric splitting around a central radial velocity value.
Because of the low SNR of the individual spectra and significant contamination
from the H$\alpha$ absorption line of the donor star (and possibly from the
white dwarf itself), we were unable to extract a radial velocity for the
accretor in this manner. However, we used the masses estimated from our model
to construct a coadded spectrum in the rest frame of the accretor, as seen in
Figure 3, and then fit the accretion disk profile to determine the splitting
between the red and blue components in the disk. We find a splitting of
$17.33\pm 0.49\rm\,\AA$, which corresponds to a velocity of $810\pm
23\rm\,km\,s^{-1}$. Assuming a Keplerian disk, this corresponds to a radius of
$0.163\pm 0.009\,R_{\rm\odot}$ around the white dwarf. This is approximately
79 percent the radius of the the white dwarf’s Roche lobe, a regime in which
material should not be able to maintain a stable orbit[21]. We verify that we
find a similar result with the He I line at $6678\,\rm\AA$. It is possible
that the splitting in the disk lines is underestimated due to underlying
contamination from the donor and WD absorption lines, but our results suggest
that the matter producing the emission in this system exists near the limit of
where a stable accretion disk can form in the white dwarf’s Roche lobe.

We also investigated the object’s ZTF power spectrum to look for signatures of
a peak associated with a superhump excess, which could allow for an
independent measurement of the mass ratio in the system. However, the only
peaks we were able to identify in the power spectrum were those associated
with the orbital period and beat frequencies of the orbital period and the
sidereal day.

Curiously, the accretion disk does not manifest itself in a deep eclipse of
the donor, which is something we would expect for an optically thick accretion
disk. This may suggest that the instantaneous mass transfer rate of the system
is lower than the longer term average mass transfer rate predicted by our MESA
models, as some models would predict an optically thick accretion disk given
the estimated average mass transfer rate. This is in contrast to some other
sub-period minimum systems such as GALEX J194419.33+491257.0[22], which
exhibits signatures of a higher mass transfer rate, though we emphasize that
the donor in ZTF J1813+4251 is significantly more evolved than any other known
sub-period minimum CV, and thus its donor properties may also be influencing
its mass-transfer behavior in a unique manner.

###  0.14 Kinematic Analysis

We performed a kinematic analysis of ZTF J1813+4251’s orbit in the Galaxy
(plotted in Extended Data Figure 7) and find that it is consistent with
residing in the Galactic thick disk. This suggests that the progenitor stars
of the system formed several billion years ago and likely had a low
metallicity (though we emphasize that we do not believe the low metallicity
value inferred from spectroscopic fitting of the donor is reliable given its
highly exotic current state). We used the galpy[23] package to compute its
trajectory around the Milky Way over 6 Gyr, using the McMillan2017
potential[24].

Extended Data Figure 7: A set of panels illustrating the orbit of ZTF
J1813+4251 around the Milky Way. The system is consistent with residing in the
Galactic thick disk, orbiting between 5 and 8 kpc from the Galactic center,
within half a kpc of the Galactic disk in height.

###  0.15 MESA models

To investigate the formation history and future evolution of ZTF J1813+4251,
we calculated a grid of binary evolution models using MESA (version r12778)
[25, 26]. These models build on those calculated by [27], who used them to
interpret observations of evolved CVs with longer periods (2-6 hours). See
also the extensive grid of CV models, including ultra short period systems
with He donors, generated with MESA in Kalomeni et al. (2016)[28].

We initialized the models with a detached, circular-orbit binary containing a
$1.1\,M_{\odot}$ zero-age main sequence star (the donor) and a
$0.55\,M_{\odot}$ point mass representing the white dwarf. Following [29], we
use the Skye + CMS equation of state. Orbital angular momentum is removed via
gravitational wave radiation and magnetic braking. Roche lobe radii are
computed using the fit from [30], and mass transfer rates during Roche lobe
overflow are determined following the prescription from [31]. We model mass
transfer as being fully non-conservative, with the mass that is transferred
from the donor to the white dwarf eventually being lost from the vicinity of
the latter. The orbital separation evolves under the assumption that the mass
that is lost has the same specific angular momentum as the accretor. Although
mass transfer in real CVs is not instantaneously non-conservative, the mass
accreted by the white dwarf is expected to be ejected from the system by
classical nova explosions, which recur on a timescale that is short compared
to the timescale on which the orbit evolves.

Magnetic braking follows the prescription from [32], with $\gamma_{\rm mb}=3$.
Following [33], we assume that the magnetic braking torque is exponentially
suppressed when the donor’s convective envelope contains less than 2 percent
of the star’s mass. This is the default implementation of magnetic braking in
MESA as of version r15140; it effectively makes magnetic braking inefficient
when the donor’s effective temperature is $T_{\rm eff}\gtrsim 6000\,\rm K$. We
turn magnetic braking off at the period minimum.

We explored models with initial periods ranging from 0.5 to 3 day. At fixed
donor mass, longer initial periods lead to donors that are older, and hence
more evolved, at the onset of mass transfer. In Extended Data Table 1 and
Figure 8 we highlight models which have similar properties to ZTF J1813+4251
at a period of 51 minutes. These models begin mass transfer just as the donor
is completing its main sequence evolution. Models with shorter initial periods
begin mass transfer before the donor has undergone significant nuclear
evolution; these evolve into normal CVs. Models with significantly longer
initial periods begin mass transfer as the donor is evolving up the giant
branch and terminate their evolution as extremely low-mass white dwarfs with
long periods (in the case of stable Roche lobe overfilling) or short periods
(in the case of common envelope evolution).

In all the models shown in Figure 8 there is an initial period of thermal-
timescale mass transfer in which the mass transfer rate reaches $\sim
10^{-7}\,M_{\odot}\,\rm yr^{-1}$ and most of the donor’s outer envelope is
lost. By the time the models reach a period of 6 hours, all that remains of
the donor is a helium core and a puffy hydrogen envelope containing a few
$\times 10^{-2}\,M_{\odot}$ of hydrogen. This envelope is steadily removed as
the orbit shrinks, still driven primarily by magnetic braking. At periods of a
few hours, the donors in these models are very similar to the evolved CVs in
the [27] sample. They heat up as their envelope is removed, with the more
evolved models (those in which mass transfer begins later) being hotter at
fixed orbital period.

By the time the models reach $P_{\rm b}=51$ minutes, their surface is
predicted to be $\sim 75$ percent helium by mass. The surface nitrogen-to-
carbon ratio is predicted to be $\sim 2,000$ times the solar value, because
the present-day surface of the donor was previously inside the CNO-burning
core while the donor was on the main sequence. The predicted mass transfer
rates and accretor temperatures at 51 minutes differ significantly between the
three models we show: the more evolved models with higher $T_{\rm eff}$ have
less efficient magnetic braking and thus lower mass transfer rates.

We calculate the effective temperature of the accretor assuming it is set by
compressional heating, as described by Equation 2 of [34]. This makes the
accretor temperature a sensitive probe of the time-averaged accretion rate.
Both the donor and accretor temperatures of ZTF J1813+4251 are best-matched by
the model with $t_{\rm RLOF}/t_{\rm MS}=0.95$. Here $t_{\rm MS}$ is the main-
sequence lifetime the donor would have if it were an isolated star, and
$t_{\rm RLOF}$ is the age of the system at the onset of Roche lobe overflow.
In this model, magnetic braking is still the dominant mode of angular momentum
loss, leading to a predicted orbital inspiral that is a factor of $\sim 3$
larger than expected from gravitational radiation alone.

The models reach a minimum period between 13 and 22 minutes. After this, their
orbits begin to widen and they transition to primarily helium mass transfer
with fully-degenerate donors. The predicted inspiral time is $\approx$ 70 Myr
for the model with $t_{\rm RLOF}/t_{\rm MS}=0.95$. Observational constraints
on the period derivative will test these models, allowing us to directly
measure the relative importance of magnetic braking and gravitational waves in
removing angular momentum from the system.

Extended Data Figure 8: MESA binary evolution models (Extended Data Table 1). Red, black and blue lines show models that overflow their Roche lobes after 94, 95 and 97 percent of the donor’s main-sequence lifetime. The dashed vertical line shows 51 minutes.  Table 1: MESA models. All models have an initial donor mass of $1.1\,M_{\odot}$ and a white dwarf mass of $0.55\,M_{\odot}$. The three models correspond to the red, black and blue lines in Extended Data Figure 8.  $P_{\mathrm{b}}^{\rm init}$ | $t_{\rm RLOF}/t_{\rm MS}$ | $P_{\mathrm{b}}^{\rm RLOF}$ | $T_{\rm eff,\,Donor}$ | $M_{\rm Donor}$ | $P_{\mathrm{b}}^{\rm min}$ | $\log(\dot{M}/M_{\odot}\,\rm yr^{-1})$ | $\rm Time\,\,since\,\,RLOF$
---|---|---|---|---|---|---|---
$[\rm day]$ |  | $[\rm day]$ | $[\rm K]$ | $M_{\odot}$ | $[\rm min]$ |  | $[\rm Gyr]$
$2.63$ | $0.94$ | $0.60$ | $5068$ | $0.100$ | $22$ | $-9.31$ | $0.55$
$2.69$ | $0.95$ | $0.63$ | $6009$ | $0.107$ | $18$ | $-9.75$ | $0.63$
$2.75$ | $0.97$ | $0.66$ | $7156$ | $0.115$ | $13$ | $-10.46$ | $1.26$

##  1 Data Availability

Reduced HiPERCAM photometric data, LRIS spectroscopic data, and MESA tracks
resulting from the models are available at
https://github.com/kburdge/ZTFJ1813-4251.git. The ZTF used data are all in the
public domain. The proprietary period for the spectroscopic data will expire
at the start of 2022, at which point the raw spectroscopic images will also be
accessible via the Keck observatory archive.

##  2 Code Availability

Upon request, the first author will provide code (primarily in Python) used to
analyze the observations, create MESA models, and any data used to generate
figures (MATLAB was used to generate most of the figures). The LCURVE
modelling code can be found at: https://github.com/trmrsh/cpp-lcurve

## References

  * [1] Zechmeister, M. & Kürster, M.  The generalised Lomb-Scargle periodogram. A new formalism for the floating-mean and Keplerian periodograms.  _A &A_ 496, 577–584 (2009).
  * [2] Chambers, K. C. _et al._ The Pan-STARRS1 Surveys.  _arXiv e-prints_ arXiv:1612.05560 (2016).
  * [3] Burdge, K. B. _et al._ A Systematic Search of Zwicky Transient Facility Data for Ultracompact Binary LISA-detectable Gravitational-wave Sources.  _ApJ_ 905, 32 (2020).
  * [4] Perley, D. A.  Fully Automated Reduction of Longslit Spectroscopy with the Low Resolution Imaging Spectrometer at the Keck Observatory.  _PASP_ 131, 084503 (2019).
  * [5] Dhillon, V. S. _et al._ HiPERCAM: A quintuple-beam, high-speed optical imager on the 10.4-m Gran Telescopio Canarias.  _MNRAS_ (2021).
  * [6] Green, G. M., Schlafly, E., Zucker, C., Speagle, J. S. & Finkbeiner, D.  A 3D Dust Map Based on Gaia, Pan-STARRS 1, and 2MASS.  _ApJ_ 887, 93 (2019).
  * [7] Allard, F., Homeier, D. & Freytag, B.  Models of very-low-mass stars, brown dwarfs and exoplanets.  _Philosophical Transactions of the Royal Society of London Series A_ 370, 2765–2777 (2012).
  * [8] Buchner, J.  A statistical test for Nested Sampling algorithms.  _arXiv e-prints_ arXiv:1407.5459 (2014).
  * [9] Buchner, J.  Collaborative Nested Sampling: Big Data vs. complex physical models.  _arXiv e-prints_ arXiv:1707.04476 (2017).
  * [10] Buchner, J.  UltraNest - a robust, general purpose Bayesian inference engine.  _The Journal of Open Source Software_ 6, 3001 (2021).
  * [11] Tremblay, P. E. & Bergeron, P.  Spectroscopic Analysis of DA White Dwarfs: Stark Broadening of Hydrogen Lines Including Nonideal Effects.  _ApJ_ 696, 1755–1770 (2009).
  * [12] Tremblay, P. E., Bergeron, P. & Gianninas, A.  An Improved Spectroscopic Analysis of DA White Dwarfs from the Sloan Digital Sky Survey Data Release 4.  _ApJ_ 730, 128 (2011).
  * [13] Hauschildt, P. H., Allard, F. & Baron, E.  The NextGen Model Atmosphere Grid for 3000¡=Teff¡=10,000 K.  _ApJ_ 512, 377–385 (1999).
  * [14] Copperwheat, C. M. _et al._ Physical properties of IP Pegasi: an eclipsing dwarf nova with an unusually cool white dwarf.  _MNRAS_ 402, 1824–1840 (2010).
  * [15] Chanan, G. A., Middleditch, J. & Nelson, J. E.  The Geometry of the Eclipse of a Pointlike Star by a Roche-Lobe Companion.  _ApJ_ 208, 512–517 (1976).
  * [16] Claret, A. & Bloemen, S.  Gravity and limb-darkening coefficients for the Kepler, CoRoT, Spitzer, uvby, UBVRIJHK, and Sloan photometric systems.  _A &A_ 529, A75 (2011).
  * [17] Claret, A., Hauschildt, P. H. & Witte, S.  New limb-darkening coefficients for Phoenix/1d model atmospheres. II. Calculations for 5000 K $\leq$ Teff $\leq$ 10 000 K Kepler, CoRot, Spitzer, uvby, UBVRIJHK, Sloan, and 2MASS photometric systems.  _A &A_ 552, A16 (2013).
  * [18] Claret, A. _et al._ Gravity and limb-darkening coefficients for compact stars: DA, DB, and DBA eclipsing white dwarfs.  _A &A_ 634, A93 (2020).
  * [19] Claret, A. _et al._ Doppler beaming factors for white dwarfs, main sequence stars, and giant stars. Limb-darkening coefficients for 3D (DA and DB) white dwarf models.  _A &A_ 641, A157 (2020).
  * [20] Soares, E. d. A.  Constraining Effective Temperature, Mass and Radius of Hot White Dwarfs.  _arXiv e-prints_ arXiv:1701.02295 (2017).  1701.02295.
  * [21] Paczyński, B.  Evolutionary Processes in Close Binary Systems.  _ARA &A_ 9, 183 (1971).
  * [22] Kato, T. & Osaki, Y.  GALEX J194419.33+491257.0: An unusually active SU UMa-type dwarf nova with a very short orbital period in the Kepler data.  _PASJ_ 66, L5 (2014).  1403.0308.
  * [23] Bovy, J.  galpy: A python Library for Galactic Dynamics.  _ApJS_ 216, 29 (2015).
  * [24] McMillan, P. J.  The mass distribution and gravitational potential of the Milky Way.  _MNRAS_ 465, 76–94 (2017).
  * [25] Paxton, B. _et al._ Modules for Experiments in Stellar Astrophysics (MESA).  _ApJS_ 192, 3 (2011).
  * [26] Paxton, B. _et al._ Modules for Experiments in Stellar Astrophysics (MESA): Binaries, Pulsations, and Explosions.  _ApJS_ 220, 15 (2015).
  * [27] El-Badry, K., Rix, H.-W., Quataert, E., Kupfer, T. & Shen, K. J.  Birth of the ELMs: a ZTF survey for evolved cataclysmic variables turning into extremely low-mass white dwarfs.  _MNRAS_ (2021).
  * [28] Kalomeni, B. _et al._ Evolution of Cataclysmic Variables and Related Binaries Containing a White Dwarf.  _ApJ_ 833, 83 (2016).
  * [29] Wong, T. L. S. & Bildsten, L.  Mass Transfer and Stellar Evolution of the White Dwarfs in AM CVn Binaries.  _ApJ_ 923, 125 (2021).
  * [30] Eggleton, P. P.  Aproximations to the radii of Roche lobes.  _ApJ_ 268, 368–369 (1983).
  * [31] Kolb, U. & Ritter, H.  A comparative study of the evolution of a close binary using a standard and an improved technique for computing mass transfer.  _A &A_ 236, 385–392 (1990).
  * [32] Rappaport, S., Verbunt, F. & Joss, P. C.  A new technique for calculations of binary stellar evolution application to magnetic braking.  _ApJ_ 275, 713–731 (1983).
  * [33] Podsiadlowski, P., Rappaport, S. & Pfahl, E. D.  Evolutionary Sequences for Low- and Intermediate-Mass X-Ray Binaries.  _ApJ_ 565, 1107–1133 (2002).
  * [34] Townsley, D. M. & Gänsicke, B. T.  Cataclysmic Variable Primary Effective Temperatures: Constraints on Binary Angular Momentum Loss.  _ApJ_ 693, 1007–1021 (2009).