Asteroseismology of massive stars with the TESS mission: the runaway Beta Cep pulsator PHL 346 = HN Aqr

TL;DR

This paper presents the first detailed asteroseismic analysis of the runaway Beta Cep pulsator PHL 346 using TESS data, revealing multiple oscillation modes and providing insights into its evolutionary status and interior structure.

Contribution

It is the first to analyze a runaway Beta Cep pulsator with TESS, identifying numerous oscillation modes and constraining its age and internal structure.

Findings

At least 34 oscillation modes detected, including g- and p-modes.

The star's age is constrained to 23 ± 1 million years.

Pulsation models suggest the star is older than initial models indicated.

Abstract

We report an analysis of the first known Beta Cep pulsator observed by the TESS mission, the runaway star PHL 346 = HN Aqr. The star, previously known as a singly-periodic pulsator, has at least 34 oscillation modes excited, 12 of those in the g-mode domain and 22 p modes. Analysis of archival data implies that the amplitude and frequency of the dominant mode and the stellar radial velocity were variable over time. A binary nature would be inconsistent with the inferred ejection velocity from the Galactic disc of 420 km/s, which is too large to be survivable by a runaway binary system. A kinematic analysis of the star results in an age constraint (23 +- 1 Myr) that can be imposed on asteroseismic modelling and that can be used to remove degeneracies in the modelling process. Our attempts to match the excitation of the observed frequency spectrum resulted in pulsation models that were too young. Hence, asteroseismic studies of runaway pulsators can become vital not only in tracing the evolutionary history of such objects, but to understand the interior structure of massive stars in general. TESS is now opening up these stars for detailed asteroseismic investigation.

Tables

Table 1. Table 1: Multifrequency solution for our TESS photometry of HN Aqr. Error estimates for the independent frequencies are given in braces in units of the last significant digit; the errors on the amplitudes are about ± 0.03 plus-or-minus 0.03 \pm 0.03  mmag.

ID	Freq.	Ampl.	S/N
	(d-1 )	(mmag)
$f_1$	0.2241(8)	0.93	5.1
$f_2$	0.2731(7)	0.88	5.0
$f_3$	0.3421(2)	4.00	24.2
$f_4$	0.3881(5)	1.55	9.7
$f_5$	0.4522(10)	0.84	5.5
$f_6$	0.4844(17)	0.75	5.0
$f_7$	0.5260(4)	2.06	14.0
$f_8$	0.6058(4)	1.67	11.8
$f_9$ aaLikely a close doublet (0.679/0.715 d-1 ).	0.697(1)	1.09	8.0
$f_{10}$	1.1461(7)	0.89	7.7
$f_{11}$	1.3589(6)	1.12	10.3
$f_{12}$	1.4919(7)	0.87	8.2
$f_{13}$	5.4562(8)	0.79	12.4
$f_{14}$	6.079(1)	0.39	6.4
$f_{15}$	6.245(2)	0.83	13.8
$f_{16}$	6.267(1)	1.06	17.6
$f_{17}$	6.371(2)	0.38	6.4
$f_{18}$ bbPossibly a close doublet (6.555/6.5666 d-1 ).	6.5653(5)	20.65	349.4
$f_{19}$	7.747(1)	0.60	10.9
$f_{20}$	7.900(1)	0.58	10.7
$f_{21}$	8.084(1)	0.69	12.7
$f_{22}$	8.273(3)	0.29	5.4
$f_{23}$	8.649(1)	0.64	12.1
$f_{24}$	8.829(1)	0.54	10.5
$f_{25}$	9.021(2)	0.26	5.0
$f_{26}$	9.194(1)	0.61	11.9
$f_{27}$	9.486(2)	0.45	8.9
$f_{28}$	9.677(2)	0.38	7.6
$f_{29}$	9.714(2)	0.36	7.3
$f_{30}$	10.043(2)	0.37	7.6
$f_{31}$	10.095(3)	0.29	6.0
$f_{32}$	10.226(2)	0.34	7.0
$f_{33}$	10.380(2)	0.44	9.0
$f_{34}$	11.267(2)	0.37	8.0
2$f_{18}$	13.1306(7)	1.13	26.1
$f_{18}+f_{19}$	14.312(1)	0.19	4.5
$f_{18}+f_{20}$	14.468(1)	0.15	3.6
$f_{18}+f_{21}$	14.649(1)	0.18	4.3

Table 2. Table 2: Determinations of basic parameters of PHL 346.

