TIC 435850195: The Second Triaxial Tidally Tilted Pulsator
Rahul Jayaraman
1
, Saul A. Rappaport
1
, Brian Powell
2
, Gerald Handler
3
, Mark Omohundro
4
, Robert Gagliano
5
,
Veselin Kostov
2
,
6
, Jim Fuller
7
, Donald W. Kurtz
8
,
9
, Valencia Zhang
10
, and George Ricker
1
1
Department of Physics and Kavli Institute for Astrophysics and Space Research, Massachusetts Institute of Technology, 77 Massachusetts Ave., Camb
ridge, MA
02139, USA;
rjayaram@mit.edu
2
NASA Goddard Space Flight Center, 8800 Greenbelt Rd., Greenbelt, MD 20771, USA
3
Nicolaus Copernicus Astronomical Center, Polish Academy of Sciences, ul. Bartycka 18, PL-00-716 Warszawa, Poland
4
Citizen Scientist, c
/
o Zooniverse, Department of Physics, University of Oxford, Denys Wilkinson Bldg., Keble Rd., Oxford OX1 3RH, UK
5
Amateur Astronomer, Glendale, AZ 85308, USA
6
SETI Institute, 189 Bernardo Ave., Ste. 200, Mountain View, CA 94043, USA
7
TAPIR, Mailcode 350-17, California Institute of Technology, Pasadena, CA 91125, USA
8
Centre for Space Research, North-West University, Mahikeng 2745, South Africa
9
Jeremiah Horrocks Institute, University of Central Lancashire, Preston PR1 2HE, UK
10
Phillips Academy, Andover, MA 01810, USA
Received 2024 August 3; revised 2024 September 3; accepted 2024 September 4; published 2024 October 29
Abstract
The Transiting Exoplanet Survey Satellite
(
TESS
)
has enabled the discovery of numerous tidally tilted pulsators
(
TTPs
)
, which are pulsating stars in close binaries where the presence of a tidal bulge has the effect of tilting the
primary star
’
s pulsation axes into the orbital plane. Recently, the modeling framework developed to analyze TTPs
has been applied to the emerging class of triaxial pulsators, which exhibit nonradial pulsations about three
perpendicular axes. In this work, we report on the identi
fi
cation of the second-ever discovered triaxial pulsator,
with 16 robustly detected pulsation multiplets, of which 14 are dipole doublets separated by 2
ν
orb
. We jointly
fi
t the
spectral energy distribution and TESS light curve of the star, and
fi
nd that the primary is slightly evolved off the
zero-age main sequence, while the less massive secondary still lies on the zero-age main sequence. Of the 14
doublets, we associate eight with
Y
10
x
modes and six with novel
Y
10
y
modes. We exclude the existence of
Y
11
x
modes in this star and show that the observed pulsation modes must be
Y
10
y
. We also present a toy model for the
triaxial pulsation framework in the context of this star. The techniques presented here can be utilized to rapidly
analyze and con
fi
rm future triaxial pulsator candidates.
Uni
fi
ed Astronomy Thesaurus concepts:
Asteroseismology
(
73
)
;
Tidal distortion
(
1697
)
;
Detached binary
stars
(
375
)
1. Introduction
The Transiting Exoplanet Survey Satellite
(
TESS;
G. R. Ricker et al.
2015
)
has revealed diverse classes of
pulsating stars, including the novel class of
“
tidally tilted
pulsators
”
(
TTPs
)
, wherein pulsating stars in tight binaries
(
with periods typically less than
∼
2 days
)
have their pulsation
axes tilted into the orbital plane by the tidal bulge induced by
their companions. Speci
fi
cally, in these systems, the pulsation
axis is aligned with the tidal bulge, rather than with the spin
axis of the star. To date, there have been a handful of TTPs
discovered
(
G. Handler et al.
2020
; D. W. Kurtz et al.
2020
;
S. A. Rappaport et al.
2021
; R. Jayaraman et al.
2022a
;
Z. Jennings et al.
2024b
; V. Zhang et al.
2024
)
, and several
candidates reported in the literature
(
see, e.g., F. Kahraman
Aliçavu
ş
et al.
2024
)
.
A key property of tidally tilted systems is that the observer is
able to study the star through a wide range of latitudinal angles
(
from 0
°
to 360
°
)
with respect to the pulsation axis. Tidally
tilted modes exhibit amplitude and phase modulations around
the orbit, thereby leading to splitting of the modes, seen in
periodograms of their light curves. Unlike typical pulsators,
where the observer
’
s viewing angle remains constant with
respect to the pulsation axis, tidally tilted pulsations yield a
unique perspective on pulsating stars as they orbit their
respective companions.
“
Tidally perturbed
”
pulsators, which also exhibit oscillation
modes whose amplitudes are modulated over the course of the
orbit, could also be interpreted as less extreme manifestations
of the tidal tilting phenomenon
(
see, e.g., J. Southworth et al.
2020
; T. Van Reeth et al.
2022
,
2023
; C. Johnston et al.
2023
;
Z. Jennings et al.
2024a
)
. However, there may also be other
factors at play in these stars causing the observed changes in
the pulsation modes, such as asynchronous rotation
(
see, e.g.,
the introduction of D. M. Bowman et al.
2019
and references
therein
)
or a tidal ampli
fi
cation mechanism, as described in
J. Fuller et al.
(
2020
)
.
Signi
fi
cant modeling of the
fi
rst three TTPs that were
discovered was performed by J. Fuller et al.
(
2020
)
; a unique
attribute of the tidal tilting phenomenon is that the pulsation
modes can be identi
fi
ed with speci
fi
c
ℓ
and
|
m
|
values,
providing valuable insights into both the star
’
s evolutionary
state and its interior structure
(
as was shown in R. Jayaraman
et al.
2022a
)
. Initially, the framework for understanding TTPs
was based on the oblique pulsator model, developed by
D. W. Kurtz
(
1982
)
for rapidly oscillating Ap
(
roAp
)
stars; in
this scheme, the gravitational
fi
eld of the companion star in
TTPs plays an analogous role to the magnetic
fi
eld intrinsic to
roAp stars. However, it remains an open question why only
certain tight binaries experience tidal tilting or, for that matter,
The Astrophysical Journal,
975:121
(
13pp
)
, 2024 November 1
https:
//
doi.org
/
10.3847
/
1538-4357
/
ad77c3
© 2024. The Author
(
s
)
. Published by the American Astronomical Society.
Original content from this work may be used under the terms
of the
Creative Commons Attribution 4.0 licence
. Any further
distribution of this work must maintain attribution to the author
(
s
)
and the title
of the work, journal citation and DOI.
1
why certain modes in a given star can be tidally tilted, while
others in the same star are not.
Recently, V. Zhang et al.
(
2024
, hereafter
Z24
)
reported the
unique discovery of TIC 184743498, the
fi
rst
“
triaxial pulsator
”
(
TAP
)
. This star was found to pulsate along three different
axes, and Section 8 of
Z24
presented a perturbative model for
how modes coupled by the tidal bulge could naturally produce
pulsations along three different axes. This work presents the
second TAP discovered with TESS. Section
2
presents the
details of the TESS data used in the analysis, and our initial
processing of the light curve. Section
3
presents our joint
fi
tto
the system
’
s light curve and spectral energy distribution
(
SED
)
,
and enumerates the system parameters resulting from that
analysis. Section
4
presents an asteroseismic analysis of the
system, and Section
5
presents a toy model of the tidal tilting
framework for this system, while exploring
(
and refuting
)
alternative explanations for this star
’
s behavior.
2. Observations and Data Processing
A team of citizen scientists, including M. Omohundro and
R. Gagliano, have been conducting a visual survey of light
curves from the TESS full-frame images; more information
about this effort can be found in M. H. K. Kristiansen et al.
(
2022
)
. One major goal of this survey is to
fi
nd unique stars in
eclipsing binary
(
EB
)
systems. To this end, B. P. Powell et al.
(
2022
)
generated a set of millions of EB light curves from the
TESS full-frame images using the
eleanor
pipeline
(
A. D. Feinstein et al.
2019
)
and a machine learning approach,
detailed in Section 2 of B. P. Powell et al.
(
2021
)
.
TIC 435850195 was initially identi
fi
ed as part of this effort,
during a review of EB light curves identi
fi
ed by the neural
network from TESS Sector 56. This star was observed for the
fi
rst time in Sector 56 at a 200 s cadence. Further details about
the system are given in Table
1
.
Follow-up analyses utilized the 200 s cadence light curve
from the MIT Quick-look Pipeline
(
QLP; C. X. Huang et al.
2020a
,
2020b
)
. The QLP light curve was downloaded using
lightkurve
(
Lightkurve Collaboration et al.
2018
)
, and we
used the orbital period calculated by the Gaia mission
(
Gaia
Collaboration et al.
2016
,
2023
)
, reported in Table
1
(
C. Siopis
et al. 2024, in preparation
)
, to analyze the eclipses. Orbital
period estimates from both the TESS data and archival data
from the All Sky Automated Survey for Supernovae
(
C. S. Kochanek et al.
2017
)
were consistent with the Gaia
value. Using the TESS data, we also calculated a value for
the ephemeris
t
0
, using the time of the
fi
rst primary eclipse. We
found that
t
0
=
BJD 2459826.36448
(
where BJD is the
“
Barycentric Julian Date
”
in the TDB system, as de
fi
ned in
J. Eastman et al.
2010
)
.
After downloading the normalized light curve, shown in
Figure
1
,we
fi
t 50 harmonics of the orbital period, along with a
constant offset, to the light curve. This procedure was used to
reconstruct the EB light curve for subsequent analysis with our
light-curve
fi
tting code. We used the residuals from this
fi
tas
the
“
pulsational
”
light curve for our asteroseismic analysis in
Section
4
. With this
“
pure
”
pulsational light curve, we then
fi
t
for the most signi
fi
cant frequencies, including their amplitudes
and phases. High-amplitude singlet frequencies, and those
frequencies that formed part of a multiplet, are enumerated,
along with their best-
fi
t amplitude and phase, in Table
4
. The
discrete Fourier transforms
(
DFTs; D. W. Kurtz
1985
)
of each
light curve
(
the raw QLP light curve, the
“
pure
”
pulsational
light curve, and the residuals after all the frequencies were
subtracted
)
are shown in Figure
2
.
The SED was downloaded from the VizieR SED viewer
11
(
F. Ochsenbein et al.
2000
)
. We used the photometric
measurements of this system available from the following
sources: Gaia
(
Gaia Collaboration et al.
2016
,
2023
)
, Pan-
STARRS1
(
K. C. Chambers et al.
2016
)
, the Wide-
fi
eld
Infrared Survey Explorer
(
E. L. Wright et al.
2010
)
, the Two
Micron All Sky Survey
(
R. M. Cutri et al.
2003
; M. F. Skruts-
kie et al.
2006
)
, the Tycho-2 catalog from the Hipparcos
survey
(
E. Høg et al.
2000
)
, the Sloan Digital Sky Survey
(
R. Ahumada et al.
2020
)
, and the Galaxy Evolution Explorer
(
GALEX
)
ultraviolet source catalog
(
L. Bianchi et al.
2017
)
.
3. Estimating Binary System Parameters
We simultaneously
fi
t the orbital light curve and the
composite SED of TIC 435850195 using a custom Markov
Chain Monte Carlo
(
MCMC
)
fi
tting code. There are 14
(
or 16
)
free parameters in the
fi
t: the two masses
(
M
1
,
M
2
)
, the system
age, the orbital inclination angle
i
, eight physically motivated
Table 1
Information about the TIC 435850195 System
Parameter
Value
References
R.A.
(
hh:mm:ss
)
22:54:51.96
(
1
)
Decl.
(
dd:mm:ss
)
20:47:52.17
(
1
)
T
mag
10.551
±
0.009
(
2
)
G
mag
10.7251
±
0.0009
(
1
)
G
BP
−
G
RP
0.3493
±
0.0043
(
1
)
K
1
(
km s
−
1
)
46.13
±
2.61
(
3
)
Distance
(
pc
)
521
±
10
(
1
)
Parallax
(
mas
)
1.801
±
0.034
(
1
)
μ
R.A.
(
mas yr
−
1
)
−
16.91
±
0.03
(
1
)
μ
decl.
(
mas yr
−
1
)
−
12.97
±
0.04
(
1
)
Period
(
days
)
1.36719
±
0.00006
(
1, 4
)
t
0
(
BJD
–
2457000
)
2826.36448
This work
Note.
The reference numbers in column
(
3
)
correspond to
(
1
)
Gaia
Collaboration et al.
(
2016
,
2023
)
,
(
2
)
K. G. Stassun et al.
(
2019
)
,
(
3
)
D. Katz
et al.
(
2023
)
, and
(
4
)
C. Siopis et al.
(
2024, in preparation
)
.
Figure 1.
A
∼
3.5 day portion of the TESS light curve of TIC 435850195 from
the QLP
(
C. X. Huang et al.
2020a
,
2020b
)
, sampled at 200 s cadence. The
amplitudes and phases of most of the pulsation modes are seen to vary
systematically with the orbital phase. There are also prominent ellipsoidal light
variations
(
discussed further in Section
3.1.2
)
and prominent primary eclipses.
11
http:
//
vizier.cds.unistra.fr
/
vizier
/
sed
/
2
The Astrophysical Journal,
975:121
(
13pp
)
, 2024 November 1
Jayaraman et al.
Fourier components of the light curve
(
and an additional
constant offset
)
, the distance to the system, and the line-of-sight
extinction
A
G
. The aforementioned eight Fourier components
correspond to the
fi
rst four orbital harmonics, and represent a
combination of the ellipsoidal light variations
(
ELVs
)
,
illumination effects, and any corotating starspots; physical
explanations for these components are given in Table
2
.
To characterize this system, we make use of the MIST stellar
evolution tracks
(
J. Choi et al.
2016
; A. Dotter
2016
)
, which
yield the stellar radii
(
R
1
,
R
2
)
, effective temperatures
(
T
eff,1
,
T
eff,2
)
, and luminosities
(
L
1
,
L
2
)
as a function of the stellar
masses and system age. As a result, at each step of the MCMC,
we are able to utilize the masses to calculate the orbital
separation
a
(
from the period and Kepler
’
s third law
)
, allowing
us to compute the eclipse geometry. We can also use the stellar
radii and effective temperatures to model the composite SED,
using the model atmospheres of F. Castelli & R. L. Kurucz
(
2003
)
. As part of our analysis, we assume that the two stars
have been evolving in a coeval manner since their formation
—
in particular, we assume that there have been no prior episodes
of mass transfer between them.
The
fi
nal possible free parameters are the distance to the
source and the line-of-sight extinction
A
G
.We
fi
rst ran an
MCMC with both the distance and the extinction as free
parameters, to ensure agreement with the best-
fi
t values from
Gaia. Then, we
fi
xed the distance to 521 pc, and the extinction
to 0.42, and reran the
fi
t. For each value in the SED, we
corrected for extinction using the J. A. Cardelli et al.
(
1989
)
extinction law and this
A
G
value.
In the following subsections, we describe in more detail the
simpli
fi
ed orbital light-curve model and the SED
fi
t. The
MCMC was run with 10 parallel walkers, for 175,000 steps.
We discarded the
fi
rst 20,000 samples in each chain as burn-in.
3.1. Light-curve Fitting
We performed a
fi
t to the folded and binned Fourier-
reconstructed light curve. This light curve was generated using
the coef
fi
cients from the
fi
t to 50 orbital harmonics described in
Section
2
, and was binned to 2 minutes. This approach was
taken in order to minimize any leakage of signals from the
pulsation frequencies into the
“
pure
”
EB light curve. Our
simpli
fi
ed light-curve
fi
tting code incorporates three parts:
(
i
)
eclipses from two spherical stars
12
;
(
ii
)
limb darkening, for
which we utilized the approximation scheme described in the
Appendix
; and
(
iii
)
a series of four sine and four cosine terms
to capture out-of-eclipse variability, including ELVs, the
illumination effect, and starspots that are corotating with the
orbit. The bolometric correction for the TESS band, which is a
function of the stellar effective temperature, was approximated
by assuming a blackbody spectrum and then calculating the
ratio of the integrated
fl
ux between 6000 and 10000
Å
to the
total integrated
fl
ux across all wavelengths. Bolometric
correction values were calculated for blackbodies having
temperatures between 4000 and 10,000 K.
3.1.1. Eclipse Modeling
At each link in the MCMC chain, we used the stellar radii
and luminosities for the eclipse modeling. We set the baseline
fl
ux as the sum of the two luminosities, multiplied by their
respective bolometric corrections. We then calculated the sky-
projected separation between the two stars
’
centers for the
current values of the orbital inclination angle
i
and semimajor
Figure 2.
The DFT
(
D. W. Kurtz
1985
)
of TIC 435850195, with sequential
subtractions of the EB signal in the raw light curve
(
visible in the top panel
)
and the
δ
Sct pulsations
(
visible in the middle panel
)
. Each set of periodic
signals was
fi
t for using least squares and then subtracted from the light curve;
the bottom panel shows the DFT of the residual light curve, after both the EB
signal and pulsations were
fi
tted and removed. Note the differing scales
between the top panel and the lower two.
Table 2
Best-
fi
t Amplitudes to Our Semiempirical Model of the Out-of-eclipse
Variability
Coef
fi
cient
Value
(
ppm
)
Physical Explanation
a
1
−
15,178
±
28
Spots
+
ELV
a
2
−
13,870
±
41
ELV
(
dominant
)
a
3
846
±
19
Spots
+
ELV
a
4
1233
±
23
Spots
b
1
−
83
±
8
Doppler boosting?
b
2
−
1229
±
10
Spots
b
3
−
97
±
10
Spots
b
4
−
15
±
12
Spots
c
0
−
3977
±
25
Offset from 0
Notes.
Two of these parameters
(
a
2
and
b
1
)
correspond to physical attributes of
the system, and are further discussed in Section
3.1.2
. Other parameters
account for corotating spots on the stellar surface.
12
The primary star in this system
fi
lls 50% of its Roche lobe, and the
secondary
∼
30%
(
for details, see Section
3.3
)
. Such stars, in the Roche
potential, are found to be spherical to within
2% in all dimensions, and to
within
1% in the dimensions relevant during eclipses
—
i.e., the
y
-axis and
z
-
axis perpendicular to the tidal axis.
3
The Astrophysical Journal,
975:121
(
13pp
)
, 2024 November 1
Jayaraman et al.
axis
a
, which then allowed us to determine the overlapping area
between the two stars as a function of orbital phase.
During the portions of overlap, we utilized a one-
dimensional integral to calculate the limb darkening
fl
ux
variations for the primary star. We used a quadratic limb
darkening model and the coef
fi
cients from A. Claret
(
2017
)
.
We numerically integrated the limb darkening over the
overlapping area, in a scheme that is a simpli
fi
cation of the
standard approach from K. Mandel & E. Agol
(
2002
)
. The
relatively simple geometry of this system allows for a large
number of the cases
(
III
–
XI
)
in Mandel & Agol to be omitted
for the case of the primary eclipse. Further details of our
approach to calculating the limb darkening can be found in the
Appendix
. For the secondary star, whose limb darkening has a
signi
fi
cantly smaller effect on the light-curve shape, we
calculated an
“
average
”
limb darkening coef
fi
cient for the
overlap region and then applied this across that entire area; this
approach provided a reasonable
fi
t
(
see residuals in Figure
3
)
.
For the primary, we used the limb darkening coef
fi
cients for
a star with
T
eff
=
7500 K,
=
g
log
3.5
, and microturbulent
velocity
ξ
=
2kms
−
1
. For the secondary, we used the
coef
fi
cients for a star with
T
eff
=
4250 K,
=
g
log
4.0
, and
microturbulent velocity
ξ
=
0kms
−
1
. We recognize that there
exists signi
fi
cant debate in the exoplanet community over
whether limb darkening coef
fi
cients should be
fi
xed or free
parameters when
fi
tting transit light curves
(
see, e.g.,
S. Csizmadia et al.
2013
; N. Espinoza & A. Jordán
2015
)
.
Tests with limb darkening coef
fi
cients as free parameters did
not yield a signi
fi
cant improvement in the residuals, so we
fi
xed
them to the Claret values.
Two orbital cycles of the Fourier-reconstructed light curve,
with our best-
fi
tting model overlaid
(
including the out-of-
eclipse
fi
t described in Section
3.1.2
)
, are shown in Figure
3
.
Our model
fi
ts the light curve rather well, with a maximum
residual of approximately 1 ppt.
3.1.2. Out-of-eclipse Variability
For the eclipse modeling, we assumed a
fl
at out-of-eclipse
baseline that was normalized to 1. However, in reality, light
curves of tight binaries exhibit signi
fi
cant out-of-eclipse
variability, for the reasons described earlier. As a result, we
added in a semiempirical model with nine free parameters that
has the following functional form, which is similar to Equation
(
8
)
from J. A. Carter et al.
(
2011
)
:
()
å
ww
++
=
cantbnt
cos
sin
.
1
n
nn
0
1
4
orb
orb
The values of these parameters and their physical signi
fi
cances
are enumerated in Table
2
. With these best-
fi
t parameter values,
we can constrain certain physical attributes of the system.
ELVs
. Stars in tight binaries are distorted into ellipsoids by
the tidal forces from their companion
(
a more detailed
discussion can be found in Section IV.2 of Z. Kopal
1959
)
.
This affects the amount of
fl
ux that we observe, due to the fact
that both the cross-sectional area from which the
fl
ux appears to
be emitted and the gravity darkening vary throughout the orbit
with a frequency of 2
ν
orb
. These contribute to the dominant
term in the ELV
(
w
t
cos 2
)
; terms with frequencies
ν
orb
and
3
ν
orb
also contribute to the ELV, but to a much lesser degree.
Equations
(
14
)
–
(
16
)
from J. A. Carter et al.
(
2011
)
relate the
amplitude of each cosine term
(
a
1
,
a
2
, and
a
3
in our notation
)
to
physical parameters.
The ELVs are analytically related to the physical properties
of the system through the linear limb darkening coef
fi
cient
u
,
and the gravity darkening exponent
y
(
H. v. Zeipel
1924
)
. For
TIC 435850195, we use the
u
and
y
values derived for the
TESS band by A. Claret
(
2017
)
; these are 0.4114 and 0.1245,
respectively. The equation relating
a
2
and the physical
attributes of the system is
⎛
⎝
⎞
⎠
()
()
-
aZq
R
a
i
2sin.
2
21
1
3
2
We present below a simpli
fi
ed version of the expression for
Z
1
given in S. L. Morris
(
1985
)
. Our expression does not include
the
k
1
term in the original equation that accounts for precession
induced by a third body:
()
()
()
()
=
+
-
+=
Z
u
u
y
2
45 3
20 3
1
1.004.
3
1
The dominant source of uncertainty in this expression arises
from the uncertainty in the coef
fi
cients
u
and
y
. We substitute
our value of
Z
1
(
2
)
into Equation
2
, along with the relevant best-
fi
t parameters from Table
3
. Using our calculated value of
a
=
6.96
R
e
,we
fi
nd that
a
2
;
9813 ppm. This is
∼
70% of the
best-
fi
t value for
a
2
in Table
2
, suggesting that the ELV
contribution is the dominant term for the component of the out-
of-eclipse variability with a frequency twice that of the orbit.
Doppler Boosting
. Of particular interest is the
w
bt
sin
1
term,
which corresponds to the Doppler boosting effect in tight
binaries
fi
rst detected by P. F. L. Maxted et al.
(
2000
)
;a
theoretical treatment of this effect is given in A. Loeb &
B. S. Gaudi
(
2003
)
. Doppler boosting is caused by the motion
of the primary inducing three key effects, all of which are
comparable in magnitude: an increased rate of photon arrivals,
an increase in the net energy of the emitted photons, and slight
relativistic beaming due to the motion of the star. Using
synphot
(
STScI Development Team
2018
)
, we calculated the
Doppler boosting coef
fi
cient
α
TESS
to be 2.89 for a 7800 K star
(
assuming a blackbody spectrum
)
; for a Vega-type star, this
coef
fi
cient increases to 2.92.
Figure 3.
Our best-
fi
t model to the binned, folded, Fourier-reconstructed light
curve of TIC 435850195, along with the residuals. The residual structure
around the times of eclipse could be partially explained by our use of
fi
xed limb
darkening coef
fi
cients, rather than letting them be free parameters.
4
The Astrophysical Journal,
975:121
(
13pp
)
, 2024 November 1
Jayaraman et al.