A tidally tilted sectoral dipole pulsation mode in the eclipsing binary TIC 63328020
S.A. Rappaport, D.W. Kurtz, G. Handler, D. Jones, L.A. Nelson, H. Saio, J. Fuller, D.L. Holdsworth, A. Vanderburg, J.Žák, M. Skarka, J. Aiken, P.F.L. Maxted, D.J. Stevens, D.L. Feliz, F. Kahraman Aliçavu?
MMon. Not. R. Astron. Soc. , 000–000 (2020) Printed 4 February 2021 (MN L A TEX style file v2.2)
A tidally tilted sectoral dipole pulsation mode in the eclipsing binaryTIC 63328020
S. A. Rappaport (cid:63) , D. W. Kurtz , , G. Handler , D. Jones , , L. A. Nelson , H. Saio ,J. Fuller , D. L. Holdsworth , A. Vanderburg , J. ˇZ´ak , , , M. Skarka , , J. Aiken ,P. F. L. Maxted , D. J. Stevens , D. L. Feliz , , and F. Kahraman Alic¸avus¸ , Department of Physics, and Kavli Institute for Astrophysics and Space Research, M.I.T., Cambridge, MA 02139, USA Centre for Space Research, Physics Department, North West University, Mahikeng 2745, South Africa Jeremiah Horrocks Institute, University of Central Lancashire, Preston PR1 2HE, UK Nicolaus Copernicus Astronomical Center, Polish Academy of Sciences, ul. Bartycka 18, 00-716, Warszawa, Poland Instituto de Astrof´ısica de Canarias, E-38205 La Laguna, Tenerife, Spain Departamento de Astrof´ısica, Universidad de La Laguna, E-38206 La Laguna, Tenerife, Spain Department of Physics and Astronomy, Bishop’s University, 2600 College St., Sherbrooke, QC J1M 1Z7 Astronomical Institute, Graduate School of Science, Tohoku University, Sendai 980-8578, Japan Division of Physics, Mathematics and Astronomy, California Institute of Technology, Pasadena, CA 91125, USA Department of Astronomy, The University of Texas at Austin, 2515 Speedway, Stop C1400, Austin, TX 78712, USA Department of Theoretical Physics and Astrophysics, Masaryk Univesity, Kotl´aˇrsk´a 2, 60200 Brno, Czech Republic ESO, Karl-Schwarzschild-str. 2, D-85748 Garching, Germany Astronomical Institute, Czech Academy of Sciences, Friˇcova 298, 25165, Ondˇrejov, Czech Republic Astrophysics Group, Keele University, Staffordshire, ST5 5BG, UK Department of Astronomy & Astrophysics and Center for Exoplanets and Habitable Worlds, PSU, 525 Davey Lab, University Park, PA 16802, USA Department of Physics and Astronomy, Vanderbilt University, Nashville, TN 37235, USA Department of Physics, Fisk University, 1000 17th Avenue North, Nashville, TN 37208, USA C¸ anakkale Onsekiz Mart University, Faculty of Sciences and Arts, Physics Department, 17100, C¸ anakkale, Turkey C¸ anakkale Onsekiz Mart University, Astrophysics Research Center and Ulup˜onar Observatory, TR-17100, C¸ anakkale, Turkey
Accepted XXX. Received YYY; in original form ZZZ
ABSTRACT
We report the discovery of the third tidally tilted pulsator, TIC 63328020. Observations withthe
TESS satellite reveal binary eclipses with an orbital period of 1.1057 d, and δ Scuti-typepulsations with a mode frequency of 21.09533 d − . This pulsation exhibits a septuplet oforbital sidelobes as well as a harmonic quintuplet. Using the oblique pulsator model, theprimary oscillation is identified as a sectoral dipole mode with l = 1 , | m | = 1 . We find thepulsating star to have M (cid:39) . (cid:12) , R (cid:39) (cid:12) , and T eff , (cid:39) K, while the secondaryhas M (cid:39) . (cid:12) , R (cid:39) (cid:12) , and T eff , (cid:39) K. Both stars appear to be close tofilling their respective Roche lobes. The properties of this binary as well as the tidally tiltedpulsations differ from the previous two tidally tilted pulsators, HD74423 and CO Cam, inimportant ways. We also study the prior history of this system with binary evolution modelsand conclude that extensive mass transfer has occurred from the current secondary to theprimary.
Key words: stars: oscillations – stars: variables – stars: individual (TIC 63328020)
The single-sided or tidally-tilted pulsators are a newly recognizedtype of pulsating star in close binary systems. In these stars the tidaldistortion caused by the companion aligns the pulsation axis of theoscillating star with the tidal axis. This has two important conse- (cid:63)
E-mail: [email protected] quences. First, the pulsation axis corotates with the orbit and there-fore the stellar pulsation modes are seen at varying aspect, leadingto amplitude and phase variations with the orbital phase. Second,the tidal distortion of the pulsating star causes an intrinsically un-even distribution of pulsation amplitude over the stellar surface.These two facts can be used to our astrophysical advantage.Viewing the stellar pulsation over the full orbital cycle allows us toconstrain the orbital inclination, i , and the obliquity of the pulsa- © 2020 RAS a r X i v : . [ a s t r o - ph . S R ] F e b S. Rappaport et al. tion axis, β , to the orbital axis. The required mathematical frame-work, the oblique pulsator model, has been developed over the pastfour decades, starting with Kurtz (1982). Even though it was de-veloped for the rapidly oscillating Ap (roAp) stars, whose pulsa-tion axes are tilted with respect to their pulsation axes due to thestellar magnetic fields, it is readily applicable to tidally tilted pul-sators. Since the pulsation axes of those stars are located in the or-bital plane, the analyses of the binary-induced variability and of thetidally tilted pulsation mutually constrain each other. The predom-inant shape of the distorted oscillations, and thus the “pulsationalquantum numbers” – the spherical degree l and the azimuthal order m – can be determined from the variation of the pulsation ampli-tude and phase over the orbit. In reality, the tidal distortion inducescoupling between modes of different spherical degree l , modify-ing the perturbed flux distribution at the surface, and sometimes“tidally trapping” the oscillations on one side of the star. The theo-retical groundwork for this type of analysis was laid by Fuller et al.(2020).Two single-sided pulsators have so far been reported in the lit-erature; both are ellipsoidal variables and each contains at least one δ Scuti pulsator. Ellipsoidal variables are close binary stars withtidally distorted components. They exhibit light variations as theprojected stellar surface area and surface gravity vary towards thedirection to a distant observer (e.g., Morris 1985) over the orbitalcycle. The δ Scuti stars, on the other hand, are a common groupof short-period pulsators located at the intersection of the classicalinstability strip with the Main Sequence (Breger 1979, 2000). Nat-urally, some δ Scuti stars are also located in binary systems (e.g.,Liakos & Niarchos 2017), with a wide range of phenomenologyoccurring, for instance pulsators in eclipsing binaries (KahramanAlicavus et al. 2017) or the so-called ‘heartbeat’ stars with tidallyexcited stellar oscillations in binary systems with eccentric orbits(e.g., Welsh et al. 2011). However, in none of these systems wasevidence for tidal effects on the pulsation axes, such as those oc-curring in the single-sided pulsators, reported.The first such discovery was HD 74423 (Handler et al. 2020),which contains two chemically peculiar stars of the λ Bootis typein a 1.58-d orbit that are close to filling their Roche lobes. Al-though the two components are almost identical, only one of themshows δ Scuti pulsation. Rather unusual for this type of pulsatingstar, there is only a single mode of oscillation, and it is not clearwhich of the two components is the pulsator. Shortly afterwards,Kurtz et al. (2020) reported the discovery of a second such system,CO Cam. The properties of this binary are different from those ofHD 74423. Its orbital period is somewhat shorter (1.27 d), the sec-ondary component is spectroscopically undetected, hence consid-erably less luminous than the pulsating primary, which is far fromfilling its Roche lobe. The pulsating star in the CO Cam systemis also chemically peculiar, but it is a marginal metallic-lined A-Fstar (often denoted with the spectral classification “Am:”), and itpulsates in at least four tidally distorted modes.Both HD 74423 and CO Cam were designated as “single-sidedpulsators” because they show enhanced pulsation amplitude on theL side of the star facing the secondary, as explained by Fuller et al.(2020).. In the present paper, we report the discovery of the thirdsingle-sided pulsator, TIC 63328020, which is different from thetwo systems studied earlier. With the discovery of a sectoral pulsa-tion mode in TIC 63328020, we now refer to these stars generallyas ‘tidally tilted pulsators’, where the more specific name, ‘single-sided pulsators’, can still be used for those stars that have stronglyenhanced pulsation amplitude on one side of the star.In Section 2 we discuss how this object was first noticed in Transiting Exoplanet Survey Satellite ( TESS ) data and we presenta detailed analysis of the pulsations. In particular, we show that thepulsations are strongly modulated in amplitude and phase aroundthe orbit, with the peak pulsation amplitudes coinciding in timewith the maxima of the ellipsoidal light variations (‘ELVs’), hencein quadrature with the eclipses. In Section 3.1 we present the re-sults of a study of archival data for the spectral energy distribution(‘SED’) as well as the long-term eclipse timing variations (‘ETVs’)for the system. Our radial velocity (RV) data for the system arepresented and analysed in Section 3.2.1, while Section 3.3 utilizesthe RV and SED data to analyse the system properties. As a com-plementary analysis of the system parameters, in Section 3.4 wecombine the RV data and the
TESS light curve via the phoebe2 code to derive the system parameters. Finally, in Section 4 we use aseries of
MESA binary evolution grids to understand the formationand evolutionary history of TIC 63328020. We conclude that therewas almost certainly a prior history of mass transfer in the systemand that the roles of the primary and secondary stars have reversed.
TIC 63328020 = NSVS 5856840 was reported as an eclipsing bi-nary by Hoffman et al. (2008), the only literature reference to thisstar. It has an apparent visual magnitude (cid:39) , and no spectral typewas given. The archival properties of TIC 6338020 are summarisedin Section 3.1.1 below.In addition to the eclipses, δ Scuti-type pulsations were dis-covered by one of us (DWK) during a visual inspection of
TESS
Sector 15 light curves. During its main two-year mission, − , TESS observed almost the whole sky in a search for transitingextrasolar planets around bright stars ( < I c < ) in a wide red-bandpass filter. The measurements were taken in partly overlapping × ◦ sectors around the ecliptic poles that were observed fortwo 13.6-d satellite orbital periods each (Ricker et al. 2015).TIC 63328020 was observed by TESS in Sectors 15 and 16in 2-min cadence. We used the pre-search data conditioned sim-ple aperture photometry (PDCSAP) data downloaded from theMikulski Archive for Space Telescopes (MAST) . The data havea time span of 51.95 d with a centre point in time of t =BJD 2458737 . , and comprise 34657 data points after someoutliers were removed with inspection by eye.Fig. 1 shows the TESS photometry for TIC 63328020. The toppanel displays the full Sectors 15 and 16 light curve, where theeclipses and ellipsoidal light variations are obvious. The bottompanel shows a short segment of the light curve where more detailscan be seen. Importantly, a careful look reveals that the pulsationalvariations are largest on the ELV humps where the star is brightest.That, of course, is at orbital quadrature, so this is distinctly differentfrom the first two single-sided pulsators, HD 74423 and CO Cam.The frequency analysis below bears out this first impression. First,we look at the light curve, which shows the ellipsoidal orbital vari- http://archive.stsci.edu/tess/all products.html This t was used to begin the analysis, but later changed to the time ofpulsation maximum to test the oblique pulsator model. For the assessmentof phase errors with nonlinear least-squares fitting, it is important that the t chosen is near to the centre of the data set. Since frequency and phase aredegenerately coupled in the fitting of sinusoids, when t is not the centreof the data set, small changes in frequency result in very large changes inphase, since phase is referenced from t .© 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 ations clearly, and the amplitude modulation of the pulsations oncareful inspection. We analysed the data using the frequency analysis package
PERIOD
04 (Lenz & Breger 2005) , a Discrete Fourier Transformprogram (Kurtz 1985) to produce amplitude spectra, and a com-bination of linear and nonlinear least-squares fitting to optimisefrequency, amplitude and phase. The derived orbital frequency is ν orb = 0 . ± . d − ( P orb = 1 . ± . d), where variance from the pulsations and from somelow frequency artefacts have been filtered for a better estimate ofthe uncertainty in the frequency. A 50-harmonic fit by least-squareswas done to see how this frequency fits the data, and to show thatthe pulsation frequencies are not orbital harmonics. In practice, theharmonics at frequencies higher than × ν orb do not have statisti-cally significant amplitudes. In particular, the 20th orbital harmonichas S/N = 7 in amplitude. All higher harmonics have S/N < inamplitude.After pre-whitening the data by the 50-harmonic fit, some lowfrequency artefacts from the data reduction were removed with ahigh-pass filter, which was a simple consecutive pre-whitening oflow frequency peaks extracted by Fourier analysis until the noiselevel was reached in the frequency range − d − . This was doneto study the pulsations with white noise for the purpose of estimat-ing the uncertainties. TIC 63328020 pulsates principally in a single mode at a frequencyof ν = 21 . ± . d − , typical of δ Sct stars. Becausethis oblique nonradial pulsation mode is observed from changingaspect with the orbit of the star and its synchronous rotation, am-plitude and phase modulation of the pulsation generate a frequencyseptuplet . There is also a harmonic at ν that generates a quin-tuplet. These are typical of oblique pulsators. There are also twolow-amplitude frequencies at 10.502 d − and 11.406 d − that areseparated by the orbital frequency. It is likely that one of these is amode frequency and the other is part of a frequency multiplet fromoblique pulsation where the signal-to-noise is too low to detect theother multiplet components. Whichever of these two frequenciesis the mode frequency, it is close to, but is not, a sub-harmonicof the principal mode frequency. These peaks at .
502 d − and .
406 d − have amplitude signal-to-noise ratios of 5 and 6, re-spectively, which are too low for further discussion of these fre-quencies here.The top panel in Fig. 2 shows the amplitude spectrum for thehigh-pass filtered data. An inspection shows that there is a septu-plet centred on ν = 21 . d − and a quintuplet centred on ν = 42 . d − . The detected frequencies are listed in Ta-ble 1. Because the orbital inclination is close to ◦ , and the tidalpulsation axis is inclined ◦ to that, the pulsation shows ampli-tude maxima twice per orbit, as will be seen in the next subsec-tion. That then generates two principal peaks in the amplitude spec- Because this oblique nonradial pulsation mode is a distorted dipole modeobserved from changing aspect with the orbit of the star and its synchronousrotation, amplitude and phase modulation of the pulsation generate a fre-quency septuplet (a pure dipole mode would generate a frequency triplet).
Table 1.
A least squares fit of the two low frequencies, the frequency sep-tuplet for ν and the frequency quintuplet for ν . The zero point for thephases, t = BJD 2458737 . , has been chosen to be a time when thetwo first orbital sidelobes of ν have equal phase.frequency amplitude phased − mmag radians ± . ν low − ν orb − . ± . ν low . ± . ν − ν orb . ± . ν − ν orb . ± . ν − ν orb . ± . ν . ± . ν + ν orb . ± . ν + 2 ν orb − . ± . ν + 3 ν orb . ± . ν − ν orb − . ± . ν − ν orb . ± . ν . ± . ν + ν orb − . ± . ν + 2 ν orb − . ± . trum at ν − ν orb and ν + ν orb , along with the other sidelobes.We determined the frequencies in the multiplets by a combinationof linear and nonlinear least-squares fitting, and determined thatall separations in the multiplets are equal to the orbital frequency ν orb = 0 . ± . d − determined in the last sectionwithin 1.5 σ . We took the average of those two highest amplitudesidelobe frequencies to determine the value of ν , and generatedthe frequency multiplets from that. Thus we only give frequencyuncertainties on these two highest amplitude frequencies, i.e., thosefor the first orbital sidelobes ν − ν orb and ν + ν orb , which are ± . d − and ± . d − , respectively.We then forced a frequency septuplet for ν and a quintu-plet for ν , all split by exactly the orbital frequency, ν orb =0 . ± . d − . The reason for choosing this ex-act splitting is that for oblique pulsators it is instructive to exam-ine the pulsation phases, and those are inextricably coupled to thefrequencies, as can be seen by examining the function we fitted,( cos 2 πf ( t − t o ) + φ ). Pre-whitening by that solution (and includ-ing in the fit the two low-frequencies 10.502 d − and 11.406 d − )leads to no variance in the data above noise, as can be seen in thebottom panel of Fig. 2. This shows that exact splitting by the orbitalfrequency about ν and 2 ν fits the data, and that there are no otherdetectable pulsations. This star has one principal oblique pulsatingmode.Finally, we chose a time zero point, t , such that the pulsa-tion phases were equal for the two dominant peaks seen in the toppanel of Fig. 2. This effectively choses a time of pulsation ampli-tude maximum with the orbit and rotation of the star, thus providesthe orbital phase for the time of pulsation maximum. Table 1 showsa least-squares fit of the determined frequencies. While the frequencies, amplitudes and phases determined byFourier analysis and least-squares fitting in the last subsection con-tain the information to study the oblique pulsation, it is instructiveand easier to see how this pulsation varies with orbital aspect by © 2020 RAS, MNRAS , 000–000
S. Rappaport et al.
Figure 1.
Top: The full Sectors 15-16 light curve of TIC 63328020 showing the orbital variations. Bottom: A section of the light curve where close inspectionshows the pulsations, particularly at orbital quadrature. The zero point of the ordinate scale is the mean.
Figure 2.
Top: The amplitude spectrum of the high pass data. There is a frequency septuplet centred on ν = 21 . c/d, and a quintuplet centred on ν .The reason there are two high amplitude first sidelobes about an almost zero amplitude mode pulsation frequency (marked by the vertical dotted red line) isthat we are seeing this pulsation inclined by ◦ to the orbital axis, which itself is near to i = 90 ◦ , thus giving two pulsation amplitude maxima per orbit.© 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 plotting the pulsation amplitude and phase as a function of orbital(rotational) phase. To do this we fitted the pulsation frequency, ν ,and its harmonic, ν , to chunks of the data that are P orb / induration. It is immediately obvious from doing this that pulsationamplitude maximum occurs in quadrature to the orbital eclipses.This is the signature of an oblique pulsation in a dipole sectoralmode.To show this, we have set the time zero point of anorbital period prior to pulsation maximum, as determined fromthe fitting in the last subsection. That zero point is t =BJD 2458736 . . We emphasise that this time has been chosenwith reference to the pulsation amplitude maximum, hence is inde-pendent of the determination of the time of primary eclipse fromthe study of the orbital variations.Fig. 3 shows the results. The mode is a sectoral dipole modewith (cid:96) = 1 , | m | = 1 . Amplitude maximum coincides with or-bital quadrature, and the phase reverses by π rad at the times of theeclipses when the line of sight aligns with the tidal axis. This isnew and currently unique among tidally tilted pulsators. Becausethe mode is sectoral, it has a symmetry with respect to the tidal dis-tortion such that the star is not strongly “single-sided”. The thirdpanel of Fig. 3 shows that the harmonic distortion of the mode isstrongest during secondary eclipse when the L side of the pulsat-ing star is closest to the observer, thus the L and L sides of thepulsator do differ and the star is mildly a “single-sided pulsator”. The standard simple relation for a toy model pulsator of P (cid:113) ρρ (cid:12) = Q can be rewritten in terms of observables as log Q = log P + log g + M bol + log T eff − . ,where Q is a pulsation “constant” that can be compared with mod-els, P is the pulsation period in days and log g is in cgs units. Tak-ing T eff (cid:39) K, log g (cid:39) . (Sect. 3.3), and M bol (cid:39) M V (cid:39) . from the Gaia parallax and V magnitude then gives for ν avalue of Q = 0.016, indicative of radial overtone around n = 3 − .The same calculation for the low frequency peak at 11.406 d − gives Q = 0 . , suggesting a second radial overtone mode. We have collected a set of archival magnitudes for TIC 63328020and report these values in Table 2. Additionally, we list the Gaiainformation about this object in Table 2. There is a fainter ( G =19 . ) neighbour star some 2.84 (cid:48)(cid:48) away, but there is insufficient Gaiainformation to tell whether that star is physically associated withTIC 63328020. The spectral energy distribution (SED) points for this object areplotted in flux units in Fig. 4, and many of them are reported as The bolometric correction near F0 is zero (see Table 2 of Morton &Adams 1968).
Figure 3. ν : Top and second panels: The pulsation amplitude and phasevariation as a function of orbital phase taking ν to be the pulsation fre-quency. The zero point in time, t = BJD 2458736 . , has been setto be of an orbital period before the time of pulsation maximum, whichwas determined by choosing the time when the phases of the first orbitalsidelobes are equal. Third panel: The pulsation amplitude as a function oforbital phase for the second harmonic 2 ν . The phase diagram for this fre-quency is not shown, since it is uninformative because the low amplituderesults in high scatter in the phase determinations. The second harmonicpeaks when the L point is closest to the observer. Bottom: the orbital lightvariations as a function of orbital phase for comparison. For this panel thedata were binned by a factor of 10. It can be seen that orbital light maximumcoincides with pulsation maximum, and that the pulsation amplitude is verysimilar at orbital phases 0.0 and 0.5. This is currently unique. magnitudes in Table 2. The SED points have all been corrected forinterstellar extinction at a level of E ( g − r ) (cid:39) . ± . (Greenet al. 2018, 2019), which we take to mean A V (cid:39) . We then usethe wavelength dependence of extinction given by Cardelli, Clay-ton & Mathis (1989), in particular, the fitting formulae in theirEqn. (2) and (3) for (cid:104) A ( λ ) /A ( V ) (cid:105) . Also shown on the figure arefitted curves to the SED based on the contributions from both starsin the binary. These will be discussed in Section 3.3.2. In addition to the new
TESS photometry on TIC 63328020, we haveutilized archival photometric data from ASAS-SN, WASP, KELT,and DASCH (for references see Table 3) to establish the long-termorbital ephemeris for this binary. The time intervals for the variousdata sets, and references to the archival data are given in Table 3.The DASCH data cover more than a century, but only about 1100 (http://argonaut.skymaps.info/query)© 2020 RAS, MNRAS , 000–000 S. Rappaport et al.
Table 2.
Properties of the TIC 63328020 SystemParameter ValueRA (J2000) (h m s) 21:20:14.41Dec (J2000) ( ◦ (cid:48) (cid:48)(cid:48) ) 51:23:41.04 T a . ± . G b . ± . G BP b . ± . G RP b . ± . B a . ± . V a . ± . J c . ± . H c . ± . K c . ± . W1 d . ± . W2 d . ± . W3 d . ± . W4 d > . R ( R (cid:12) ) b . +0 . − . L ( L (cid:12) ) b . ± . Orbital Period (d) e . K (km s − ) e . ± . γ (km s − ) e − . ± . Distance (pc) b ± µ α (mas yr − ) b +2 . ± . µ δ (mas yr − ) b +1 . ± . Notes. (a) ExoFOP (exofop.ipac.caltech.edu/tess/index.php). (b) Gaia DR2(Lindegren et al. 2018; Luri et al. 2018; Gaia Collaboration et al. 2018). (c)2MASS catalog (Skrutskie et al. 2006). (d) WISE point source catalog(Cutri & et al. 2013). (e) This work; see Table 4 for details on the orbitalperiod and Table 5 and Fig. 6 for the RV results.
Figure 4.
SED data and model for TIC 63328020 spanning the blue to 20micron region for the non-coeval case described in Section 3.3.2. The con-tinuous green, blue and red curves represent the contributions of the sec-ondary, the primary, and the total system flux, respectively. scanned plates for this object were available; we divided these upinto three roughly 30-yr long intervals which are denoted “1”, “2”and “3”. The same was done for the KELT data which spannedthree observing seasons and was divided into three ∼ Figure 5.
Long-term eclipse timing variations for TIC 63328020. The datasets used are marked with labels and include
TESS , ASAS-SN, WASP,KELT, and DASCH (see Table 3 for references). The latter two data setshave each been divided into three subsets of 30 years and ∼
10 months, re-spectively. The time intervals covered by each data set and the referencesare given in Table 3. The error bars plotted here and used in the long-term fitare taken from Eqn. (2) with σ J found from the fit to be 0.0057 d (see textand Eqn. 1). We did this so as not to allow the tiny error bars of the modern-epoch points to totally dominate the determination of the long-term periodor its derivative. The red curve is the best fit to the linear plus quadraticterms. See Table 4 for results. The multiple curves in faint grey are 100random draws from the posteriors of the MCMC fit. Table 3.
Eclipse Timing Variation Data for TIC 63328020Source Start Date a End Data a ETV b TESS c +0 . ± . ASAS-SN d +0 . ± . KELT 3 e +0 . ± . KELT 2 e +0 . ± . KELT 1 e +0 . ± . WASP f − . ± . DASCH 3 g − . ± . DASCH 2 g − . ± . DASCH 1 g +0 . ± . Notes. (a) HJD-2 400 000. (b) The fold period is 1.105 769 8 days with anepoch of HJD 2 444 000.2795. The ETV value is expressed in days. (c)This work. (d) Shappee et al. (2014) and Kochanek et al. (2017). (e)Pepper et al. (2007) and Pepper et al. (2012). (f) Collier Cameron et al.(2006); Pollacco et al. (2006). (g) Grindlay et al. (2009). term orbital period is well defined to a about a part per million with P orb = 1 .
105 769 8(3) d, it is also apparent that there are signifi-cant non-linear ETVs present. At the moment, there is insufficientinformation to quantify whether these are due to orbital motion in-duced by a distant companion or some other effect causing jitterin P orb (e.g., Applegate 1992). To get a handle on the long-termtrend in the orbital period, we modelled the measured ETVs as aquadratic function. We fitted simultaneously for a “jitter” term thatwe added in quadrature to the measured statistical uncertainties.The jitter term models an independent noise term in our ETV mea-surements that is not captured by our formal uncertainties (perhapsa systematic uncertainty due to the way we measured the ETVs).Our log likelihood function was: log L = − (cid:88) i (cid:20) ( y i − m i ) σ i + log σ i (cid:21) (1) © 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 Table 4.
Orbital Period Determinations for TIC 633280200Parameter Value Uncertainty
TESS
Only P orb a [days] 1.105 749 0.000 001 P orb b [days] 1.105 9 0.000 2 P orb c [days] 1.105 751 0.000 001 P orb d [days] 1.105 754 0.000 001 P orb / ˙ P orb e [years] − P orb f [days] 1.105 769 8 0.000 000 3 P orb g [days] 1.105 770 3 0.000 000 3 P orb / ˙ P orb h [years] +12 . × . × Jitter Noise, σ J i [days] 0.0057 +0 . − . Notes. (a) Based on the frequency analysis. (b) Spacing between the ν pulsation sextuplet. The corresponding epoch time is JD 2 458 736.8729.(c) Eclipse timing analysis (‘ETV’) assuming no period changes. Thereference fold epoch is JD 2 458 736.8720. (d) ETV analysis allowing for ˙ P orb . (e) P orb / ˙ P orb from the ETV analysis. (f) From the long-termphotometric ETV analysis assuming P orb is a constant. The fold epochtime is JD 2 444 000.2795. (g) From the long-term photometric ETVanalysis allowing for ˙ P orb . (h) P orb / ˙ P orb from the long-termphotometric ETV analysis (i) See Eqn. 2 for definition. where y i are the measured ETVs, m i are evaluations of thequadratic model, and σ i = (cid:113) σ ,i + σ J , (2)where σ ,i are the formal uncertainties on the ETVs and σ J is theextra uncertainty added in quadrature. The four free parameters,constant, linear, quadratic and σ J were found via an MCMC fit(see, e.g., Ford 2005).In Table 4 we summarize all the information that we haveabout the orbital period, and its derivative, derived in several dif-ferent ways from the various available data sets. TIC 63328020 was observed with the Intermediate DispersionSpectrograph (IDS) on the 2.5-m Isaac Newton Telescope (INT)between 28 November and 1 December 2019. The blue-sensitiveEEV10 detector was used along with a (cid:48)(cid:48) wide slit and the R1200Bgrating centred at 4000 ˚A for an unvignetted spectral coverageof ∼ STARLINK routines.
We extracted the K -velocity of the primary star via two differ-ent approaches. In the first we simply fit the deep Ca II K line (at3934 ˚A) and thereby estimated RVs with the corresponding uncer-tainties. In the second, we did a cross-correlation analysis. For the Table 5.
Radial Velocity Data for TIC 63328020Epoch (BJD) RV (Ca II K) a RV (CCF) b RV (average) c ( − − km s − km s − . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . . ± . Epoch of φ K d γ d cycles after TESS km s − km s − . ± . − . ± . . ± . Notes. (a) Based on measurements of the Ca II K line only. (b) Based on across correlation analysis of the spectrum over the region of 3700 to4700 ˚A. (c) Formal statistically weighted average of the two RV analyses.(d) Based on the average of the RVs, and a simple sinusoidal fit that doesnot take into account the Rossiter-McLaughlin effect (but see Section 3.4for a more complete analysis).
Figure 6.
Radial velocity curve for the primary in TIC 63328020. Becausethe orbital period is relatively close to a day, the orbital phase sampled overthe four nights of observations did not cover the lower half of the RV curve.The red curve is a model fit using
Phoebe2 , yielding ± km/s (seeSect. 3.4). latter, we removed the spectrograph blaze function from the spectraby breaking the spectra into 4-nm wide bins, identifying the highest10 per cent of flux measurements within each bin (most of whichare in the continuum), fitting a basis spline to the continuum points,and dividing the spectrum by the best-fit spline. We then measuredradial velocities by cross-correlating each blaze-corrected spectrumwith the highest signal-to-noise observation. The results are verysimilar to those obtained from the Ca II K line alone, but the un-certainties for the cross-correlation result are empirically somewhatsmaller. We list both sets of RVs in Table 5.We do not see any lines from the companion star. We estimatethat the luminosity of the secondary is (cid:46)
10 per cent that of theprimary.The RV data were taken over four consecutive nights, andsince the orbital period is 1.1 d, there was insufficient time for theorbital phases during the observations to drift more than ∼
40% ofan orbital cycle. Nonetheless, we were able to measure the K ve-locity of the primary star to be . ± . km s − (see Table 5). © 2020 RAS, MNRAS , 000–000 S. Rappaport et al.
Table 6.
Spectrally Determined Stellar Properties of the PrimaryParameter Value T eff [K] ± g [cgs] . ± . / H] / [Fe / H] (cid:12) [dex] − . ± . v sin i [km/sec] ± The radial velocities are plotted in Fig. 6 along with a superposedmodel fit which is discussed in Section 3.4.The RV orbital phase zero comes out to be . ± . or-bital cycles after the TESS reference phase zero, and thus is consis-tent to within an uncertainty of ∼ We also analyzed these same spectra to extract some basic prop-erties of the primary star. Because there is no sign of the coolerbinary component in the available spectra, we simply assumed thatthe primary star dominates the spectral signatures. The combinedspectrum was used to determine the stellar atmospheric parame-ters ( T eff , surface gravity log g , metallicity Fe/H) and also the pro-jected rotational velocity ( v sin i ). During the spectral analysis, theKurucz line list was considered and also the ATLAS9 theoreticalmodel atmospheres (Kurucz 1993) were generated with the SYN-THE code (Kurucz & Avrett 1981). The synthetic spectra werecompared with the combined spectrum to obtain the final atmo-spheric parameters with a χ minimization method. By using thisapproach, first we determined the T eff value from the H γ line whichis very sensitive to T eff . The T eff value was searched over the rangeof − K and for T eff (cid:54) log g was fixed to be4.0 (cgs) because the hydrogen lines weakly depend on log g forstars having T eff values over 8000 K (Smalley et al. 2002). The re-sult was a derived T eff of 7900 ±
150 K. The comparison of thesynthetic and the observed H γ line is shown in the lower panel ofFig. 7.By then fixing the derived T eff and also the microturbulencevalue to be 2 km s − , we also determined the log g , Fe/H and the v sin i parameters of the system by applying the same method in thespectral range of − nm. The resulting parameters are givenin Table 6. The uncertainties of the parameters were calculated bythe procedure used by Kahraman Alic¸avus¸ et al. (2020). The best fitto the observed spectrum is illustrated in the upper panel of Fig. 7.These spectrally inferred stellar parameters are summarized inTable 6. In order to evaluate the binary system parameters of TIC 63328020,we utilised two essentially independent approaches to the analy-sis (see also Kurtz et al. 2020). In the first, we find the stellarmasses, inclination and system age that best yield a match to the kurucz.harvard.edu/linelists.html Figure 7.
Comparison of the synthetic (dashed line) and observed (solidline) spectrum. The best fits to the H γ and other lines are shown in thelower and upper panels, respectively. existing measurements of the spectral energy distribution (SED)and the measured radial velocity of the primary star. In the secondapproach, we model the TESS light curve and simultaneously theradial velocity curve with the phoebe2 binary light curve emula-tor (Prˇsa et al. 2016). Both methods utilise a Markov chain MonteCarlo (MCMC) approach to evaluate the uncertainties in the pa-rameters.
The first method for finding the system parameters utilises threebasic ingredients: (1) the known K velocity for the primary star(see Fig. 6); (2) the measured SED points between 0.4 and 22 µ m;and (3) the Gaia distance (Lindegren et al. 2018).In this part of the analysis we also make use of the MIST ( MESA
Isochrones & Stellar Tracks; Dotter 2016; Choi et al. 2016;Paxton et al. 2011; Paxton et al. 2015; Paxton et al. 2019) evolu-tion tracks for stellar masses between 0.7 and 3.0 M (cid:12) with solarcomposition , in steps of 0.1 M (cid:12) . Both here and in Section 3.3.2, http://viz-beta.u-strasbg.fr/vizier/sed/doc/; see also Table 2. We have chosen solar metallicity for this part of the analysis due to (i)© 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 we utilise the Castelli & Kurucz (2003) model stellar atmospheresfor < T eff < , K in steps of 250 K. A solar chemicalcomposition is assumed.Our approach follows that of Kurtz et al. (2020; and refer-ences therein), but we briefly describe our procedure here for com-pleteness. We use an MCMC code (see, e.g., Ford 2005) that evalu-ates four parameters: the primary mass, M , secondary mass, M ,system inclination angle, i , and the MIST equivalent evolutionaryphase (EEP) of the primary star. The use of EEPs as a fitted param-eter are described in detail in Kurtz et al. (2020).For each step in the MCMC analysis we use the value of M and the EEP value for the primary to find R and T eff , fromthe corresponding MIST tracks, using interpolation for masses be-tween those that are tabulated. That also automatically provides anage, τ , for the star. Since, in this first step of the analysis, we as-sume that the two stars in the binary are coeval and have experi-enced no mass exchange, we use the value of τ to find the EEP forthe secondary. In turn, that yields the values of R and T eff , .We check to see that neither star overfills its Roche lobe, andif one does, then that step in the MCMC chain is rejected.The two masses and the orbital inclination angle determinewhat the K velocity of the primary should be. This is then com-pared to the measured value of ± km s − , and determines thecontribution to χ due to the RV evaluation.Finally, we use R and T eff , , as well as R and T eff , , alongwith interpolated Castelli & Kurucz (2003) model spectra, to fitthe 26 available SED points. Here log g is simply fixed at 4.0. Thevalue of χ for this part of the analysis is added to the contributionfrom the RV match, and a decision is made in the usual way viathe Metropolis-Hastings jump condition (Metropolis et al. 1953;Hastings 1970) as to whether to accept the new step or not.This is done times and the posterior system parame-ters are collected. The parameter posterior distributions are fur-ther weighted according to the derivative of the age with respectto the primary EEP number: dτ /d (EEP) as described in Kurtz etal. (2020). This corrects for the unevenly spaced EEP points withina larger evolutionary category, and across their boundaries.The results of this analysis are summarised in a single plot ofdistributions in Fig. 8. We show the posterior distributions for M , M , R and R , in solar units, while T and T are in units of K, and R /R L is dimensionless ( R L is the radius of the primary’sRoche lobe). For this coeval and no mass exchange scenario, the ra-dius of the lower mass secondary, R is much smaller than for theprimary, R . This results from the fact that if the more massive pri-mary is only somewhat evolved off the ZAMS, then the secondarywith a much lower mass cannot be hardly evolved at all.These results are summarised in the second column of Table7. Here we utilised the same information as in Section 3.3.1 (viz, the K velocity of the primary and the SED points, but we relax the the uncertainty in the interior vs surface composition of the primary starin TIC 63328020, (ii) the weak spectroscopic determination of [ Z/H ] = − . ± . (Table 6), and (iii) the fact that [ Z/H ] = 0 turns out toyield the largest number of acceptable models within the wide range ofplausible Z values that we explored (see Table 7). Solar composition for the MIST tracks we used was defined by Choi et al. (2016) as: X (cid:12) = 0 . , Y (cid:12) = 0 . , and Z (cid:12) = 0 . (taken from Asplund et al. 2009). Figure 8.
Posterior distributions from the MCMC analysis for the casewhere the stars are assumed to be coeval and without any prior mass ex-change. For M , M , R , and R , the x axis is labeled in solar units, while T and T are in units of K, and R /R L is dimensionless. Note, inparticular, that R is small compared to R . Figure 9.
Posterior distributions from the MCMC analysis for the casewhere the stars are not assumed to be coeval and without any prior massexchange. The axis labeling is the same as in Fig. 9. Note, that the distribu-tions are broader than in Fig. 8, and M has gotten larger, as has R . constraint that the two stars must be coeval and have undergone nomass exchange.The MCMC system parameter evaluation results for this caseare summarised by the distributions in Fig. 9. We find three majordifferences from this removal of the coeval constraint: (1) the dis-tributions are considerably broader than in Fig. 8; (2) the mass ofthe primary star has shifted considerably to higher values; and (3)the radius of secondary has nearly doubled.The system parameter results for the case where the no-prior-mass-exchange assumption has been relaxed are summarised in thethird column of Table 7. In Section 3.3.1 above we analysed the basic system parametersfrom an MCMC evaluation of the two masses, the inclination an-gle, and the evolutionary phases (EEP) of the two stars. The fittedquantities were K and 26 SED points, coupled with the Gaia dis- © 2020 RAS, MNRAS , 000–000 S. Rappaport et al.
Table 7.
Derived Parameters for the TIC 63328020 SystemInput Constraints SED + RV a SED + RV b Light curve + RV c Period (days) 1.1057 1.1057 1.1057 K (km s − ) d ± ± ± v sin i (km s − ) e ± ± ...Spectral 26 SED points f
26 SED points f ...Stellar evolution tracks MIST g ... ...Light curve modeling ... ... TESS h Distance (pc) i ±
20 1054 ±
20 1054 ± A V a SED + RV b Light curve + RV c M (M (cid:12) ) . ± .
03 2 . ± .
35 2 . ± . M (M (cid:12) ) . ± .
04 1 . ± .
10 0 . ± . R (R (cid:12) ) . ± .
06 2 . ± .
09 3 . ± . R (R (cid:12) ) . ± .
05 1 . ± .
13 2 . ± . T eff , (K) ±
135 8040 ±
150 8300 ± T eff , (K) ±
140 5660 ±
460 5650 ± i (deg) ± ± . ± . a (R (cid:12) ) . ± . . ± . . ± . R /R L . ± .
02 0 . ± . (cid:38) . K (km s − ) j ± ±
17 204 ± age (Myr) ± ... ... β k ... ... . ± . A l ... ... . ± . Notes. (a) MCMC fits to the measured RV amplitude plus the SED points. The assumption is made that the two stars are coeval in their evolution, and havenot exchanged any mass. We give the same weight to the K ‘data point’ as to any one SED point. We have also tested other weightings (e.g., weighting theone K value several times higher) and it does not change the results significantly. If we had allowed the stellar metallicity of the primary to vary freelyinstead of fixing it at solar, then we could have at best constrained Z (cid:12) / (cid:46) Z (cid:46) Z (cid:12) . The corresponding uncertainties in M , M , R , and system agewould have increased to ± . (cid:12) , ± .
15 M (cid:12) , ± .
15 R (cid:12) and ± Myr, respectively, while R and T eff , would remain unchanged. (b) Same as (a)except that the assumption of no prior mass exchange has been dropped. (c) phoebe2 fit to the TESS orbital light curve plus the RV amplitude. (d) Thiswork (see Sect. 3.4). (e) Determined from the observed spectra (see Table 6). (f) See Fig. 4. (g) Dotter (2016) and Choi et al. (2016). (h) Modelled with phoebe2 . (i) Gaia DR2 (Lindegren et al. 2018). (j) Predicted from the MCMC parameter evaluations. (k) Gravity brightening exponent. (l) Bond bolometricalbedo of the secondary. tance. In Section 3.3.2 we relaxed the coeval constraint on the twostars and fit independently for their masses and radii.We now proceed to analyse the system parameters via simul-taneous fitting of the
TESS orbital light curve as well as the radialvelocity curve using the next-generation Wilson-Devinney code phoebe2 (Prˇsa et al. 2016; Horvat et al. 2018; Jones et al. 2020;Conroy et al. 2020). First, we removed the pulsations from the lightcurve. In addition, a visual inspection of the orbital light curve inFig. 1 shows that the ELV peak following the primary eclipses islower than the preceding ELV peak. This difference can be em-pirically removed by subtracting a simple sinusoid at the orbitalfrequency and of amplitude 4110 ppm. The phasing of this sinu-soid would be correct for the Doppler boosting (DB) effect (Loeb& Gaudi 2003; van Kerkwijk et al. 2010) if the lower luminositysecondary were the source of the DB. From the orbital solutionsalready in hand (see columns 2 and 3 of Table 7), the expected DBamplitude would be (cid:46) ppm for the primary and (cid:46) ppmfor the secondary, and with opposite signs . Thus, the observedamplitude is far too large to be the DB effect (which should be only The Doppler boosting amplitude for this system should be A DB =( α K L − α K L ) /cL tot (van Kerkwijk et al. 2010), where the L ’sare the luminosities and α ’s are the Doppler boosting coefficients in the TESS band. If we take α (cid:39) . ± . , and α (cid:39) . ± . , and L /L tot = 0.9, then we find A DB (cid:39) ppm. ∼
400 ppm net), and we attribute it to spots on the secondary thatare corotating with the orbit. Therefore, we elected to subtract offa sinusoidal component with amplitude 4110 ppm from the lightcurve before carrying out the fitting with phoebe2 .The light-curve fitting procedure utilized the MCMC method-ology outlined in Boffin et al. (2018) and Jones et al. (2019).The component masses, radii and temperatures, and the orbitalinclination were allowed to vary freely over ranges consistentwith the observed SED. The only additional free parameterswere the gravity brightening exponent, β (where T , local = T , pole ( g local /g pole ) β ), of the primary and the (Bond) bolomet-ric albedo (Horvat et al. 2019), A , of the secondary, which arecritical for constraining the ELV and irradiation effect amplitudes,respectively.The best-fitting phoebe2 orbital light curve is shown inFig. 10 while the corresponding model fit to the radial velocitieswas presented in Fig. 6. The resultant model parameters for thesystem are listed in the last column of Table 7.It is clear that the phoebe2 model provides a remarkablygood fit to both the observed light and radial velocity curves, withall model variables extremely well constrained. The model vari-ables all present with strongly Gaussian posteriors, however sev-eral are strongly correlated. For example, due to the use of a singlephotometric band, the posteriors of the primary and secondary tem-peratures show a weak positive correlation. Likewise, the primary’s © 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 Figure 10.
The model phoebe2 light curve (red) on the observations(black). A sinusoidal component of amplitude 0.4 per cent was removedfrom the
TESS light curve to equalize the two ELV maxima before per-forming the fit (see text for details). radius is positively correlated with its mass, with larger primarymasses necessitating larger primary radii in order to maintain thesame Roche lobe filling factor and thus the same amplitude of ELV.Ultimately, further observations are required to break these corre-lations but, nonetheless, the current data are sufficient to stronglyconstrain the properties of the system (see Kurtz et al. 2020, for fur-ther discussion of the fitting of a single-band ELV light curve for asimilar case, albeit without the observed eclipses of TIC 63328020which provide additional strong constraints).The results for the system parameters derived from the
Phoebe fit to the
TESS are summarised in the fourth column ofTable 7.
In order to determine the formation history of TIC 63328020, wemust first consider whether or not mass transfer between the binarycomponents has had a significant effect on its evolution. The secondcolumn of Table 7 lists the inferred properties of TIC 63328020 thatwere derived based on the SED and RV data under the assumptionof no mass transfer during the binary’s evolution. With this latterconstraint relaxed, the properties of the binary deduced using (i)the SED and RV data, and (ii) the RV data in conjunction with a phoebe2 fit to the light curve, are shown in the third and fourthcolumns, of Table 7, respectively.The most glaring discrepancy between any of the predictedproperties of the two stars occurs for the radius of the secondary( R ). The inference for the radius made under the assumption ofno mass transfer disagrees with the other two inferences by nearlya factor of two. Given that the last two inferences were derivedindependently (although they do share the same RV data) andgiven that no constraint on mass transfer was imposed, the rela-tively good agreement between these two cases seems to imply thatTIC 63328020 very likely experienced mass transfer in the past.Moreover, an orbital period on the order of days is typical of many‘Algol-like’ binaries for which mass transfer/loss occurred duringtheir prior evolution (see, e.g., Batten 1989; Eggleton 2000).In analysing the evolution of the progenitor binary we willtherefore assume that mass transfer occurred. The next questionto address is whether or not the current primary was the originalprimary of the progenitor system (i.e., the more massive one) orwhether a mass-ratio reversal occurred (i.e., an Algol-like evolu- tion). If a mass-ratio reversal did not occur, then we are forced toconclude that the original primary could not have lost much masssimply because the original mass ratio ( M /M ) would have beenso large ( (cid:38) ) that the binary would have undergone a dynamicalinstability leading to the presumed merger of the two stars (see,e.g., Webbink 1976; Soberman, Phinney & van den Heuvel 1997for a discussion of the conditions leading to dynamical instability).For the other scenario, the masses of the two primordial com-ponents (we will refer to them as ‘Star 1’ and ‘Star 2’) mayhave changed significantly during the course of the evolution. Thelargest uncertainty concerns the degree to which mass transfer be-tween the components is conservative. For a fully conservativetransfer, all of the mass that is lost by the more massive star (Star1) is subsequently accreted by Star 2. On the other hand, for com-pletely non-conservative mass transfer, the mass of Star 2 wouldnot change as the mass of Star 1 decreased. Because we do notknow how non-conservative mass transfer could have been (or theamount of angular momentum transported out of the binary), wehave investigated a realistic range of possibilities. To determine the properties of putative primordial binaries thatcould evolve to approximately match the observational propertiesinferred for TIC 63328020, we have computed evolutionary tracksfor an extremely wide range of initial conditions and assumed inputphysics. To optimise the numerical computations, we were guidedby a number of grids that had been previously generated to solvefor the evolution of other types of interacting binaries. Specifically,we have used the grids generated for post-Algol binaries such asMWC882 (Zhou et al. 2018) and wide, hot subdwarf binaries (Nel-son & Senhadji 2019) to try to constrain the range of possible initialconditions. Once this was accomplished, additional (more precise)grids were successively computed until we were able to enumeratea reasonably precise set of primordial binaries that could evolve toproduce reasonable facsimiles of TIC 63328020.The evolutionary tracks were calculated using the binary ver-sion of
MESA for which the evolution of both the donor and ac-cretor stars are computed simultaneously (see Paxton et al. 2011;Paxton et al. 2015; Paxton et al. 2019). The grids cover a range ofinitial conditions describing the properties of the primordial bina-ries. Specifically, we created grids for primary masses (i.e., Star 1,the more massive component) in the range of (cid:54) M , / M (cid:12) (cid:54) .The mass of the secondary (Star 2) was expressed in terms of themass ratio ( q ) of the primary’s mass to the secondary’s mass. Weexplored the range of . (cid:54) q (cid:54) . Finally, the primordial orbitalperiod was expressed in terms of the critical period ( P c ) for whichthe primordial primary would just be on the verge of overflowingits Roche lobe. Orbital periods in the range of (cid:54) P orb /P c (cid:54) were computed. In all, some 1800 new binary evolution modelswere generated, in addition to the original (cid:39) that we alreadyhad in our library of computations for Algol-like systems.We also investigated the effects of metallicity. Given the bi-nary’s proximity to the mid-plane of the Galaxy, we chose valuesin the range of . (cid:54) Z (cid:54) . which is a reasonable rangefor Population I stars. We found that a metallicity of Z = 0 . The results presented in this paper were computed with release 10108. This value of Z does not have to be same as used to compute the 2ndcolumn in Table 7 since those results, which assume a coeval evolution, turnout to be invalid, regardless of Z .© 2020 RAS, MNRAS , 000–000 S. Rappaport et al.
Figure 11.
Representative tracks of three distinctly different primordial binary evolutions leading to the formation of TIC 63328020 are shown in the HR di-agram. The evolution of the more massive component of the primordial binary (Star 1) is denoted by a solid curve and the evolution of its less massivecompanion (Star 2) is represented by a dashed curve. Arrows superposed on the curves denote the direction of increasing age (i.e., the direction in which thecomponent is evolving). Different colours are used to denote each pair of the components of the primordial binary. Each component initially starts evolvingfrom the ZAMS and those masses for all three cases are labelled. The red tracks correspond to a binary that experiences highly non-conservative mass-transfer( β = 0 . ), while the blue and green tracks assume fully conservative mass-transfer. The solid dots of the same colour as the track indicate the points in theevolution when the orbital period matches the observed value (1.1057 d). The black dots indicate the points for which (rapid) thermal-timescale mass-transferis first initiated. All of the tracks terminate when both stars have simultaneously fill their Roche lobes. The error bars mark the current observed luminosityand T eff for both components of the system. best matched the observations of the effective temperatures, and forthis reason we adopted that value when computing the final grid ofmodels (with X = 0 . and Y = 0 . ) . We will return to thisissue later and discuss our choice for the metallicity.In terms of the input physics, the degree to which mass trans-fer is non-conservative and the mechanism describing the systemicloss of orbital angular momentum is very uncertain. This uncer-tainty can be parametrized in terms of the quantities α and β (fordetails see Tauris & van den Heuvel 2006). In the MESA code, α is the fraction of the mass lost by the donor star that gets ejectedfrom the binary such that the mass carries away the specific angu-lar momentum of the donor star (i.e., fast Jeans’ ejection), and β is the fraction of the mass that is ejected from the accretor and isassumed to carry away the specific angular momentum correspond-ing to that star. Thus the mass gained by the accretor (secondary)can be written as: δM = − (1 − α − β ) δM . (3)We further assume that none of the mass that is lost from the bi-nary forms a circumbinary torus that can extract additional orbitalangular momentum during the binary’s evolution.Given the uncertainty in the values of α and β , our evo-lutionary tracks were computed for a range of values such that (cid:54) α (cid:54) . and (cid:54) β (cid:54) , under the constraint that α + β (cid:54) .We draw the qualitative conclusion that the sum of α + β has a These
MESA values are in excellent agreement with those of Coelhoet al. (2007) as interpolated from their Table 1: X = 0 . , Y = 0 . , Z = 0 . . much greater effect on the evolutionary outcomes than the combi-nation of individual choices of α and β that give the same sum. Thus to minimise numerical computations, our final set of mod-els has been computed with α = 0 . Finally, orbital angular mo-mentum dissipation was calculated based on the torques associatedwith gravitational radiation and magnetic braking as described inGoliasch & Nelson (2015) and Kalomeni et al. (2016), with themagnetic braking index set equal to 3. The magnetic braking for-mula (Verbunt-Zwaan law; Verbunt & Zwaan 1981) was inferredfrom observations of low-mass main-sequence stars and thus mustsometimes be extrapolated to stars that are either evolved, verylow-mass, or rapidly rotating. Although the magnitude of magneticbraking torques remains uncertain, it has relatively little effect onthe evolutionary tracks until after thermal timescale mass transferhas occurred.After generating our grid of binary evolution tracks, we foundthat both conservative and non-conservative evolutions could pro-duce the desired results given the appropriate choices of the pri-mordial masses and the primordial period. Thus we conclude thata fine-tuning of the initial conditions is not required in order to re-produce the observations. Possible evolutionary scenarios can bedivided into two separate classes: (1) the more massive star (Star1) loses a relatively small fraction of its initial mass while thecompanion (Star 2) gains some portion of that mass; or, (2) themore massive primordial star loses a large fraction of its massleading to a mass-ratio reversal (the mass ratio being defined as Note that the choices of the primordial component masses and orbitalperiod will have a profound effect on the evolution.© 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 q = M Star1 /M Star2 ), thus implying that the accretor becomesmore massive than the donor. For either scenario, the evolutioncan be fully conservative ( β = 0 ) or highly non-conservative( β = 0 . ).We find that the first scenario (mass ratio does not change sig-nificantly) never fully reproduces the observationally inferred prop-erties listed in columns 2 and 3 of Table 7. Although this grouping(class) of evolutionary tracks can reproduce most of the propertiesof TIC 63328020, we did not find a primordial binary that could ul-timately produce a secondary star with such a large radius ( ≈ (cid:12) )while simultaneously matching all of the other inferred propertiesof both stars. The problem stems from the following physical con-straints: (i) in order to bloat the accretor sufficiently, mass-transferrates in excess of ∼ × − M (cid:12) yr − are required for extendedperiods of time; and, (ii) binaries with large mass ratios tend toexperience dynamical instabilities when Roche lobe overflow firstcommences. With respect to the latter issue, the very large massratio inferred for TIC 63328020 necessarily implies that the accre-tor could only have gained a few tenths of a solar mass during theevolution (otherwise the initial evolution would have been dynam-ically unstable). And given the required high mass-transfer ratesand the small net accretion, this implies that mass transfer wouldhave occurred over an extremely short interval ( (cid:46) Myr), mak-ing the whole scenario less likely. Moreover, mass-transfer ratesof ∼ − M (cid:12) yr − are expected at the current epoch and thereis little evidence to support such a high value (see the discussionbelow).According to the second scenario, the original primary of theprimordial binary (i.e., the donor) loses so much mass to its accret-ing companion that a mass-ratio reversal occurs (in other words,the observed low mass secondary of TIC 63328020 was originallythe higher mass star). As discussed above, the largest uncertaintyconcerns the choice of β and we attempt to mitigate the effects ofthis uncertainty by creating a grid of models with the variable β taken to be one of the dimensions of parameter space. It is gen-erally expected that the evolution of Algol-like binaries will be atleast mildly non-conservative (see, e.g., Eggleton 2000) and that iswhy we chose to investigate the range (cid:54) β (cid:54) . .The evolution of both binary stars in the Hertzsprung-Russelldiagram for three representative systems is shown in Fig. 11. Theblue curves illustrate the first scenario, while two sets of tracks rep-resent the second scenario—corresponding to extreme values of β ,i.e., β = 0 and 0.8, green and red curves, respectively. For the firstscenario (see the solid and dashed blue curves for the evolution ofthe donor and accretor, respectively), we chose a primordial binaryconsisting of 2.7 and 0.9 M (cid:12) components with an orbital periodof 3.0 d. Possible solutions for TIC 63328020 are denoted by thesolid blue dots. For the second scenario, this fully conservative casehas components initially consisting of 2.25 and 1.85 M (cid:12) stars in a2.35-d orbit (see the green solid and dashed curves, respectively).The highly non-conservative evolution with β = 0 . is denoted bythe red curves.The primordial binary consisted of a 3.5 M (cid:12) donor in a 1.706-d orbit with a 2.09 M (cid:12) accretor. For these latter two sets of evolu-tionary tracks there are two sets of solid dots (green and red) thatdenote possible solutions at the observed orbital period of 1.1057 d.In each case, the latter set of dots (corresponding to a later age)better fits the inferred properties of TIC 63328020 enumerated inTable 7. The solid black dots indicate the onset of (rapid) thermaltimescale mass transfer. For all cases, the evolutionary tracks areseen to abruptly change their trajectories in the HR diagram oncemass transfer commences. The donor stars all tend to evolve to- Figure 12.
The evolution of the masses of the two components of the pri-mordial binary as a function of orbital period is illustrated. Arrows super-posed on the curves denote the direction of increasing age (i.e., the directionin which the component is evolving). Star 1 is always losing mass while Star2 is accreting mass. The colour scheme for the three cases shown in Fig. 11is repeated here. The initial orbital period for each case is also labeled (up-per left) and the observed orbital period (1.1057 d) is denoted by a verticalblack bar. The solid dots of the same colour as the curves indicate the pointsthat could be possible solutions for the properties of TIC 63328020. Thecyan-shaded region represents the approximate range of possible masses ofthe current primary based on the analysis presented in Table 7 (columns3 and 4). The yellow-shaded region corresponds to the possible range ofmasses for the current secondary. All of the curves terminate when bothstars have simultaneously filled their Roche lobes. wards lower luminosities and effective temperatures while the ac-cretors immediately evolve towards higher temperatures and lumi-nosities.Each of the three tracks terminates once the accretor has ex-panded sufficiently to fill its Roche lobe. A summary of the threetracks and the best fit to the inferred data enumerated in columns 3and 4 of Table 7 is presented in Table 8. Note that the subscripts 1and 2 denote the properties of the primary and secondary, respec-tively, for TIC 63328020 at the current epoch. The age is measuredfrom the formation of the primordial binary and log ˙ M indicatesthe (current) mass-transfer rate from the donor star.In order to further elucidate the scenarios associated with ourthree representative tracks, the evolution of the masses of eachcomponent is shown as a function of the orbital period ( P orb ) inFig. 12. The colour coding and the use of solid and dashed lines inaddition to the solid dots have the same meaning as that describedfor Fig. 11. For track © 2020 RAS, MNRAS , 000–000 S. Rappaport et al.
Table 8.
MESA
Model Parameters for the TIC 63328020 SystemModel Parameter Track M Star1 , (M (cid:12) ) a M Star2 , (M (cid:12) ) a P orb , (d) a β b q ( ≡ M , /M , ) c M (M (cid:12) ) d M (M (cid:12) ) d R (R (cid:12) ) d R (R (cid:12) ) d T eff , (K) d T eff , (K) d L bol , (L (cid:12) ) d L bol , (L (cid:12) ) d R /R L , d R /R L , d log ˙ M (M (cid:12) yr − ) e − − − ˙ P orb /P orb ( yr − ) − . E-6 +5 . E-10 − . E-10system age (Myr) 485 935 545
Notes . (a) ‘Star 1’ and ‘Star 2’ refer to the original primordial primary andsecondary, respectively. The subscript ‘0’ indicates the initial system pa-rameters. (b) β is the fraction of mass that is transferred to the accretorbut is ejected from the system with the specific angular momentum of thatstar. The parameter α (not in the Table) is fixed at 0.0 and is the fraction ofmass lost by the donor star that is ejected from the system with the specificangular momentum of the donor star. (c) Initial mass ratio of the primor-dial binary ( M Star , /M Star , ). (d) These are the model parameters for thecurrent-epoch TIC 63328020 system. (e) The total rate of mass lost by thedonor star. produced the inferred radius of the secondary. Instead, both tracks P orb as the donor stars readjust thermally. As the donorstars approach quasi-thermal equilibrium (with a much reducedmass transfer rate), they contract forcing the orbit to shrink . Theaccretors for both cases continuously gain mass with a resulting in-crease in their radii. The tracks terminate once the accretors alsofill their Roche lobes. Based on the analysis of an extensive grid of models, we concludethat there is a wide range of initial conditions that can replicate thecurrently observed properties of TIC 63328020. The simplest typeof evolution wherein the primordial primary loses a few tenths ofa solar mass to a much less massive accretor, although appealing,cannot reproduce all of the inferred properties. But a wide rangeof evolutionary scenarios for which a mass-ratio reversal occurs Note that the green track of the donor star subsequently experiences asecond orbital period minimum. This is due to the thermal re-adjustment ofthe core from purely convective energy transport to a fully radiative mode. (i.e., the primordial primary becomes the less massive secondary)can be accommodated. In particular, if the evolution is highly non-conservative, then the total mass of the primordial binary wouldhave to be considerably more massive than the presently inferredvalue with an initial mass ratio of q (cid:38) . . For more conservativeevolutions, the initial total mass can be much smaller and the massratio much closer to unity. Thus there is a wide range of initial pa-rameters for the primordial binary that can produce robust modelsof TIC 63328020.Although there are significant uncertainties associated withsystemic mass loss and the orbital angular momentum dissipatedas a result of this non-conservative mass transfer, we find thatthe individual choices of the parameters α and β are not nearlyas important as the contribution from their sum. For this reason,we parametrized the effects of non-conservative mass transfer interms of β ( α = 0 ). We conclude that the values of β in therange of 0 to 0.8 can lead to plausible solutions for the proper-ties of TIC 63328020 (see columns 3 and 4 of Table 7). How-ever, the highly conservative models tend to produce primaries withhigher effective temperatures ( (cid:38) K higher). For this reason,we somewhat prefer models for which β (cid:39) . .Another important result to note is that the ‘simple’ evolution-ary scenario implies that mass-transfer rates at the present epochshould be on the order of − M (cid:12) yr − . By way of contrast, the‘mass-reversal’ scenario requires mass-transfer rates that are typi-cally three orders of magnitude smaller (Table 8). We examined thespectrum of TIC 63328020 for P Cygni profiles and found no ev-idence for that feature. This would seem to imply a relatively lowmass transfer rate. We also examined four WISE band observationslooking for any evidence of nebulosity that might be expected dueto a significant wind emanating from the binary. We could not findany hint of nebulosity in that region. There was also no sign of anyNIR nebulosity from the PanSTARRS images. Although not con-clusive, these results seem to hint at a relatively low mass transferrate or one that is not highly non-conservative.One of the very intriguing features of many of our evolution-ary tracks that reproduce robust models of TIC 63328020 is that theaccretor is very close to filling its Roche Lobe. Since it is relativelyunlikely that we would find such a configuration based solely on theobserved pulsational properties, the question arises as to whetherthe binary had already evolved to a point where both stars tem-porarily over/filled their Roche lobes before one of them contractedleading to the currently observed configuration. If one of the starscontracted then it is possible that it could remain in a detached statefor at least a Kelvin time. On the other hand, it is quite possible thatthe binary would have merged were contact to have occurred. Weare currently trying to address this and questions related to the for-mation and evolution of WU Ma binaries using smoothed particlehydrodynamics (SPH; S. Tripathi, L. Nelson, & T. S. Tricco 2020[in preparation]).Finally, we comment on our choice of generating the evolutiontracks with a metallicity of Z = 0 . . In the process of selectingan appropriate Z , we have explored the effects of metallicity on theevolution of representative models describing the observed proper-ties of TIC 63328020. Specifically, we generated grids of modelsfor the mass fraction of metals in the range of . < Z < . corresponding to between 60 per cent and 170 per cent of the esti-mated solar value. We conclude that in order to reproduce the typeof evolution described by Track Z (cid:46) . ) produce secondary massesthat are too massive by factors of 25 per cent compared to what isexpected based on the results presented in Table 7. For Tracks © 2020 RAS, MNRAS , 000–000 sectoral dipole mode in TIC 63328020 and Z = 0 . , we find that T eff increases by ∼
400 Kfor Track ∼ ◦ ,its distance above the galactic mid-plane is only about 20 pc (thescale height of the thin disk being about 300 pc). Based on Gaia’sestimate of the tangential velocity and using our radial velocity ofthe binary’s centre of mass ( γ = 8 . km s − in Table 5), we esti-mate the spatial velocity to be between ∼
15 and 20 km s − . Theseproperties suggest that TIC 63328020 could well be a young, highmetallicity Population I system. It is also worth noting that all ofour evolutionary models – including the low-metallicity ones de-scribed above – suggest a relatively young age of < Z = 0 . as a reasonablevalue for the metallicity. In this work we report the discovery of a short-period binary withtidally-tilted pulsations at ν = 21 . d − . The pulsation ampli-tude varies with orbital phase and is a maximum at orbital quadra-ture, i.e., when the ellipsoidal light variations are at a maximum.The phase of the pulsations rapidly change by more than π radiansaround the time of the primary eclipse, and there is a smaller jumpin phase at the secondary eclipse by about half that amount in theopposite direction. We note that the phase is not a pure π -radianjump because the mode is distorted from a pure sectoral dipolemode.In order to help visualize how the tidally tilted pulsationswould appear to an observer on an circumbinary planet orbitingTIC 63328020, we include a simulation in the form of an MP4video (‘TIC63328020.mp4’). The video is supplied as SupportingInformation for the paper. This same video is also presented inFuller et al. (2020).The pulsating star has M (cid:39) . (cid:12) , R (cid:39) (cid:12) , and T eff , (cid:39) K, while the secondary has M (cid:39) . (cid:12) , R (cid:39) (cid:12) , and T eff , (cid:39) K. Both stars appear to be close to fillingtheir respective Roche lobes. The orbital period is constant to a partin ∼ over the last century. However, the period appears to varyerratically on timescales of weeks to decades. At present we haveno firm explanation for this behavior.We have carried out an investigation of the history of this sys-tem with an extensive set of binary evolution models. We concludethat the most likely scenario is that there has been a prior epoch ofmass transfer which has reduced the mass of the original primary sothat it is currently the low-mass secondary. By contrast, the originalsecondary is now the pulsating primary star. The mass transfer maystill be ongoing with a low mass-transfer rate of ∼ × − M (cid:12) yr − .Although the architecture and evolutionary histories of thethree known tidally tilted pulsators (HD 74423, CO Cam, andTIC 63328020) are unique, they all feature tidally distorted δ Sct pulsators in short-period orbits. Whereas HD 74423 and CO Camfeature axisymmetric tidally tilted pulsations trapped on the L side of the star, the tidally tilted mode in TIC 63328020 exhibitsvery different phase and amplitude modulation, indicative of a non-axisymmetric ( | m | = 1 ) mode that is not completely trapped on ei-ther side of the star (Fuller et al. 2020). Future discoveries of tidallytilted pulsators will likely reveal more diversity amongst this newclass of stars. ACKNOWLEDGEMENTS
We are grateful to an anonymous referee whose comments and sug-gestions helped clarify the presentations in this paper.This paper includes data collected by the TESS mission. Fund-ing for the TESS mission is provided by the NASA Science Missiondirectorate. Some of the data presented in this paper were obtainedfrom the Mikulski Archive for Space Telescopes (MAST). STScI isoperated by the Association of Universities for Research in Astron-omy, Inc., under NASA contract NAS5-26555. Support for MASTfor non-HST data is provided by the NASA Office of Space Sciencevia grant NNX09AF08G and by other grants and contracts.Based on observations made with the Isaac Newton Telescopeoperated by the Isaac Newton Group of Telescopes, which resideson the island of La Palma at the Spanish Observatorio del Roque delos Muchachos of the Instituto de Astrof´ısica de Canarias. The au-thors thankfully acknowledge the technical expertise and assistanceprovided by the Spanish Supercomputing Network (Red Espa˜nolade Supercomputaci´on), as well as the computer resources used: theLaPalma Supercomputer, located at the Instituto de Astrof´ısica deCanarias.G. H. acknowledges financial support from the Polish NationalScience Center (NCN), grant no. 2015/18/A/ST9/00578. D. J. ac-knowledges support from the State Research Agency (AEI) of theSpanish Ministry of Science, Innovation and Universities (MCIU)and the European Regional Development Fund (FEDER) undergrant AYA2017-83383-P. DJ also acknowledges support undergrant P/308614 financed by funds transferred from the SpanishMinistry of Science, Innovation and Universities, charged to theGeneral State Budgets and with funds transferred from the Gen-eral Budgets of the Autonomous Community of the Canary Islandsby the Ministry of Economy, Industry, Trade and Knowledge. Thisresearch was supported by the Erasmus+ programme of the Eu-ropean Union under grant number 2017-1-CZ01-KA203-035562.L. N. thanks the Natural Sciences and Engineering Research Coun-cil (Canada) for financial support through the Discovery Grantsprogram. Some computations were carried out on the supercom-puters managed by Calcul Qu´ebec and Compute Canada. The op-eration of these supercomputers is funded by the Canada Foun-dation for Innovation (CFI), NanoQu´ebec, R´eseau de M´edecineG´en´etique Appliqu´ee, and the Fonds de recherche du Qu´ebec –Nature et technologies (FRQNT). J. A. thanks NSERC (Canada)for an Undergraduate Student Research Award (USRA). D. J. S.acknowledges funding support from the Eberly Research Fellow-ship from The Pennsylvania State University Eberly College ofScience. The Center for Exoplanets and Habitable Worlds is sup-ported by the Pennsylvania State University, the Eberly College ofScience, and the Pennsylvania Space Grant Consortium. M. S. ac-knowledges the financial support of the Operational Program Re-search, Development and Education – Project Postdoc@MUNI(No. CZ.02.2.69/0.0/0.0/16 027/0008360).This project utilized data from the Digital Access to a Sky © 2020 RAS, MNRAS , 000–000 S. Rappaport et al.
Century@Harvard (‘DASCH’) project at Harvard that is partiallysupport from NSF grants AST-0407380, AST-0909073, and AST-1313370. This paper also makes use of the WASP data set asprovided by the WASP consortium and services at the NASAExoplanet Archive, which is operated by the California Instituteof Technology, under contract with the National Aeronautics andSpace Administration under the Exoplanet Exploration Program(DOI 10.26133/NEA9).This project also makes use of data from the Kilodegree Ex-tremely Little Telescope (‘KELT’) survey, including support fromThe Ohio State University, Vanderbilt University, and Lehigh Uni-versity, along with the KELT follow-up collaboration.
Data availability
The
TESS data used in this paper are available on MAST. All otherdata used are reported in tables within the paper. The
MESA binaryevolution ‘inlists’ are available on the
MESA
Marketplace: http://cococubed.asu.edu/mesa_market/inlists.html . REFERENCES
Applegate J. H., 1992, ApJ, 385, 621Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ARA&A,47, 481Batten A. H., 1989, Space Sci. Rev., 50Boffin H. M. J. et al., 2018, A&A, 619, A84Breger M., 1979, PASP, 91, 5—, 2000, in Astronomical Society of the Pacific Conference Se-ries, Vol. 210, Delta Scuti and Related Stars, Breger M., Mont-gomery M., eds., p. 3Cardelli J. A., Clayton G. C., Mathis J. S., 1989, ApJ, 345, 245Castelli F., Kurucz R. L., 2003, in IAU Symposium, Vol. 210,Modelling of Stellar Atmospheres, Piskunov N., Weiss W. W.,Gray D. F., eds., p. A20Choi J., Dotter A., Conroy C., Cantiello M., Paxton B., JohnsonB. D., 2016, ApJ, 823, 102Coelho P., Bruzual G., Charlot S., Weiss A., Barbuy B., FergusonJ. W., 2007, MNRAS, 382, 498Collier Cameron A. et al., 2006, MNRAS, 373, 799Conroy K. E. et al., 2020, ApJS, 250, 34Cutri R. M., et al., 2013, VizieR Online Data Catalog, II/328Dotter A., 2016, ApJS, 222, 8Eggleton P. P., 2000, NewAR, 44, 111Ford E. B., 2005, AJ, 129, 1706Fuller J., Kurtz D. W., Handler G., Rappaport S., 2020, MNRASGaia Collaboration et al., 2018, A&A, 616, A1Goliasch J., Nelson L., 2015, ApJ, 809, 80Green G. M., Schlafly E., Zucker C., Speagle J. S., Finkbeiner D.,2019, ApJ, 887, 93Green G. M. et al., 2018, MNRAS, 478, 651Grindlay J., Tang S., Simcoe R., Laycock S., Los E., Mink D.,Doane A., Champine G., 2009, Astronomical Society of the Pa-cific Conference Series, Vol. 410, DASCH to Measure (and pre-serve) the Harvard Plates: Opening the 100-year Time DomainAstronomy Window, Osborn W., Robbins L., eds., p. 101Handler G. et al., 2020, Nature AstronomyHastings W., 1970, Biometrica, 57, 97Hoffman D. I., Harrison T. E., Coughlin J. L., McNamara B. J.,Holtzman J. A., Taylor G. E., Vestrand W. T., 2008, AJ, 136,1067 Horne K., 1986, PASP, 98, 609Horvat M., Conroy K. E., Jones D., Prˇsa A., 2019, ApJS, 240, 36Horvat M., Conroy K. E., Pablo H., Hambleton K. M., KochoskaA., Giammarco J., Prˇsa A., 2018, ApJS, 237, 26Jones D., Boffin H. M. J., Sowicka P., Miszalski B., Rodr´ıguez-Gil P., Santander-Garc´ıa M., Corradi R. L. M., 2019, MNRAS,482, L75Jones D. et al., 2020, ApJS, 247, 63Kahraman Alicavus F., Soydugan E., Smalley B., Kub´at J., 2017,MNRAS, 470, 915Kahraman Alic¸avus¸ F., Poretti E., Catanzaro G., Smalley B.,Niemczura E., Rainer M., Handler G., 2020, MNRAS, 493, 4518Kalomeni B., Nelson L., Rappaport S., Molnar M., Quintin J.,Yakut K., 2016, ApJ, 833, 83Kochanek C. S. et al., 2017, PASP, 129, 104502Kurtz D. W., 1982, MNRAS, 200, 807—, 1985, MNRAS, 213, 773Kurtz D. W. et al., 2020, MNRAS, 494, 5118Kurucz R., 1993, ATLAS9 Stellar Atmosphere Programs and 2km/s grid. Kurucz CD-ROM No. 13. Cambridge, 13Kurucz R. L., Avrett E. H., 1981, SAO Special Report, 391Lenz P., Breger M., 2005, Communications in Asteroseismology,146, 53Liakos A., Niarchos P., 2017, MNRAS, 465, 1181Lindegren L. et al., 2018, A&A, 616, A2Loeb A., Gaudi B. S., 2003, ApJ, 588, L117Luri X. et al., 2018, A&A, 616, A9Metropolis N., Rosenbluth A. W., Rosenbluth M. N., Teller A. H.,Teller E., 1953, J. Chem. Phys., 21, 1087Morris S. L., 1985, ApJ, 295, 143Morton D. C., Adams T. F., 1968, ApJ, 151, 611Nelson L. A., Senhadji A., 2019, in American Astronomical So-ciety Meeting Abstracts, Vol. 234, American Astronomical So-ciety Meeting Abstracts © 2020 RAS, MNRAS000