$T_{\mathrm{eff}}$ (K)	log $g$	Ref.	Method
$22000\pm 900$	$3.4\pm 0.2$ aa$M_v=-3.3\pm 0.4$, $E\left(b-y\right)=0.037$, calibration by Napiwotzki et al. (1993), $H_{\beta }$ from Handler (2011)	Kilkenny et al. (1977)	Strömgren photometry
$21000\pm 1500$	$3.6\pm 0.3$ bb$v\sin i=75\pm 25$ km s-1 , $E\left(b-y\right)=0.030$	Kilkenny & Lydon (1986)	Optical spectroscopy
$22600\pm 1000$	$3.6\pm 0.2$ ccMacroturbulence $\zeta =12\pm 3$ km s-1	Keenan et al. (1986)	Optical spectroscopy
$22900$	$3.882$ dd$M_{bol}=-5.05$	Heynderickx et al. (1994)	Walraven photometry
$22600$	$3.890$ ee$M_{bol}=-4.97$, calibration by North & Nicolet (1990)	Heynderickx et al. (1994)	Geneva photometry
$22300\pm 1000$	$3.7\pm 0.2$ ffMacroturbulence $\zeta =16\pm 5$ km s-1	Ryans et al. (1996)	Optical spectroscopy
$20700\pm 1000$	$3.58\pm 0.10$ gg$v\sin i=45$ km s-1	Ramspeck et al. (2001)	Optical spectroscopy
$21500\pm 900$	4.1 hhlog $g$ derived from photometry, $\left[M/H\right]=0.21\pm 0.09$, $E\left(B-V\right)=0.068\pm 0.010$	Niemczura & Daszyńska-Daszkiewicz (2005)	UV spectroscopy
$23800\pm 900$	$3.6\pm 0.2$ jj$M_v=-3.3\pm 0.4$, $E\left(b-y\right)=0.025$, calibration by Napiwotzki et al. (1993)	Handler (2011)	Strömgren photometry
$22290\pm 450$	$3.84\pm 0.10$ kk$v\sin i=30.3\pm 0.3$ km s-1 , macroturbulence $\zeta =18\pm 2$ km s-1 , microturbulence $\xi =8\pm 1$ km s-1 , $\left[M/H\right]\approx 0.3$, $RV=60\pm 5$ km s-1 (heliocentric)	this work	Optical spectroscopy
$22400\pm 300$	$3.80\pm 0.06$ ll$v\sin i=26\pm 4$ km s-1 , $\zeta =20\pm 7$ km s-1 , $\xi =14\pm 1$ km s-1 ; abundances $C=8.78\pm 0.07$, $N=8.20\pm 0.05$, $O=8.90\pm 0.08$, $Si=7.76\pm 0.07$, and $Fe=7.66\pm 0.07$	this work	Optical spectroscopy

Equations

$$a(\nu )=\frac{a_0}{1+(\frac{\nu }{{\nu }_0}{)}^{\lambda }}+P_0,$$

Full text

stars: early-type — stars: individual (HN Aqr) — stars: interiors — stars: kinematics and dynamics — stars: massive — stars: oscillations (including pulsations)

††journal: ApJL††facilities: TESS, ESO VLT-UT2, UVES, SAAO††software: MIARMA (Pascual-Granado et al., 2015), SigSpec (Reegen, 2007)

1 Introduction

The Transiting Exoplanet Survey Satellite (TESS) is a NASA mission whose primary objective is to discover hundreds of transiting planets smaller than Neptune with host stars bright enough for spectroscopic follow-up to measure planetary masses and atmospheric compositions (Ricker et al., 2015). TESS has commenced its almost-all-sky survey of bright stars ($4<I_{c}<13$) in a wide red-bandpass filter. In the first two years of operation, precision time series photometry is obtained for 200 000 pre-selected targets with 2-min cadence, for about 30 million stars every 30 minutes.

The characterization of extrasolar planets requires information about their host stars. One of the methods that yields this information is asteroseismology (e.g., Lundkvist et al., 2018), the study of stellar interiors by utilizing their pulsations as seismic waves. Asteroseismology comes at no additional cost to planet-search photometric missions as the observational technique is essentially the same: high-accuracy time-resolved photometry. Consequently, asteroseismology and exoplanet space missions often are combined (e.g., Michel et al., 2006; Gilliland et al., 2010), just like TESS does.

Because TESS will be the first precision photometry mission that surveys almost the whole sky, some types of stars that were not prime targets for searches for extrasolar planets will now be observed in large amounts. In particular, the asteroseismic potential of hot massive stars does not appear to have been fully exploited yet. Pedersen et al. (2019) give a first overview of what TESS can do for OB star astrophysics. In this Letter, we report a study of the first known $\beta$ Cep pulsator observed with the TESS mission, PHL 346 = HN Aqr = TIC 69925250.

1.1 PHL 346 = HN Aqr

PHL 346 is a star of $V=11.44$ located at a Galactic latitude of $b\approx 58^{o}$. Kilkenny et al. (1977) classified it as spectral type B1, and Keenan et al. (1986) reported a surface gravity consistent with an evolved main-sequence evolutionary status and Pop. I metal abundances. Given the radial velocity they measured ($+66\pm 10$ km s*-1* ), Keenan et al. (1986) had to conclude that PHL 346 would not have had enough time to attain such a high Galactic latitude within its lifetime and thus must have been formed far from the Galactic plane. This puzzling result was amended by Ramspeck et al. (2001) who, based on a new spectroscopic analysis and the first proper motion measurement of PHL 346, reconciled the stellar flight time with its lifetime, meaning the star could have been formed in and ejected from the Galactic plane.

PHL 346 was the subject of several studies that resulted in determinations of its effective temperature and surface gravity, summarized in Table 2 and available online only. $\beta$ Cep-type pulsations of PHL 346 were discovered by Waelkens & Rufener (1988) and confirmed by Kilkenny & van Wyk (1990). The star was also observed during the ASAS survey (Pigulski & Pojmański, 2008) and by Handler & Shobbrook (2008). All these authors detected the same single oscillation frequency near 6.566 d*-1* . Heynderickx et al. (1994) and Cugier et al. (1994) suggested this pulsation is due to a nonradial $l=1$ mode. On the other hand, Handler & Shobbrook (2008) derived $l=2$ or 4 for this oscillation, but noted a possible problem with the $U$ filter data their identification critically hinged upon.

2 Observations

PHL 346 was observed with the TESS mission in Sector 2, from August 23, 2018 to September 20, 2018 in 2-min cadence. 18317 brightness measurements were secured over a time base of $\Delta T=27.4$ d, with a duty cycle of 92.8%. The photometry was downloaded from MAST111https://mast.stsci.edu/portal/Mashup/Clients/Mast/Portal.html, and the PDC_SAP fluxes were converted into magnitude. No further manipulations with the data were made. The second part of the light curve is shown in Fig. 1. Clearly, HN Aqr is not a singly-periodic variable. Variations in the mean light level that are not singly periodic as well as changes in the amplitude of the dominant short-period signal are visible.

3 Analysis

3.1 Periodicities in the TESS light curve

Given that HN Aqr is the first $\beta$ Cep pulsator observed with TESS, given that its light variations are multiperiodic, and given that the time base of the observations may not allow all pulsational signals to be resolved in frequency, we analyzed the data with various methods. Method 1 comprised classical single-frequency power spectrum analysis and simultaneous multi-frequency sine-wave fitting with input parameter optimization. The sine-wave fits are subtracted from the data and the residuals examined for the presence of further periodicities (“prewhitening”). The noise level was assessed in apparently signal-free frequency regions and the $S/N>4$ criterion by Breger et al. (1993) was adopted to evaluate the significance of signal detection.

Method 2 used essentially the same approach, but the background noise in the Fourier spectrum was modelled in log-log space following Pablo et al. (2017) using the amplitude spectrum from 0.036 to 30 d*-1* prewhitened with the strongest oscillation. Again, the $S/N>4$ criterion was used. Method 3 applied the MIARMA gap-filling method (Pascual-Granado et al., 2015) to improve the spectral window function of the data. This modified light curve was then frequency analyzed with the SigSpec algorithm (Reegen, 2007) that also involves pre-whitening.

Method 4 (Kallinger & Weiss, 2016) employed a Bayesian algorithm allowing a probabilistic assessment of the significance of signal detection. It was run twice, once allowing only a single frequency to represent a peak in the amplitude spectrum (equivalent to multifrequency sine-wave fitting), and the other time allowing close frequency doublets (spacing $\leq 3/\Delta T=0.11$ d*-1* ). Finally, Method 5 used a Morley wavelet transform (Torrence & Compo, 1998). Output from this approach reflected the amplitude modulation of the short-period pulsations (Fig. 1). In the low-frequency region ($f<2$ d*-1* ), considerable changes in the amplitudes of individual signals were visible, in most cases with some degree of regularity indicating multifrequency beating.

To reconcile the outcome of the different methods we have compared the results of the first four techniques that resulted in lists of frequencies and amplitudes. They yielded fairly consistent results for strong, well-resolved signals, but diverged in frequency regions where densely spaced signals were present. The number of detected frequencies strongly depended on the adopted signal detection threshold.

To provide a set of pulsation frequencies that can be reasonably safely used for asteroseismic investigations, we only accepted signals detected by at least three of the methods independently. The amplitudes of these signals had to exceed the noise level $a(\nu)$ (Eq. 1) by factors of 5 (independent signals) and 3.5 (combination frequencies), respectively. The more conservative S/N threshold was chosen to avoid picking up spurious frequencies in a data set of the present size and sampling (Baran, Koen, & Pokrzywka, 2015); the noise level was computed according to Method 2:

[TABLE]

where $P_{0}=1.4388\times 10^{-5}$ mag is a constant Gaussian noise term, $\lambda=2.0$, $a_{0}=7.76\times 10^{-5}$ mag and $\nu_{0}=(2\pi\tau)^{-1},\tau=0.076$ d are the amplitude and characteristic frequency describing the red noise, respectively.

The resulting list of frequencies is given in Table 1. The frequency values and error estimates were adopted from the Bayesian method as it does not rely on prewhitening. The error bars correspond well to the least-squares errors (Montgomery & O’Donoghue, 1999) for well-separated signals, but take systematic uncertainties of closely-spaced frequencies into account. As it is not yet well known how much data processing affects the low-frequency domain, periods longer than 4 d should be treated with caution. Some steps of prewhitening of the TESS data are illustrated in Fig. 2. Our multifrequency fit leaves a residual scatter of 3.4 mmag per point (617 ppm/hr) containing residual systematic (presumably mostly stellar) variability contributing some 15% to the total scatter.

\startlongtable

The low-frequency domain of detected signals extends up to 1.49 d*-1* . As the large number of frequencies cannot be explained by effects of binarity or rotation most, if not all, of these must be due to pulsation. There is no evidence of rotational or binary-induced modulation, or of a regular period spacing. A re-occurring frequency spacing is present within the signals in the p-mode domain ($\Delta f\approx 0.187$ d*-1* ), which may be due to rotational splitting. No multiplets were identified with confidence; the best candidate consists of frequencies $f_{22}-f_{25}$ which could be part of an $l=2$ quintuplet with the $m=0$ mode missing.

3.2 Stability of the period and amplitude

Thirty years have passed since pulsation in PHL 346 was discovered. Thus the long-term stability of the period and amplitude of the dominant pulsation mode can be studied. We gathered all available time-series photometry of the star (Waelkens & Rufener, 1988; Kilkenny & van Wyk, 1990; Handler & Shobbrook, 2008; Pigulski & Pojmański, 2008). Additionally, we include some previously unpublished observations from 1989 and 1990 by one of us (DK), IUE spectrophotometry (Dufton et al., 1996), as well as NSVS (Woźniak et al., 2004), and ASAS-SN (Shappee et al., 2014) data. Each of the individual data sets used spans at least 13 d.

The resulting O – C diagram and amplitudes over time are shown in Fig. 3. Taking into account passband differences between the various individual data sets, the amplitude of the dominant mode remained relatively stable over 30 years at a level of about 22 mmag and dropped only in the late 1980s and early 1990s, down to 11.9 $\pm$ 0.7 mmag in 1990. (As expected, the amplitude in the UV was higher than in the visual.) During this drop, the O – C diagram shows a significant period change: the 1990 value of O – C is $\approx$ 25 minutes off the ephemeris. We checked whether the shape of the O – C diagram can be explained in terms of the light-time effect in a binary system, with no satisfactory solution. In any case, such an orbit needs to be highly eccentric (e\ {\raise-2.15277pt\hbox{\buildrel>\over{\sim}}}\ 0.7) with periastron passage in the mid 1990’s, when unfortunately a gap of almost a decade in the photometric data occurs. This would be also the time when large radial-velocity (RV) changes would occur. We therefore gathered all available radial-velocity data for HN Aqr (Keenan et al., 1986; Kilkenny & Muller, 1989; Hambly et al., 1996; Ramspeck et al., 2001; Lynn et al., 2002), and plot them in the lower panel of Fig. 3. Although there is no large change of RV at the predicted time of the periastron passage, the first two RVs (Kilkenny & Muller, 1989), are significantly different.

Assuming a period of 30 yr for a hypothetical binary orbit of HN Aqr and $M=10$ M*⊙* for the primary, a companion of $M_{2}\sin i=1.15$ M*⊙* would be able to explain the 0.013 d semi-amplitude of the (O–C) variations. However, such a companion would cause an orbital RV semi-amplitude of only 2.2 km s*-1* , more than an order of magnitude less than the observed RV change. To summarize, in view of the large uncertainties of RVs, the presence of pulsations that increase the scatter of RVs, the large gap in photometric data which causes cycle count ambiguities in the 1990s, the inconsistency of the hypothetical light time effect and the radial velocity change, the results as to the binarity of the star are inconclusive. What one can conclude, however, is that the observed pulsational (O–C) variations and RV changes cannot be caused by binarity alone.

3.3 Spectroscopy and kinematics

High-resolution ($R\approx 35\,000$) spectra of HN Aqr were taken with the UVES spectrograph (Dekker et al., 2000) attached to the VLT-UT2 at the European Southern Observatory (ESO) on the night of August 27, 2009. Integration times of 1150 s at both the blue and the red arm led to S/N ratios of about $150-190$ in the range of $304-669$ nm, and around 70 up to 1043 nm. The parameters derived from these spectra using two analysis strategies (Irrgang et al., 2014; Hubeny & Lanz, 2017) are listed in online Table 2.

Summarizing all these results and using bolometric corrections from Flower (1996) suggests that a range of basic parameters $T_{\rm eff}=22300\pm 900$ K, log $g=3.75\pm 0.15$, $M_{\mathrm{bol}}=-5.2\pm 0.3$ and $[M/H]=0.2\pm 0.1$ should contain the parameter space in which a seismic model for HN Aqr ought to be located. Furthermore, a projected rotational velocity $v\sin i=30\pm 5$ km s*-1* can be adopted. The Gaia DR2 parallax $\pi=0.10\pm 0.08$ mas for HN Aqr (Gaia Collaboration et al., 2016, 2018; Luri et al., 2018) cannot be used to better constrain the stellar luminosity, but is consistent with the results above.

The increased metallicity agrees very well with a kinematic investigation based on the radial velocity, Gaia DR2 proper motions, and a derived spectrophotometric distance of $6.7\pm 0.8$ kpc (99% confidence interval; assuming a mass of $9.5\pm 0.4\,M_{\odot}$). It suggests that the star stems from the inner part of the Galaxy (Galactocentric radius at disk intersection $2.3\pm 0.4$ kpc with 68% confidence interval). After taking the Galactic abundance gradients into consideration, the abundance pattern appears to be perfectly normal except for an underabundance of $\sim 0.3$ dex in carbon. From the calculation of stellar trajectories (for details, see Irrgang et al., 2013), we infer a flight time to the Galactic disk of $23\pm 1$ Myr (68% confidence interval).

3.4 Pulsational modeling

We follow Daszyńska-Daszkiewicz et al. (2017) to model the pulsation spectrum of HN Aqr preliminarily. The basic stellar parameters determined earlier are represented best by models of 11 $M_{\sun}$ that we take for the purpose of example. We attempt to reproduce the frequency domains of the observed p and g-modes (Sect. 3.1) and search for required modifications of the stellar models. Given the low projected rotational velocity of the star, and that the suspected first-order p-mode splitting from Sect. 3.1 would suggest $v_{\rm rot}\approx 64$ km s*-1* , the observed frequency ranges should not differ much from those in the stellar frame of rest. In addition to the frequencies, we consider the age constraint of $23\pm 1$ Myr from the kinematic analysis. The most important influence on mode stability is provided by the overall metallicity, as shown in Fig. 4.

However, increasing $Z$ alone, consistent with the spectroscopically determined abundances (Table 2), which gives a better match to the p-mode domain, is insufficient to reproduce the observed g-mode region. Also, increasing $Z$ requires more massive models that evolve too rapidly: the model with $Z=0.025$ has an age of 13.5 Myr and a main sequence lifetime of 18.0 Myr. Possible ways out of this problem would be the inclusion of convective core overshooting prolonging the main sequence lifetime (a 11 $M_{\sun}$, $Z=0.025$, $\alpha_{\mathrm{o}v}=0.2$, $\log T_{\rm eff}=4.35$ model has an age of 15.5 Myr), increasing the He abundance, leading to models of lower mass, and modifications of the input opacities, in particular near $\log T=5.46$, corresponding to an enhancement of the nickel opacity (see, e.g. Daszyńska-Daszkiewicz et al., 2017).

4 Summary and conclusions

TESS photometry of the pulsating runaway star HN Aqr provided the detection of 38 frequencies of variability (34 independent modes within them) with evidence for more. The star has rich p- and g-mode pulsation spectra. Some of the oscillation frequencies are formally unresolved during the 27.4-d time base of the observations, a problem expected to affect the analysis of other $\beta$ Cep and related types of pulsating stars as well. Hence we applied five different frequency analysis techniques whose combination resulted in a reliable solution.

The signals in the low-frequency domain are dominated by g-mode pulsation. The frequency spectrum is denser at the low-frequency end, which is expected for opacity-driven g modes, but the occurrence of internal gravity waves (e.g., Aerts & Rogers, 2015; Bowman et al., 2019) (or residual instrumental effects) may also be suspected. The p-mode frequency region is surprisingly wide and spans some five radial overtones (cf. Fig. 4); we detected 22 independent frequencies in this domain.

An analysis of archival data showed that the frequency and amplitude of the dominant mode, as well as the radial velocity of the star were not stable over the last 30 years. We could speculate about the presence of a binary companion, but this would be rather surprising because such systems should not survive the ejection of the star from the Galactic disc (Perets & Šubr, 2012), and capture of a companion by a fast-moving runaway star is unfeasible. Furthermore, our Bayesian frequency analysis provided evidence that the dominant pulsation frequency of HN Aqr may be a close doublet, which would provide an alternative interpretation for its amplitude and frequency variations. Therefore, the star should be included in a long-term spectroscopic and photometric observing program.

Pulsational modeling shows that the observed pulsation spectrum cannot be reproduced by increasing the metallicity only; an increase in the opacities in certain stellar interior regions is required. Interestingly, and perhaps most importantly, the nature of HN Aqr as a runaway star provides additional constraints on the modeling, as the age of the models must be reconciled with the stellar flight time ($23\pm 1$ Myr). Our initial attempts in this direction were unsuccessful as we could not obtain models older than 18 Myr that would give a reasonable match to the observed pulsation spectrum.

In asteroseismic modelling of massive stars there are degeneracies between opacity, metallicity, age, mass, and overshooting (Aerts et al., 2018). Having a tight constraint on stellar age obviously will remove at least part of these degeneracies. Therefore, asteroseismic studies of runaway pulsators with precise age determinations may become as important as the studies of pulsators in double-lined eclipsing binaries and can become vital in tracing the evolutionary history of such objects. HN Aqr may not be the only star that can be studied that way (e.g., see Table 1 of Perets 2009).

The analysis of the first known $\beta$ Cep pulsator observed with TESS already demonstrates its potential for massive star asteroseismology. HN Aqr, over the last 30 years believed to pulsate in a single frequency, exhibits at least 34 independent modes. This is the level of progress that space photometry of lower-mass stars has already achieved thanks to the Kepler mission; TESS is now opening the domain of massive stars for in-depth asteroseismology as well.

This paper includes data collected by the TESS mission. Funding for the TESS mission is provided by the NASA Explorer Program. Funding for the TESS Asteroseismic Science Operations Centre is provided by the Danish National Research Foundation (Grant agreement no.: DNRF106), ESA PRODEX (PEA 4000119301) and Stellar Astrophysics Centre (SAC) at Aarhus University. We thank the TESS team and staff and TASC/TASOC for their support of the present work. This work is also based on observations collected at the European Southern Observatory under ESO programme 383.D-0909(A). Funding through the Polish NCN grants 2015/18/A/ST9/00578, 2016/21/B/ST9/01126, 2015/17/B/ST9/02082 and 2014/13/B/ST9/00902 is gratefully acknowledged. GMM acknowledges funding by the STFC consolidated grant ST/R000603/1. The research leading to these results has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement no. 670519: MAMSIE). SS-D acknowledges funding by the Spanish MCIU (projects AYA2015-68012-C2-1-P and SEV2015-0548) and the Gobierno de Canarias (project ProID2017010115). GH thanks Daniel Heynderickx for supplying the photometry by Waelkens & Rufener (1988), David Jones for help in retrieving archival data and Andrzej Baran for helpful comments on the manuscript.

Appendix A Online table