Ellipsoidal Oscillations Induced by Substellar Companions: A Prospect
  for the Kepler Mission
DRAFT VERSION NOVEMBER 1, 2018
Preprint typeset using LATEX style emulateapj v. 08/22/09
ELLIPSOIDAL OSCILLATIONS INDUCED BY SUBSTELLAR COMPANIONS:
A PROSPECT FOR THE KEPLER MISSION
ERIC PFAHL1, PHIL ARRAS1,2 , AND BILL PAXTON1
Draft version November 1, 2018
ABSTRACT
Hundreds of substellar companions to solar-type stars will be discovered with the Kepler satellite. Kepler’s
extreme photometric precision gives access to low-amplitude stellar variability contributed by a variety of
physical processes. We discuss in detail the periodic flux modulations arising from the tidal force on the star
due to a substellar companion. An analytic expression for the variability is derived in the equilibrium-tide
approximation. We demonstrate analytically and through numerical solutions of the linear, nonadiabatic stellar
oscillation equations that the equilibrium-tide formula works extremely well for stars of mass <1.4M⊙ with
thick surface convection zones. More massive stars with largely radiative envelopes do not conform to the
equilibrium-tide approximation and can exhibit flux variations &10 times larger than naive estimates. Over the
full range of stellar masses considered, we treat the oscillatory response of the convection zone by adapting a
prescription that A. J. Brickhill developed for pulsating white dwarfs. Compared to other sources of periodic
variability, the ellipsoidal lightcurve has a distinct dependence on time and system parameters. We suggest that
ellipsoidal oscillations induced by giant planets may be detectable from as many as ∼100 of the 105 Kepler
target stars. For the subset of these stars that show transits and have radial-velocity measurements, all system
parameters are well constrained, and measurement of ellipsoidal variation provides a consistency check, as
well as a test of the theory of forced stellar oscillations in a challenging regime.
Subject headings: planetary systems — stars: oscillations — techniques: photometric
1. INTRODUCTION
The upcoming Kepler3 satellite will continuously monitor
∼105 main-sequence stars of mass ≃0.5–1.5M⊙ over 4–6
years with fractional photometric precisions of ∼10−5. Such
high sensitivity, which is unattainable from the ground, will
allow for the robust detection of Earth-size planets that transit
their host stars, and the measurement of asteroseismic oscilla-
tions as a probe of stellar structure (e.g., Borucki et al. 2004;
Basri et al. 2005). These missions will also discover hundreds
of “hot Jupiters” with orbital periods of <10 days, revealed
by their transits or reflected starlight (e.g., Jenkins & Doyle
2003). Continuous observations of these systems are likely to
show a myriad of novel physical effects, including Doppler
flux variability of the host stars (Loeb & Gaudi 2003), pho-
tometric dips due to moons or rings around the planets
(Sartoretti & Schneider 1999; Brown et al. 2001), and the
impact of additional perturbing planets on transit timing
(Miralda-Escudé 2002; Agol et al. 2005; Holman & Murray
2005). The same ideas apply if the companion is a more mas-
sive brown dwarf, but these are rarely found in close orbits
around solar-type stars (e.g., Grether & Lineweaver 2006).
Here we scrutinize another mechanism for generating peri-
odic variability of a star closely orbited by a giant planet or
brown dwarf. A star subject to the tidal gravity of a binary
companion has a nonspherical shape and surface-brightness
distribution. In the simplest approximation, the stellar surface
is a prolate ellipsoid with its long axis on the line connecting
the two objects. As the tidal bulge tracks the orbital motion,
differing amounts of light reach the observer. For a solar-type
1 Kavli Institute for Theoretical Physics, University of California, Santa
Barbara, CA 93106; pfahl@kitp.ucsb.edu, paxton@kitp.ucsb.edu
2 Department of Astronomy, University of Virginia, P.O. Box 400325,
Charlottesville, VA 22904-4325; arras@virginia.edu
3 http://kepler.nasa.gov
star orbited by a perturbing companion of mass Mp with pe-
riod Porb, the expected fractional amplitude of this ellipsoidal
variability is ∼10−2(Mp/M⊙)(1day/Porb)
2. This effect has a
long history in the study of eclipsing binary stars (see the re-
view by Wilson 1994), but was mentioned only recently in the
exoplanet context.
Udalski et al. (2002), Drake (2003), and Sirko & Paczyński
(2003) noted that if ellipsoidal light variations are detected
from the ground, where the fractional photometric precision
is &10−3, then the perturber must be fairly massive (e.g.,
&0.1M⊙). They offered this idea as a test to distinguish be-
tween planetary transits and eclipses by low-mass stars. The
superior sensitivity of Kepler offers the possibility of mea-
suring ellipsoidal variability induced by giant planets (Mp ∼
10−3–10−2 M⊙) with orbital periods of .10 days.
Loeb & Gaudi (2003) compare the ellipsoidal variability in-
duced by a planetary companion to flux modulations arising
from reflected starlight and the Doppler effect. The three am-
plitudes are similar when the companion has an orbital period
of .3 days and an optical albedo of .0.1. In a sufficiently
long observation it should be possible to separately extract
each of the signals, since their Fourier decompositions are dis-
tinct. Precise physical modeling of the ellipsoidal lightcurve
could provide an independent constraint on the mass of the
companion, as well as important clues regarding stellar tidal
interactions.
Ellipsoidal variability is typically modeled under the as-
sumption that the distorted star maintains hydrostatic balance
and precisely fills a level surface of an appropriate poten-
tial (e.g., the Roche potential). The measured flux is then
just an integral of the intensity over the visible stellar sur-
face, where the intensity includes the effects of limb darken-
ing and gravity darkening (e.g., Kopal 1942). This approach
is strictly valid only when the orbit is circular and the star
rotates at the orbital frequency, so that a stationary configura-
http://arxiv.org/abs/0704.1910v1
http://kepler.nasa.gov
2 PFAHL, ARRAS, & PAXTON
tion exists in the coorbital frame. These conditions may not
be satisfied when the companion has a low mass or long pe-
riod, because of the weak tidal interaction. In fact, a state of
tidal equilibrium may not be attainable in the case of a plan-
etary companion (e.g., Rasio et al. 1996). Equilibrium mod-
els of ellipsoidal lightcurves do have a realm of validity for
noncircular orbits and asynchronously rototating stars, and
have been applied successfully to somewhat eccentric bina-
ries (e.g., Soszynski et al. 2004). However, by construction,
such models ignore fluid inertia and the possibility exciting
normal modes of oscillation, effects that may be of critical
importance in a wide range of observationally relevant cir-
cumstances. Here we apply the machinery of linear stellar
oscillation theory to the weak tidal forcing of stars by sub-
stellar companions. Conceptually, our investigation bridges
Kepler’s planetary and astroseismology programs.
Section 2 describes the geometry of the problem, provides
quantitative measures for the strength of the tidal interaction,
discusses our simplifying assumptions, and presents the math-
ematical framework for calculating ellipsoidal variability. In
§ 3, we consider the equilibrium-tide approximation and de-
rive an analytic expression for the ellipsoidal lightcurve. A
brief review of von Zeipel’s theorem and its limitations is
given in § 4. Tidally forced, nonadiabatic stellar oscillations
are addressed in § 5, where we argue for a simple treatment
of perturbed surface convection zones, use this prescription to
calculate the ellipsoidal variability of deeply convective stars,
estimate analytically the surface flux perturbation in mainly
radiative stars, and show select numerical results. Our main
conclusions are summarized in § 6. We conclude in § 7 with
remarks on the measurement of ellipsoidal oscillations in the
presence of other sources of periodic variability.
2. PRELIMINARIES
Consider a star of mass M and radius R is orbited by a sub-
stellar companion of mass Mp and radius Rp. We work in
spherical coordinates (r,θ,φ) with the origin at the star’s cen-
ter and the pole direction (θ = 0) parallel to the orbital angular
momentum vector. The orbit is then described by (d,π/2,φp),
where d and φp are, respectively, the time-dependent or-
bital separation and true anomaly; φp = 0 marks the phase
of periastron. We assume that the orbit is strictly Keple-
rian with fixed semimajor axis a and eccentricity e, such that
d = a(1 − e2)/(1 + ecosφp). The direction to the observer from
the center of the star is (θo,φo), so that the conventional orbital
inclination is I = π − θo.
We imagine that the gravity of the companion raises nearly
symmetrical tidal bulges on opposite sides of the star that ro-
tate at the orbital frequency. A rough measure of both the
height of the tides relative to the unperturbed stellar radius
and the fractional amplitude of the ellipsoidal variability is
given by the ratio of the tidal acceleration to the star’s surface
gravity:
∼ 10−5
2.8 hr
1 day
, (1)
where MJ ≃ 10
−3 M⊙ is the mass of Jupiter, P∗ =
2π(R3/GM)1/2 = 2.8[(R/R⊙)
3(M⊙/M)]
1/2 hr is the dynami-
cal time of the star. For main-sequence stars with R/R⊙ ≃
M/M⊙, we see that ε ∝ MpMP
orb. The maximum value of
ε is attained when the companion fills its Roche lobe at an
orbital separation of a ≃ 2Rp(M/Mp)
1/3, which gives
εmax ≃
≃ 10−4
0.1R⊙
, (2)
where we have applied a fixed value of Rp = 0.1R⊙, appro-
priate for both giant planets and old brown dwarfs. Note that
εmax ∼ 1 for massive brown dwarfs (Mp/MJ ∼ 80). Hereafter,
we consider only cases with ε≪ 1.
For orbital periods as short as ≃1 day, tidal torques on the
star from a planetary companion are rather ineffective at alter-
ing the stellar rotation rate (e.g., Rasio et al. 1996). Therefore,
as already mentioned in § 1, we should not generally expect
the star to rotate synchronously with the orbit, and so there
is no frame in which the star appears static. This holds when
the orbit is circular, and is obviously true when the there is a
finite eccentricity. In fact, ≃30% of the known exoplanets4
with Porb < 10 days have eccentricities of >0.1. Small vari-
able distortions of the star from its equilibrium state, due to
a combination of asynchronous rotation and orbital eccentric-
ity, should be viewed as waves excited by the tidal force of
the companion. Our task is to study such tidally forced stel-
lar oscillations in the linear domain in order to understand the
corresponding lightcurves.
In order to greatly simplify the mathematical description of
the stellar oscillations, we assume that the star is nonrotating
in the inertial frame. When the stellar rotation frequency is
nonzero, but much smaller than the tidal forcing frequency,
the effect of rotation is to introduce fine structure into the os-
cillation frequency spectrum, and cause the oscillation eigen-
functions to be slightly modified as a result of the Coriolis
force (for a discussion, see Unno et al. 1989). Tidal pump-
ing of a slowly rotating star by an orbiting companion has a
dominant period of Porb/2—a few days in the cases of inter-
est. By contrast, single solar-type stars with ages>1 Gyr tend
to have rotation periods of >10 days (e.g., Skumanich 1972;
Pace & Pasquini 2004); the Sun has an equatorial rotation pe-
riod of ≃25 days. Slowly rotating stars with masses of ≃1M⊙
are prime targets for Kepler, since they exhibit low intrinsic
variability. Based on this selection effect, and the inability of
tidal torques to spin up the star, our assumption of vanishing
stellar rotation seems generally justified.
The general framework for calculating the measurable flux
modulations associated with ellipsoidal stellar oscillations is
as follows. We consider small perturbations to a spherical,
nonrotating background stellar model, such that fluid ele-
ments at equilibrium position x are displaced in a Lagrangian
fashion to position x +ξ. Variations in the measured flux from
an oscillating star arise from two physically distinct contribu-
tions (e.g., Dziembowski 1977): (1) changes in the shape of
the star due to radial fluid displacements ξr = ξ · er, where er
is the radial unit vector, and (2) hot and cold spots generated
by local Lagrangian perturbations ∆F to the heat flux. Our
main task in §§ 3 and 5 is to compute ξr and ∆F according to
the relevant physics.
Given the dependences of ξr and ∆F on (r,θ,φ), it is
straightforward to compute the time varying component of the
measured flux. The flux5 received from a star at distance D is
4 http://vo.obspm.fr/exoplanetes/encyclo/encycl.html
5 Our calculations concern the bolometric flux, although is relatively
http://vo.obspm.fr/exoplanetes/encyclo/encycl.html
ELLIPSOIDAL OSCILLATIONS 3
(e.g., Robinson et al. 1982)
dS n ·no F h(n ·no) (3)
where dS is an area element at the stellar photosphere, F is
the net flux of radiation out of the surface element, h is the
limb-darkening function, n and no are unit vectors normal to
the surface and toward the observer, respectively, and the in-
tegration is over the visible stellar disk. Vertical displacement
at the surface yields changes in F through changes in surface
area and n · no. Following Dziembowski (1977), we expand
ξr and ∆F in spherical harmonics,
ξr(r,θ,φ, t) =
ξr,ℓm(r, t)Yℓm(θ,φ) , (4)
∆F(r,θ,φ, t) =
∆Fℓm(r, t)Yℓm(θ,φ) , (5)
and carry out the appropriate linear expansions to obtain the
fractional variability
(2bℓ − cℓ)
ξor,ℓ
. (6)
Here ξor,ℓ and ∆F
are components evaluated at the surface
(r = R) and in the direction of the observer:
ξor,ℓ
ξr,ℓm(R, t)
Yℓm(θo,φo) , (7)
∆Fℓm(R, t)
Yℓm(θo,φo) . (8)
The terms bℓ and cℓ are given by
dµµPℓ h , cℓ =
dµ(1 −µ2)
where µ = n ·no, the Pℓ(µ) are ordinary Legendre polynomials,
and h(µ) is normalized such that
dµµh = 1. The linear
limb-darkening function is
h(µ) =
(3 −γ)
1 −γ(1 −µ)
; (10)
more general nonlinear functions of µ (e.g., Claret 2000)
will not be considered here. The classical Eddington limb-
darkening function is h = 1 + 3µ/2 (γ = 3/5; e.g., Mihalas
1970). Table 1 shows shows functional forms and particular
values of bℓ and cℓ for ℓ = 2 and 3.
3. EQUILIBRIUM TIDE
Vertical displacement of the stellar surface is often accu-
rately modeled by assuming that the tidally perturbed fluid
remains in hydrostatic balance. The cause and magnitude of
the surface flux perturbation is a more complicated affair. In
this section, we apply a simple parameterization of the flux
perturbation and obtain a complete set of formulae for com-
puting the ellipsoidal lightcurve. Subsequent sections provide
more detailed calculations. In particular, we show in § 5.2
straightforward to modify the analysis for narrow-band measurements.
TABLE 1
LIMB DARKENING PARAMETERS
bℓ cℓ
ℓ General γ γ = 3/5 General γ γ = 3/5
2 (1 +γ)/[20(3 −γ)] 13/40 3(1 + 3γ)[10(3 −γ)] 39/20
3 γ/[4(3 −γ)] 1/16 3γ/(3 −γ) 3/4
that stars with deep convective envelopes (the majority of Ke-
pler targets) have surface flux variations that conform to the
equilibrium-tide approximation.
When the tidal forces on the stellar fluid change sufficiently
slowly, the star can stay very nearly in hydrostatic equilib-
rium. If the net acceleration required to balance the pressure
gradient is derivable from a potential, then equilibrium im-
plies that a fluid element remains on an equipotential surface.
Since we neglect stellar rotation, there is no centrifugal force,
and the total potential is the sum of the gravitational potential
ϕ from the spherical background stellar model and the per-
turbing tidal potential U ∼ εϕ≪ ϕ. For our analytic work,
we neglect the modification of ϕ due to the tide. In general,
the Eulerian variation δϕ should be added to U , as we do in
our numerical models (see § 5.4 and the Appendix); we find
that |δϕ|/|U | ∼ 10−2.
In the absence of tidal forces, a given fluid element sits at
equilibrium position x with total potential ϕ(x). Gentle inclu-
sion of the tidal potential causes the fluid element to move to
position x +ξ while preserving the value of the total potential.
This is expressed mathematically by
ϕ(x) =ϕ(x +ξ) +U(x +ξ, t)
=ϕ(x) +ξ ·∇ϕ+U(x, t) +O(ξ2, ξU) . (11)
We see that ξ ·∇ϕ = ξrg, where g = GMr/r
2 is the background
gravitational acceleration at mass coordinate Mr. To first or-
der, the radial displacement of the equilibrium tide is (see also
Goldreich & Nicholson 1989)
ξr(x, t)≃−U(x, t)/g , (12)
which tells us the geometry of the star as a function of time.
The tidal potential within the star can be expanded as
U(r,θ,φ, t) = −
Pℓ(cosψ) , (13)
where cosψ = sinθ cos(φp −φ). There is no ℓ = 1 term, since
this would give the acceleration of the star’s center of mass,
which is already incorporated into the orbital dynamics. The
angular expansion of ξr follows immediately from eq. (12):
ξr(r,θ,φ, t)
Pℓ(cosψ) . (14)
In order to express U and ξr in spherical harmonics, we utilize
the addition theorem,
Pℓ(cosψ) =
2ℓ+ 1
ℓm(π/2,φp) Yℓm(θ,φ) , (15)
where “∗” denotes the complex conjugate. Note that
Yℓm(π/2,φp) is nonzero only when ℓ−m is even. For the dom-
4 PFAHL, ARRAS, & PAXTON
inant ℓ = 2 components of U and ξr, the surface values of U/ϕ
and ξr/R are ∼ε, as expected.
From eqs. (7), (14), and (15), the components ξor,ℓ/R of the
surface radial displacement toward the observer are immedi-
ately apparent. As we will see in § 5, the computation of
∆F/F is, in general, rather technical. However, in the special
case where the stellar fluid responds adiabatically to a slowly
varying tidal potential, ∆Fℓ/F varies in phase with and in pro-
portion to ξr,ℓ/r in the linear approximation of the equilibrium
tide. Making this assumption, we write ∆Fℓ/F = −λℓξr,ℓ/R at
the surface, where the λℓ are real constants that depend on the
stellar structure (see § 5.1). We will see in § 4 that λℓ = ℓ+ 2
is a good first guess for radiative stars, and so we might gen-
erally expect λℓ to be positive and O(ℓ).
We now have the ingredients for the fractional variability
(eq. [6]), and we obtain
)ℓ−2(
fℓ Pℓ(cosψo) , (16)
where fℓ = (2−λℓ)bℓ−cℓ, and cosψo = sinθo cos(φp −φo). The
ℓ = 2 and 3 Legendre polynomials can be expanded as
P2(cosψo) =
− (3cos2 I − 1)
+ 3sin2 I cos2(φp −φo)
, (17)
P3(cosψo) =
sin I
− 3(5cos2 I − 1)cos(φp −φo)
+ 5sin2 I cos3(φp −φo)
, (18)
where we have substituted θo = π − I. The Eddington limb-
darkening formula gives (see Table 1)
f2 = −
, f3 = −
. (19)
It is important to note that f2 < 0 when λ2 ≥ 0 (see below).
In eq. (16), the orbital dynamics are described by the evolu-
tion of d and φp (see § 2). For a circular orbit, we have d = a
and φp =Ωt, where Ω = 2π/Porb, and t is the time since perias-
tron (modulo Porb). Example lightcurves with e = 0, γ = 3/5,
λℓ = 0, and I = π/2 are shown in Fig. 1 for a/R = {2,4,8,16}.
When R/a ≪ 1, the ℓ = 2 piece of δF/F is a good approxi-
mation, and the temporal flux variation approaches a pure co-
sine with angular frequency 2Ω (see eq. [17]). Because f2 < 0,
the dominant ℓ = 2 component of the ellipsoidal variability has
minimum light when tidal bulge is aligned with the direction
to the observer. As R/a increases, so does the importance
of ℓ > 2 terms and their extra harmonic content, as seen in
eq. (18) and Fig. 1.
Additional harmonics in δF/F also result from a finite
eccentricity. At the O(e) level, signals with frequenciesΩ and
3Ω, and amplitudes of ∼εe, are present in the ℓ = 2 component
of δF/F , which compete with the ℓ = 3 piece when e ∼ R/a.
Notice that when e> 0 the flux is variable even when the orbit
is viewed face-on (I = 0 or π), by virtue of changes in d−3 = 1+
3ecosΩt +O(e2). For I = 0, we see that P3(cosψo) vanishes,
leaving the largest contribution δF/F ≃ −1.5εe f2 cos(Ωt).
4. AN ASIDE ON VON ZEIPEL’S THEOREM
FIG. 1.— Disk-averaged flux variation (eq. [16]) for an edge-on circu-
lar orbit under the equilibrium-tide approximation (eqs. [6] and [14]) with
∆F/F = 0 at the surface. The four curves correspond to a/R = 2 (black), 4
(red), 8 (blue) and 16 (green). In order to compare the shapes of the curves,
δF/F has been multiplied by (a/R)3(M/Mp). As a/R increases, higher har-
monics decrease in strength and the lightcurve approaches a pure cosine with
frequency 2/Porb. The tidal bulge closest to the companion points toward the
observer at integer values of t/Porb.
Our equilibrium calculation in the last section used the sim-
ple prescription ∆Fℓ/F = −λℓξr,ℓ/R. There remains the ques-
tion of what physics determines ∆F/F . A common practice
in empirical studies of close eclipsing binaries—systems that
tend to be nearly in tidal equilibrium—is to use some vari-
ant of the von Zeipel (1924) theorem, which was originally
formulated for purely radiative, strictly hydrostatic stars. In
equilibrium, all the thermodynamic variables depend only on
the local value of the total potential Φ. Thus, the radiative flux
can be written as (e.g., Hansen & Kawaler 1994)
∇Φ , (20)
where ρ is the mass density, T is the effective temperature,
and κ(ρ,T ) is the opacity. Equation (20) is the essence of
von Zeipel’s theorem, which says that the magnitude F of the
radiative flux is proportional to the magnitude of the net ac-
celeration A = |∇Φ|. When Φ = ϕ + U (see § 3), we obtain
A = g + ∂U/∂r +O(ξ2), so that the Lagrangian flux perturba-
tion about equilibrium is
, (21)
where ∆g/g = −2ξr/r, due to the change in radius at approxi-
mately constant enclosed mass. Substituting the equilibrium-
tide result U = −ξrg into eq. (21), we obtain the compact ex-
pression ∆F/F = −∂ξr/∂r. Using eq. (14), we find
= −(ℓ+ 2)
, (22)
from which we identify λℓ = ℓ+ 2.
Although the application of von Zeipel’s theorem is instruc-
tive, the underlying physical assumptions are inaccurate for
ELLIPSOIDAL OSCILLATIONS 5
FIG. 2.— Important oscillations frequencies and time scales as a function
of pressure for a 1M⊙ main-sequence star. The four curves show the Brunt-
Väisällä frequency N (black), Lamb frequency Sl (red; ℓ = 2 is shown), in-
verse thermal time t−1th (blue), and inverse eddy turnover time t
ed (green).
Large, real values of N occur in the radiative core and very near the photo-
sphere, while N2 < 0 in the convective envelope. Gravity waves propagate
only where the angular frequency is below both N and Lℓ. The two horizontal
lines delimit the range of tidal forcing frequencies of interest here.
slowly rotating main-sequence stars of mass 1.0–1.6M⊙ with
tidal forcing periods of days. We are now led to investigate the
general problem of forced nonadiabatic stellar oscillations.
5. FORCED NONADIABATIC OSCILLATIONS
The equilibrium analysis ignores fluid inertia and the exci-
tation of the star’s natural oscillation modes. While this as-
sumption may be valid near the surface of the star, it does
not hold deeper in the interior. Gravity waves (g-modes; re-
stored by buoyancy) can propagate in the radiative interiors
of Sun-like stars with a range of oscillation periods that in-
cludes the tidal forcing periods of interest (.3 days). Tidal
forcing of radiative regions may produce substantial devia-
tions from hydrostatic balance, as well as large surface am-
plitudes of ∆F/F , in particular if resonant oscillations are
excited. This is especially relevant for main-sequence stars of
mass M & 1.4–1.5M⊙ with mainly radiative envelopes. Less
massive stars (M . 1.3–1.4M⊙) have rather deep convective
envelopes that can block information about the dynamic inte-
rior from being conveyed to the surface. Here we investigate
each of these regimes with both analytic estimates and numer-
ical models of oscillating stars.
Our calculations employ realistic models of 0.9–1.6M⊙
main-sequence stars, constructed with the EZ stellar evolu-
tion code (Paxton 2004), a distilled and rewritten version of
the program originally created by Peter Eggleton. We adopt
Solar metallicity and a convective mixing length of 1.6 times
the pressure scale height. All stars are evolved to an age when
the core hydrogen abundance has the Solar value of XH = 0.35.
Models with 199 radial grid points are interpolated to yield
&104 points in which the g-mode radial wavelength is well
resolved in the core.
Figures 2 and 3 illustrate some of the differences between
FIG. 3.— Same as Fig. 2, but for a 1.6M⊙ main-sequence star. Note two
geometrically thin, relatively inefficient (tth ∼ ted) convection zones near the
surface. The spike in N near the center is at the edge of the convective core,
and signals a steep gradient in the mean molecular weight.
1M⊙ and 1.6M⊙ stars, and serve to introduce several impor-
tant physical quantities used in the remainder of this section.
The Lamb frequency,
Sℓ = [ℓ(ℓ+ 1)]1/2
, (23)
is the inverse of the horizontal sound-crossing time scale,
where cs is the sound speed, and [ℓ(ℓ + 1)]1/2/r ≡ kh is the
horizontal wavenumber of the oscillation. For fixed chemical
composition, the squared Brunt-Väisällä frequency is
(∇ad −∇) , (24)
where Hp = −(d ln p/dr)−1 = p/(ρg) is the pressure scale
height, and ∇ = d lnT/d ln p is the temperature gradient6 (∇ad
is the adiabatic value). Radiative regions have ∇ad −∇ > 0
(N2 > 0), and N represents the frequency of buoyancy oscilla-
tions. In convection zones, ∇ad −∇< 0 and N2 < 0, indicating
that g-modes are evanescent. When N2 < 0, the time scale
ted ∼ |N|
−1 , (25)
approximates the turnover time of convective motions (for de-
tails and modifications for radiative losses, see, e.g., Kippen-
hahn & Weigert 1990). A shell of radius r, thickness Hp (size
of the largest convective eddies), and radiative luminosity L
cools on the thermal time scale
tth =
4πr2HpρCpT
, (26)
where Cp is the specific heat at constant pressure.
The 1M⊙ model (Fig. 2) has one deep convection zone with
ted ≪ tth over most of the region, indicating that convection
6 Do not confuse the temperature gradient ∇ with the spatial gradient ∇
used in § 3.
6 PFAHL, ARRAS, & PAXTON
FIG. 4.— Eddy turnover time at base of the convective envelope versus
stellar age for a range of stellar masses.
very efficiently transports energy and causes the zone to be
essentially isentropic. By consrast, the 1.6M⊙ star (Fig. 3)
has two thin surface convection zones with ted ∼ tth, and thus
the radiative and convective fluxes are comparable. Gravity
waves with frequency ω propagate only in radiative regions
where ω < N and ω < Sℓ. For the 1M⊙ star, heat and en-
tropy generated by g-modes in the radiative interior may be
strongly mitigated owing to the long thermal time at the base
of the deep convection zone. On the other hand, g-modes in
a 1.6M⊙ star can propagate very near the surface, producing
qualitatively different results.
We now go on to elucidate the physics of the flux per-
turbations. All the analytic and numerical work that fol-
lows assumes that the tidal potential has the generic form
U ∝ rℓYℓm(θ,φ)exp(−iωt) with forcing frequency ω.
5.1. Heat Transfer in a Convective Envelope
Calculation of the perturbed convective flux in oscillating
stars is a thorny issue. For the purposes of our study, we argue
for an especially simple treatment that draws from previous
work on this subject. Specifically, we modify the prescription
of Brickhill (1983, 1990; see also Goldreich & Wu 1999a,b),
which was originally applied to white-dwarf pulsations, into
a form appropriate for the tidal flow problem.
In the mixing-length theory of convection, heat is trans-
ported by eddies with a spectrum of sizes l . Hp, speeds vl ,
and turnover times ted(l) = l/vl . The Kolmogorov scalings for
turbulent motions give vl ∝ l
1/3, ted ∝ l
2/3, and an energy den-
sity per unit mixing length interval ∝l−1/3. We see that in the
unperturbed star most of the convective energy flux (∝v3l at
scale l) is carried by the largest eddies (l ∼ Hp). Convection
efficiently transports energy when the radiative thermal time
scale associated with the dominant eddies is much longer than
ted. Alternatively, efficient convection implies that the gradi-
ent of the specific entropy s is small; i.e., d lns/d ln p ≪ 1. If
all the convective energy flux F is carried by eddies with mix-
ing length l, the flux and entropy gradient are related by (e.g.,
FIG. 5.— Thermal time at the base of the convection envelope versus stellar
age for a range of stellar masses.
Kippenhahn & Weigert 1990)
d ln p
(l) = (∇−∇ad) ∼
. (27)
Efficient convection enforces ∇ − ∇ad ≪ 1, which implies
d lns/d ln p ≪ 1, since s & Cp in the convective regions of
our background models.
Gravity waves with the tidal forcing frequency ω are ex-
cited in the radiative region below the convection zone. Con-
vective eddies can transport heat during a forcing period only
if ted < 2π/ω (e.g., Brickhill 1990; Goldreich & Wu 1999b).
Inspection of Fig. 2 shows that in the 1M⊙ model, the
largest eddies have ted ≃ 20(p/pbcz)
0.5 days, where pbcz ≃
1013.5 dyne cm−2 is the pressure at the base of convection
zone. Using the Kolmogorov scaling, the “resonant” length
lres for which ωted/2π = 1 is
∼ 10−2
1 day
)3/2(
, (28)
which is >1 for all periods 2π/ω > 1 day when p .
1012 dyne cm−2, which still encompasses much of the convec-
tion zone. Now imagine the situation where all the convective
flux is carried by eddies of size .lres. The entropy gradient
for this range of mixing lengths is
d ln p
(lres)∼ 10
1 day
, (29)
where we have adopted F/pcs ∼ 10
−8 at the base of the con-
vection zone, as indicated by our 1M⊙ stellar model.
These arguments suggest that convection is efficient in a
1M⊙ star at the forcing periods of interest even if small “res-
onant” eddies carry all the energy flux near the base of the
convection zone. At larger radii, but not too near the photo-
sphere, convection is both efficient and rapid (ωted/2π < 1)
ELLIPSOIDAL OSCILLATIONS 7
over the full spectrum of eddies. Rapid convection on all
scales l . Hp enforces isentropy in the convection zone, such
that s and its Lagrangian perturbation ∆s are nearly constant,
as in the Brickhill (1983, 1990) picture. While convection at
the base is rapid only on small scales, it is still highly efficient,
which yields s ≃ constant and further indicates that ∆s/Cp is
small in magnitude, as we demonstrate in § 5.2.
As the stellar mass increases, the convection zone thins and
ted at the base decreases (see Figs. 3 and 4). Rapid convection
holds over the bulk of the convection zone for masses &1M⊙.
However, the assumption that the convection is efficient starts
to break down at 1.4–1.5M⊙, since tth ∼ ted at the base (see
Figs. 3 and 5). For the full range of stellar masses consid-
ered here, we assume that s and ∆s are constant in convection
zones.
5.2. Analytic Result for Thick Convection Zones
In a fully convective star, the emergent luminosity is de-
termined entirely by the surface boundary conditions. Under
our assumption that ∆s is constant in the convection zone, the
perturbed luminosity is likewise a function only of the bound-
ary conditions. Stars of mass .1.3–1.4M⊙ have long thermal
times (ωtth ≫ 1) at the top of the interior radiative region (see
Fig. 5), so that the flux perturbation ∆F is approximately the
“quasi-adiabatic” value, derived by ignoring ∆s ∝ (ωtth)
−1 in
eq. (A5). We assume efficient convection continues to just
below the photosphere.
At the photosphere, we adopt the usual Stefan-Boltzmann
relation, F = σT 4, and the hydrostatic condition, pκ/A = 2/3,
where A is the total acceleration defined in § 4, and 2/3 is the
photospheric optical depth. Taking the photosphere to define
the stellar surface, we compute the Lagrangian perturbations,
, (30)
= 0 . (31)
Using s and p as our independent thermodynamic variables,
we write ∆κ/κ = κad∆p/p +κs∆s/Cp and ∆T/T = ∆s/Cp +
∇ad∆p/p. In our numerical work (see § 5.4), we self-
consistently compute the perturbation ∆A to the effective sur-
face gravity, in order to follow resonant oscillations, where the
equilibrium-tide result fails. However, we are now addressing
non-resonant forcing, for which we use the equilibrium-tide
approximation at the surface, giving ∆A/g = −∂ξr/∂r (see
§ 4). We now have
κs∆s/Cp +∂ξr/∂r
1 +κad
, (32)
and upon substitution,
1 +κad −∇adκs
1 +κad
1 +κad
. (33)
Equation (33) differs from Goldreich & Wu (1999a) in that
we retain the gravity perturbation in eq. (31), whereas they
consider a constant-gravity atmosphere (and no tidal pertur-
bation). For g-modes in white dwarfs, the interesting region
is near the surface and the motion is mainly horizontal, so that
∆A = ∆g = 0 is a good approximation. Since the equilibrium
tide has large vertical motions, the ∆A term must be retained.
The luminosity change across the convection zone is de-
rived from the entropy equation (eq. [A6]). If we ignore
horizontal flux perturbations (set ℓ = 0 in eq. [A6]) and en-
ergy generation, the equation for the luminosity perturbation
∆L/L = 2ξr/r +∆F/F is
d(∆L/L)
= iωT∆s/L . (34)
Integrating over the convection zone with constant ∆s, we
obtain
∆Lbcz
= iω∆s
dMrT/L , (35)
where the subscript “ph” refers to the photosphere. We de-
fine tcz = Cp,ph
dMrT/L to be the mean thermal time of
the convection zone, so that the right-hand side of eq. (35)
is iωtcz∆s/Cp,ph.
Figure 5 shows that the thermal time at the base of the con-
vection zone (of order tcz) for M . 1.3M⊙ is orders of magni-
tude longer than the forcing periods of 1–10 days. Insofar as
|∆L|/L ∼ |ξr|/r at any location in the star (i.e., if resonances
are neglected), we see that (|ξr|/r)
−1|∆s|/Cp ∼ (ωtcz)
−1 ≪ 1
in stars with deep convective envelopes. In this limit, eq. (33)
becomes
1 +κad
. (36)
If we had set ∆A = 0, the amplitude of the photospheric flux
perturbation would have been ∼|∆s|/Cp rather than the much
larger value ∼|ξr|/R.
Photospheric flux perturbations in tidally forced solar-
type stars with thick convective envelopes arise mainly from
changes in the local effective gravity. This statement is remi-
niscent of, but physically distinct from, von Zeipel’s theorem
(eqs. [21] and [22]). We have recovered our equilibrium-tide
scaling, ∆Fℓ/F = −λℓξr,ℓ/R, where eq. (36) gives
λℓ = 4(ℓ+ 2)
1 +κad
. (37)
For M = 1.0–1.4M⊙, we find λ2 ≃ 1.9–1.1. These estimates
neglect resonant excitation of g-modes, a point addressed in
§ 5.4.
5.3. Analytic Result for Radiative Envelopes
As the stellar mass increases beyond 1.4M⊙, the outer
convective region thins and sits close to the surface, where
ted ∼ tth. Figure 3 shows that the 1.6M⊙ model has two thin,
inefficient surface convection zones, as well as a convective
core. Radiative energy transport is important throughout the
envelopes of these more massive stars. We now consider the
idealized case of a completely radiative envelope, and obtain
an analytic approximation for ∆L/L at the surface.
Near the surface of a radiative star, we have Hp/r ≪ 1,
4πr3ρ/Mr ≪ 1, and ω
2r/g ≪ 1 for 2π/ω = 1–10 days. Under
these conditions, the quasi-adiabatic luminosity perturbation
becomes (e.g., Unno et al. 1989)
∆Lqad
gk2hHp
ξr − ξr,eq
8 PFAHL, ARRAS, & PAXTON
where
ζ = κad − 4∇ad +
d lnT
, (39)
and ξr,eq is the equilibrium-tide radial displacement (eq. [12]).
Nonzero values of ∆p/p and (ξr − ξr,eq)/Hp indicate devia-
tions from hydrostatic equilibrium. Care must be taken with
these terms, because the denominators p and Hp become very
small close to the surface.
With the help of the Appendix, we define the variables
α= y1 − y2 + y3 = −
, (40)
β = y2 +
ξr − ξr,eq
, (41)
which satisfy the differential equations
d lnρ
β + (ℓ+ 4)
, (42)
. (43)
When ω2 ≪ gk2hHp, these equations produce the g-mode
dispersion relation k2r = k
2/ω2 (in the limit k2r/k
h ≪ 1)
for radial wavenumber kr. For these propagating waves,
the surface amplitudes of α and β are determined at the
core radiative-convective boundary, where g-modes are driven
(e.g., Goldreich & Nicholson 1989). On the other hand, when
ω2 ≫ gk2hHp, the g-modes are evanescent (see Unno et al.
1989) and we neglect the term gk2hβ/ω
2 in eq. (42). This
limit yields the approximate solution α≃ −(ℓ+ 4)(HpU/gR2),
or ∆p/p ≃ (4 + ℓ)U/gR. In this case, ∆p/p is not small
compared to the fractional fluid displacement, and thus the
equilibrium-tide approximation loses validity.
From our stellar models, we find that the evanescent regime
corresponds to forcing periods of .4–8 days for M = 1–
1.6M⊙, most of the range of interest. The high-frequency
limit of eq. (38) is
∆Lqad
≃ −ζ(ℓ+ 4)
. (44)
This relation should be evaluated at the layer where ωtth ≃ 1,
above which the luminosity effectively “freezes out.” Fig-
ure 6 shows the quasi-adiabic flux perturbation ∆Fqad/F =
∆Lqad/L−2ξr/R, evaluated whereωtth = 1, for a range of forc-
ing periods and M = 1.5–1.7M⊙. Note that |∆F/F| can be an
order of magnitude larger than |U |/gR, because of the rather
large values of |ζ|(ℓ+ 4) for ℓ ≥ 2. Much larger perturbations
are possible when g-modes are resonantly excited in a radia-
tive star, as we discuss in the next section.
We must point out that the quasi-adiabatic approximation is
technically inappropriate when ωtth ∼ 1. Equation (44) should
be viewed as an estimate of the modulus of the luminosity per-
turbation at the surface. If, for instance, |∆s|/Cp & |U |/gR
where ωtth ∼ 1, then ∆L/L at the surface will have a substan-
tial imaginary part (see eq. [34]). This is what we find in the
numerical calculations summarized in the next section.
5.4. Numerical Examples
Here we show solutions of the perturbed mass, momentum,
and energy equations that describe linear, nonadiabatic oscil-
lations of a star subject to a varying tidal force. The equations
FIG. 6.— Ratio of surface Lagrangian flux perturbation ∆F/F to
equilibrium-tide displacement −U/gR for a range of forcing periods in the
limit where surface g-modes are evanescent. The flux is evaluated at the loca-
tion where ωtth = 1. Dashed, solid, and dotted curves correspond to M = 1.5,
1.6, and 1.7M⊙ , respectively.
listed in the Appendix are the same as in Unno et al. (1989)
for radiative regions, but augmented to include the tidal ac-
celeration. In convection zones, we apply the prescription
∆s = constant based on our conclusions in § 5.1. Figure 7
summarizes how the interiors of 1M⊙ and 1.6M⊙ stars re-
spond to resonant and non-resonant tidal forcing. The tidal
potential has been scaled so that ξr/R = 1 corresponds to the
equilibrium-tide surface displacement.
For our 1M⊙ model, the non-resonant response to a forc-
ing period of ≃3 days is shown in Fig. 7a. We see that ξr/R
matches the equilibrium-tide result at the surface; the imag-
inary piece is completely negligible. We also find that our
approximation for ∆F/F at the surface (eq. [36]) works very
well. A factor of ∼10 decay in |∆L|/L occurred in order for
|∆F |/F ∼ ξr/R at the surface. Variation of ∆s/Cp in the con-
vection zone (log[p/(GM2/R4)]> −2.5) is due to changes in
Cp. In the radiative interior, the oscillations are caused by
most nearly resonant g-modes, whose amplitudes rise rapidly
as the core is approached, due to conservation of wave lumi-
nosity. We have checked that the quasi-adiabatic approxima-
tion of ∆L/L is valid in the radiative region; the ratio of the
real and imaginary parts is found to be roughly constant for
the ingoing gravity-wave (see also Zahn 1975).
In order to model the resonant response of a 1M⊙ star, we
tuned the forcing period to ≃1 day (see Figs. 7b and 8). At the
surface, both ξr and ∆L have dominant imaginary parts, due
to the short radial wavelength of the g-mode compared to the
equilibrium-tide fluid displacement. The entropy at the base
of the convection zone is very strongly perturbed in compari-
son to the non-resonant case, but ∆L is still damped by orders
of magnitude as the surface is approached.
Figure 8 shows the surface values of the complex modu-
lus and phase of ξr/R and ∆F/F versus forcing period. The
phase is tan−1(Imaginary/Real) ∈ (−π,π). Solid lines connect
points halfway between g-mode resonant periods. We find
ELLIPSOIDAL OSCILLATIONS 9
FIG. 7.— Responses of tidally forced 1M⊙ and 1.6M⊙ main-sequence stars. Black, red, blue, and green curves denote, respectively, the logarithms of ξr/r,
(δp + ρδϕ)/ρgr, ∆s/Cp, and ∆L/L. Solid (dashed) curves show the real (imaginary) parts. The tidal potential has been scaled so that ξr/R = 1 corresponds
to the equilibrium-tide value. The four panels show the following: (a) non-resonant response of a 1M⊙ star tidally forced at a period of 2π/ω ≃ 2.91day; (b)
resonant response of a 1M⊙ star with 2π/ω ≃ 1.00day; (c) non-resonant response of a 1.6M⊙ star with 2π/ω ≃ 3.00day; (d) resonant response of a 1.6M⊙
star with 2π/ω ≃ 1.02day
that the equilibrium-tide approximation given by eqs. (12) and
(36) is excellent for non-resonant forcing. Dashed curves give
the maximum and minimum values that occur on resonance.
One example of a resonance is shown in the insets. Resonant
forcing at periods of <2 days yields surface values of ξr/R
and ∆F/F that differ substantially from the equilibrium-tide
results. However, the ratio of resonance width to the spac-
ing between adjacent resonances is ∼10−4, making resonant
forcing very unlikely. It is noteworthy that at forcing periods
of >2 days, the equilibrium-tide result holds extremely well
even when precisely on a resonance. As explained by Zahn
(1975), the resonant response can be considered as the sum of
the equilibrium tide and the most nearly resonant wave. As
the period increases, the g-mode radial wavelength decreases,
resulting in a reduction of the overlap integral for the mode
and the tidal force, which in turn gives a decreased amplitude
10 PFAHL, ARRAS, & PAXTON
FIG. 8.— Surface radial displacement and Lagrangian flux perturbation
versus forcing period for a 1.0M⊙ star. Solid lines connect points halfway be-
tween resonant g-mode periods, while dashed curves give the maximum and
minimum values found on resonance. The equilibrium-tide approximation is
extremely good, except when the forcing period is <2 days and resonant.
of the wave component relative to the equilibrium tide.
The non-resonant response of the 1.6M⊙ star is shown
in Fig. 7c. We see that the equilibrium-tide result pro-
vides a good match to ξr/R. Our estimate for the modu-
lus of the radiative luminosity perturbation in the evanes-
cent limit (eq. [44]) agrees reasonably well with what is in
Fig. 7c. We also see that ∆L/L does roughly “freeze-out”
when ωtth ≃ 1, just below the base of the convection zone at
log p/(GM2/R4) ≃ −9 (see Fig. 3). Our expectations in § 5.3
regarding the imaginary part of ∆L/L are borne out in Fig. 7c
A resonantly excited 1.6M⊙ star exhibits huge surface flux
perturbations, radial displacements, and phase lags, as seen in
Fig. 7d. In Fig. 9, surface values of |ξr|/R and |∆F |/F are
plotted as a function of forcing period, where we have taken
care to resolve resonances. Resonant amplitudes vary non-
monotonically with period, in contrast to the smooth behavior
of the 1M⊙ star (Fig. 8). Although we do not show the re-
sults here, similar plots for masses between 1M⊙ and 1.6M⊙
show progressively more structure as the mass increases. The
cause of this irregularity is not clear, but may have to do with
the two thin surface convection zones changing the overlap of
successive g-modes with the tidal force.
6. SUMMARY
We have investigated in detail the ellipsoidal oscillations of
0.9–1.6M⊙ main-sequence stars induced by substellar com-
panions. Classical models of ellipsoidal variability (e.g.,
Wilson 1994) are built on the assumption of hydrostatic bal-
ance in a frame corotating with the binary orbit. This ap-
proach is justified in the context of short-period (Porb .
10 days) binaries containing two stars of comparable mass,
where tidal dissipation circularizes the orbits and synchro-
nizes the stellar spins with the orbital frequency. However,
when the companion has a very low mass, we cannot assume
that the binary is in complete tidal equilibrium; in fact, this
FIG. 9.— Surface radial displacement and Lagrangian flux perturbation
versus forcing period for a 1.6M⊙ star. Curves connect evenly spaced points
away from resonances, with finer spacing near resonance periods.
state may be unattainable (see § 2). In this case, one must, in
general, appeal to a dynamical description of the tidal inter-
action. A substellar companion with Porb & 1 day raises tides
on the star that are a small fraction of the stellar radius (see
eq. [1]), permitting a linear analysis of the stellar oscillations.
While the root of our study is a dynamical treatment of
stellar tidal perturbations, the equilibrium-tide approximation
does have an important realm of validity (see below). For
this reason, we derived in § 3 a general expression (eq. [16])
for the measurable flux variation of a star that remains in
hydrostatic equilibrium under the influence of a small exter-
nal tidal force. This formula (1) assumes that the local per-
turbation to the energy flux at the stellar surface is propor-
tional to and in phase with the equilibrium-tide radial fluid
displacement at each angular order ℓ (eq. [12]), (2) neglects
stellar rotation, and (3) applies to inclined and eccentric or-
bits. As expected, the fractional amplitude of the modulation
is ∼ε≡ (Mp/M)(R/a)
3 for small eccentricities and I = 90◦, or
∼10−5(Mp/MJ)(Porb/1day)
−2 for a star like the Sun (see § 2).
A common practice is to use von Zeipel’s theorem when
computing the surface radiative flux from a tidally distorted
star (see § 4). The theorem assumes that the star is in hydro-
static equilibrium and that the energy transport in the outer
layers is purely by radiative diffusion. As already mentioned,
the hydrostatic assumption is technically unjustified for sub-
stellar perturbers. Moreover, the majority of Kepler targets
will be main-sequence stars with masses of <1.4M⊙, which
have substantial surface convection zones. Evidently, von
Zeipel’s theorem is an inappropriate starting point for the con-
ditions of interest.
Section 5.1 discusses heat transport in perturbed stars with
convective envelopes. Heuristic arguments are used to de-
velop a simple treatment of the perturbed convection zone in
main-sequence stars of mass<1.6M⊙ with forcing periods of
1–10 days. We suggest that both the specific entropy s and its
Lagrangian perturbation ∆s are spatially constant in convec-
ELLIPSOIDAL OSCILLATIONS 11
tive regions, a model partly inspired by the ideas of Brickhill
(1983, 1990).
Using this prescription, we analytically compute in § 5.2 the
perturbed flux at the photosphere of deeply convective stars
(M . 1.4M⊙), where the thermal time scale at the base of the
convection zone is much longer than the forcing period. We
find that ∆s/Cp is negligible near the top of the convection
zone, and that the photospheric flux perturbation is propor-
tional to changes in the effective surface gravity. Thus, we re-
cover the equilibrium-tide result, ∆F/F = −λℓξr/R, at the sur-
face, where λℓ depends on the adiabatic derivatives of opacity
and temperature with respect to pressure (see eq. [37]). Nu-
merical solutions of the equations of linear, nonadiabatic stel-
lar oscillations (see § 5.4 and Fig. 7a) corroborate our analytic
estimates in the non-resonant regime. Resonant excitations of
g-modes in the radiative stellar interior cause large departures
from the equilibrium-tide approximation when the forcing pe-
riod is <2 days (Figs. 7b and 8). However, the likelihood of
being on a resonance is small, and at periods of >2 day the
equilibrium-tide result holds for M ≃ 1M⊙ even with reso-
nant forcing.
Stars of mass &1.4M⊙ have thin, relatively inefficient sur-
face convection zones. Thus, g-modes can propagate very
close to the surface and produce large flux perturbations and
fluid displacements. Analytic arguments in § 5.3 indicate that
the surface flux perturbations in these stars have non-resonant
amplitudes of ∼10ε (eq. [44] and Fig. 6), in rough agree-
ment with our numerical calculations (Fig. 7c). As seen in
Figs. 7d and 9, a resonantly forced 1.6M⊙ star can exhibit
flux perturbation amplitudes of >100ε at forcing periods of
≃1 day. While the amplitudes are not as extreme at longer
periods, their dependence on period is rather erratic (Fig. 9),
an issue that deserves further study. It will be difficult to de-
rive physical interpretations from the ellipsoidal variability of
these more massive stars.
7. DETECTION PROSPECTS
The dominant sources of periodic variability of a star with a
substellar companion are transit occultations (when |cos I| <
[R + Rp]/a), Doppler flux modulations, reflection of starlight
from the companion, and ellipsoidal oscillations. For each
of these signals, Table 2 lists the characteristic amplitude,
period with the largest power in the Fourier spectrum, and or-
bital phase(s) at which the light is a maximum or minimum.
The transit contribution is included for completeness, but its
duration is sufficiently short—a fraction ≃(R + Rp)/(πa) of
Porb—that it should often be possible to excise it from the
data (see Sirko & Paczyński 2003). Of the remaining signals,
the Doppler variability is the simplest, being purely sinusoidal
with period Porb when the orbit is circular. The dominant ℓ = 2
piece of the equilibrium-tide approximation to the ellipsoidal
variability (see eqs. [16] and [17]) is also sinusoidal when
e = 0, but with period Porb/2. Reflection is more problem-
atic, as its time dependence is generally not sinusoidal and
not known a priori.
If the companion scatters light as a Lambert sphere (e.g.,
Seager et al. 2000), the Fourier spectrum of the reflection
variability has finite amplitude at all harmonics of the orbital
frequency Ω, but the amplitude at 2Ω is roughly 1/5 of the
amplitude at Ω. Therefore, the reflection and ellipsoidal vari-
ability amplitudes may be similar at a frequency of 2Ω when
α = 0.1, Mp ∼ MJ , and Porb ≃ 1 day. Also, the orbital phase
at which the reflected light is a maximum is distinct from
both the Doppler and ellipsoidal cases, further distinguishing
the signals. However, Lambert scattering is probably never
appropriate in real planetary atmospheres. Infrared reemis-
sion of absorbed optical light, multiple photon scattering, and
anisotropic scattering typically conspire to narrow the peak in
the reflection lightcurve and lower the albedo, decreasing the
prominence of the reflection signal. These issues are sensi-
tive to the atmospheric chemistry and the uncertain details in
models of irradiated giant planets. For reasonable choices re-
garding the atmospheric composition, calculated optical albe-
dos of Jovian planets range from<0.01 to ≃0.5 (Seager et al.
2000; Sudarsky et al. 2000). Recent photometric observations
of HD 209458, the star hosting the first-detected transiting gi-
ant planet (Porb ≃ 3.5 days), constrain the planetary albedo to
be <0.25 (Rowe et al. 2006).
Detailed lightcurve simulations will be required to say how
well the different periodic signals can be extracted from the
data. This is beyond the scope of the current study. We now do
the simpler exercise of isolating the ellipsoidal modulations
and assessing when this effect alone should be detectable. For
a star of apparent visual magnitude V and an integration time
of T = 6hr, Kepler’s photon shot noise is7
∼ 10−5100.2(V −12)
)−1/2
. (45)
Instrumental noise should contribute at a similar level (e.g.,
Koch et al. 2006). If the data is folded at the orbital period
and binned in time intervals T ≪ Porb, the shot noise is sup-
pressed by a factor of ∼n
orb , where norb is the number of
folded cycles. After folding 1 year of continuous photomet-
ric data using T = 6 hr, a star with V < 12 orbited by a giant
planet with Porb . 3 days may have a fractional shot noise per
time bin of .10−6. This is less than the ellipsoidal amplitude,
(δF/F )ell, when I is not too small.
The actual situation is not so simple when the data spans of
weeks or months, because the intrinsic stochastic variability
of the star will not have a white-noise power spectrum. Over
times of .1 day, the Sun shows variability of (δF/F )int ∼
10−5, but the amplitude rises steeply between ∼1 and 10 days
to ∼10−3. Intrinsic variability tends to be large near the rota-
tion period of the star, due mainly to starspots. Low-frequency
variability may not too damaging for the study of ellipsoidal
oscillations induced by planets with Porb . 3 days, but more
study is needed.
Kepler’s target list will contain ≃105 main-sequence FGK
stars with V = 8–14. The statistics of known exoplanets in-
dicate that 1–2% of all such stars host a giant planet (Mp &
MJ) with Porb < 10 days (e.g., Marcy et al. 2005). Of these
“hot Jupiters,” ≃30% have Porb = 1–3 days. It seems that a
maximum of ∼103 Kepler stars will have detectable ellip-
soidal modulations. If we neglect intrinsic stellar variabil-
ity and consider only shot noise, then many systems with
Porb . 3 days and V < 14 will have signal-to-noise S/N > 1
after ∼100 cycles are monitored; this may amount to >100
stars. Obviously, the number drops when we place higher de-
mands on S/N and include the intrinsic variability. The results
depend critically on the distributions of Mp and Porb.
In order to better estimate the number of stars with poten-
tially detectable ellipsoidal oscillation, we perform a simple
population synthesis calculation. Denote the set of star-planet
7 An integration time of T = 6 hr is chosen for convenience; Kepler’s nom-
inal exposure time is 30 min. Here we use the V -band flux as a reference, but,
in fact, the Kepler bandpass is 430–890 nm, which spans B, V , and R colors.
12 PFAHL, ARRAS, & PAXTON
TABLE 2
PERIODIC FLUX MODULATIONS
Variability Dominant Phase at
Source Amplitudea,b Harmonic Maximum/Minimumc References
Ellipsoidald . . . . . . . 2× 10−5mpmP−21 sin
2 I Porb/2 0.25(0.75)/0.00(0.50) · · ·
Dopplere . . . . . . . . . . 3× 10−6mpm−2/3P
1 sin I Porb 0.25/0.75 1
Reflectionf . . . . . . . . 6× 10−5(α/0.1)m−2/3P
1 sin I Porb 0.50/0.00 2,3
Transit . . . . . . . . . . . 10−2m−2 Porb · · · /0.00 4
REFERENCES. — (1) Loeb & Gaudi 2003; (2) Seager, Whitney, & Sasselov 2000; (3) Sudarsky, Burrows,
& Pinto 2000; (4) Seager & Mallén-Ornelas 2003
a We assume that the orbit is circular in our estimates.
b The dimensionless variables used are mp = Mp/(10−3 M⊙), m = M/M⊙ and P1 = Porb/1day. We have as-
sumed that the star and companion have respective radii of R/R⊙ = m and 0.1R⊙.
c The phase is in the range 0–1, where at phase 0 the planet is closest to the observer.
d Only the ℓ = 2 component of eq. (16), with λ2 = 2, is considered here.
e We approximate the amplitude as 4vr/c, where vr is the reflex speed of the star along the line of sight, and
the factor of 4 is approximately what one obtains for a V -band spectrum similar to the Sun.
f Here α is the geometric albedo of the companion. The inclination dependence is an approximation for I ≃ 90◦
and the Lambert phase function.
system parameters by P = {M,Mp,Porb, I}, and let f (P)dP be
the probability of having a system in the 4-dimensional vol-
ume dP. We assume that the planetary orbits are circular and
obtain (δF/F )ell from the equilibrium-tide estimate in Ta-
ble 2. Given the mass of the star, we compute its absolute
V magnitude on the main-sequence using the approximation
(see also Henry & McCarthy 1993)
MV = 4.8 − 10.3log(M/M⊙) , (46)
which is in accord with the usual mass-luminosity relation
log(L/L⊙) ≃ 4log(M/M⊙) for M ≃ 1M⊙. With a maximum
apparent magnitude of Vmax = 14 for the Kepler targets, the
maximum distance of the star is
Dmax = 10
1+0.2(14−MV ) pc . (47)
With a certain signal-to-noise threshold (S/N)min, there is a
maximum distance Dd < Dmax to which the ellipsoidal vari-
ability is detectable. For a spatially uniform population, the
detectable fraction of systems is (Dd/Dmax)
3. Thus, the net
detectable fraction among all systems is
dP f (P)
, (48)
an integral over all relevant P space.
When the only noise is intrinsic to the star, N = (δF/F )int
and S/N is independent of distance, so that Dd,int = Dmax when
S/N > (S/N)min, and Dd,int = 0 otherwise. In the case of pure
shot noise, there is a maximum magnitude Vd for which the
ellipsoidal oscillations are detectable:
Vd = 5log
(δF/F )ell
χ (S/N)min
, (49)
where χ ∼ 10−8.4(T6n100)
−1/2 is the value of (δF/F )shot for
V = 0, T = 6T6 hr, and norb = 100n100. The corresponding dis-
tance is given by log[Dd,shot/10pc] = 0.2(Vd − MV ) if Vd <
Vmax, and is Dd,shot = Dmax when Vd >Vmax. We take the max-
imum detectable distance to be Dd = min{Dd,int,Dd,shot}.
At this point the simplest approach is to assume that the
parameters {M,Mp,Porb, I} are statistically independent and
carry out a Monte Carlo integration to obtain E . To this
end, we draw M from the Kroupa et al. (1993) initial mass
TABLE 3
NUMBER OF Kepler STARS WITH DETECTABLE
ELLIPSOIDAL OSCILLATIONS
(S/N)min = 1 (S/N)min = 3 (S/N)min = 5
y x = 1 2 x = 1 2 x = 1 2
1 240 166 76 35 33 13
0 99 62 26 12 11 4
-1 33 19 7 3 2 1
function in the range of 0.5–1.5M⊙. The planetary mass is
chosen from the distribution f (Mp) ∝ M
p for Mp = 1–10MJ.
Marcy et al. (2005) find that x ≃ 1 when considering all de-
tected planets; the shape of f (Mp) is not well constrained at
Porb < 10 days. We let x = 1 and 2. We adopt f (Porb) ∝ P
over 1–10 days. Multiplying the resulting value of E by 1000
provides a crude estimate of the actual number of Kepler tar-
gets with detectable ellipsoidal variability. No single value of
y is consistent with the data, and so we consider the reason-
able range y = −1, 0, and +1. Inclinations are chosen under
the assumption that the orbits are randomly oriented, such that
f (cos I) = 1/2 for I ∈ (0,π). Our calculations use fixed values
of (δF/F )int = 10
−5 and T6 = n100 = 1.
Results of our Monte Carlo integrations are shown in Ta-
ble 3 as actual numbers of Kepler targets. The largest number
of detectable systems is obtained when x = y = 1, parameters
that yield the largest proportions short periods and massive
planets. We expect that ∼10–100 Kepler stars may exhibit
ellipsoidal oscillations with S/N & 5. A handful of systems
might have S/N & 10. Higher harmonics from the ℓ = 3 com-
ponents of eq. (16) or modest eccentricities might be accessi-
ble for at most a few stars.
Our integrations also check for cases where the planet is
transiting. As (S/N)min increases from 1 to 5, the fraction
of systems in Table 3 with |cos I| < (R + Rp)/a runs from
≃30% to ≃50%, with a weak dependence on x and y. Such
significant fractions stand to reason, since systems with the
shortest periods have the highest ellipsoidal amplitudes and
transit probabilities. Transit measurements directly give Porb,
sin I & 0.95 (for Porb & 1 day), and (Rp/R)
2. The planet mass
Mp can be determined with the addition of spectroscopic ra-
ELLIPSOIDAL OSCILLATIONS 13
dial velocity measurements, which should be possible for
most of the Kepler targets with detectable ellipsoidal oscil-
lations. The ellipsoidal amplitude then depends on the un-
measured stellar mass and radius via ε∝ R3/M2 (eq. [1]), as
well as the stellar photospheric conditions (eq. [36]). If M
and R are obtained from stellar models, ellipsoidal variabil-
ity may provide an interesting consistency check on all the
system parameters, as well as test the theory of forced stellar
oscillations.
As a last point, we emphasize that stars of mass &1.4M⊙
may have typical ellipsoidal amplitudes of ∼10ε. However,
such stars will also be younger than most Kepler targets and
probably have intrinsic variability ≫10−5. We carried out
Monte Carlo integrations with M = 1.4–1.6M⊙, (δF/F )ell =
10εsin2 I, and x = y = 1. As we vary (δF/F )int from 10
to 10−4, E decreases from large values of ≃0.4 to a small
fraction of ≃0.03 for (S/N)min = 10. Unfortunately, we do
not know how many such stars will be included in the Kepler
target list. Also, there has not yet been a discovery of a giant
planet with Porb < 10 days around a star of mass ≥1.4M⊙, but
exoplanet surveys tend to exclude these more massive stars.
We thank Tim Brown for general discussions and address-
ing Kepler questions, Jørgen Christensen-Dalsgaard for guid-
ance on stellar luminosity perturbations, and Mike Muno for
advice on signal processing. This work was supported by NSF
grant PHY05-51164.
APPENDIX
OSCILLATION EQUATIONS
Here we list the nonadiabatic, linearized fluid equations that we solve numerically. The reader is referred to Unno et al.
(1989) for a complete discussion. Scalar and vector quantities are expanded in spherical harmonics Yℓm and poloidal vector
harmonics, respectively. The momentum, mass, and energy equations are written in terms of the dimensionless variables y1 = ξr/r,
y2 = (δp/ρ + δϕ)/gr, y3 = δϕ/gr, y4 = g−1dδϕ/dr, y5 = ∆s/Cp, and y6 = ∆L/L. Here L is the total (radiative plus convective)
luminosity. The radial flux perturbation is ∆F/F = ∆L/L − 2ξr/r. In radiative zones, the nonadiabatic equations are
d lnr
gk2hr
− y5ρs +
U , (A1)
d lnr
ω2 − N2
1 − η+
−ρsy5 −
, (A2)
d lnr
= y3 (1 − η)+ y4 , (A3)
d lnr
= y1η
+ y2η
ℓ(ℓ+ 1) − η
− y4η + y5ρsη , (A4)
d lnr
+ 4(∇−∇ad) + c2
(∇ad −∇)
gk2hr
c2 + y4
∇ad + y5
∇ (4 −κs) − y6
dU/dr
, (A5)
d lnr
= y1ℓ(ℓ+ 1)
− y2ℓ(ℓ+ 1)
+ y3ℓ(ℓ+ 1)
4πr3ρCpT
ℓ(ℓ+ 1)
, (A6)
where cs is the sound speed, η = d lnMr/d lnr, c2 = (r/Hp)∇(κad −4∇ad)+∇ad(d ln∇ad/d lnr+r/Hp), and we have ignored energy
generation terms. Note that the tidal acceleration −∇U has been added to the momentum equations. In convection zones, we
ignore turbulent viscosity effects and replace the radiative diffusion equation (eq. [A5]) with the prescription ∆s = constant (see
§ 5.1), or more precisely
= 0 . (A7)
Equation (A6) still involves the total (convective plus radiative) luminosity. We ignore energy generation and horizontal flux
perturbation terms, i.e. we ignore all terms with spherical harmonic index ℓ in eq. (A6) in convection zones.
At the center of the star, we require the solutions to be finite, and also set ∆s = 0. At the surface, we set δp = ρgξr and we
require δϕ to decrease outward. This boundary condition is only approximate, as g-modes may propagate above the convection
zone for wave periods of &4days in our 1M⊙ model. The final surface boundary condition is given by eq. (31). Care must be
used in the radiative zone just below the photosphere, since the entropy perturbation is far from the quasi-adiabatic value. If we
solve the radiative diffusion equation in this region, we find that the entropy increases by ∼10 orders of magnitude in just a few
grid points. However, we regard this behavior as unphysical, because the region at the top of the convection zone is optically thin.
To eliminate this unphysical behavior, we set ∆s to a constant at such low optical depths.
REFERENCES
Agol, E., Steffen, J., Sari, R., & Clarkson, W. 2005, MNRAS, 359, 567
Basri, G., Borucki, W. J., & Koch, D. 2005, New Astronomy Review, 49,
Borucki, W., et al. 2004, ESA SP-538: Stellar Structure and Habitable Planet
Finding, 177
Brickhill, A. J. 1983, MNRAS, 204, 537
14 PFAHL, ARRAS, & PAXTON
Brickhill, A. J. 1990, MNRAS, 246, 510
Brown, T. M., Charbonneau, D., Gilliland, R. L., Noyes, R. W., & Burrows,
A. 2001, ApJ, 552, 699
Claret, A. 2000, A&A, 363, 1081
Drake, A. J. 2003, ApJ, 589, 1020
Dziembowski, W. 1977, Acta Astronomica, 27, 203
Goldreich, P., & Nicholson, P. D. 1989, ApJ, 342, 1079
Goldreich, P., & Wu, Y. 1999, ApJ, 511, 904
Goldreich, P., & Wu, Y. 1999, ApJ, 523, 805
Grether, D., & Lineweaver, C. H. 2006, ApJ, 640, 1051
Hansen, C. J., & Kawaler, S. D. 1994, Stellar Interiors: Principles, Structure,
and Evolution, ( New York: Springer)
Henry, T. J., & McCarthy, D. W., Jr. 1993, AJ, 106, 773
Holman, M. J., & Murray, N. W. 2005, Science, 307, 1288
Jenkins, J. M., & Doyle, L. R. 2003, ApJ, 595, 429
Kippenhahn, R., & Weigert, A. 1990, Stellar Structure and Evolution
(Berlin: Springer)
Koch, D., et al. 2006, Ap&SS, 304, 391
Kopal, Z. 1942, ApJ, 96, 20
Kroupa, P., Tout, C. A., & Gilmore, G. 1993, MNRAS, 262, 545
Loeb, A., & Gaudi, B. S. 2003, ApJ, 588, L117
Marcy, G., Butler, R. P., Fischer, D., Vogt, S., Wright, J. T., Tinney, C. G., &
Jones, H. R. A. 2005, Prog. Theor. Phys. Supp., 158, 24
Mihalas, D. 1970, Stellar Atmospheres (San Francisco: Freeman)
Miralda-Escudé, J. 2002, ApJ, 564, 1019
Pace, G., & Pasquini, L. 2004, A&A, 426, 1021
Paxton, B. 2004, PASP, 116, 699
Rasio, F. A., Tout, C. A., Lubow, S. H., & Livio, M. 1996, ApJ, 470, 1187
Robinson, E. L., Kepler, S. O., & Nather, R. E. 1982, ApJ, 259, 219
Rowe, J. F., et al. 2006, ApJ, 646, 1241
Sartoretti, P., & Schneider, J. 1999, A&AS, 134, 553
Seager, S., Whitney, B. A., & Sasselov, D. D. 2000, ApJ, 540, 504
Sirko, E., & Paczyński, B. 2003, ApJ, 592, 1217
Skumanich, A. 1972, ApJ,171, 565
Soszynski, I., Udalski, A., Kubiak, M., Szymanski, M. K., Pietrzynski, G.,
Zebrun, K., Wyrzykowski, O. S. L., & Dziembowski, W. A. 2004, Acta
Astronomica, 54, 347
Sudarsky, D., Burrows, A., & Pinto, P. 2000, ApJ, 538, 885
Udalski, A., et al. 2002, Acta Astron., 52, 1
Unno, W., Osaki, Y., Ando, H., Saio, H., & Shibahashi, H. 1989, Nonradial
Oscillations of Stars, Second Edition (Tokyo: Univ. Tokyo Press)
von Zeipel, H. 1924, MNRAS, 84, 665
Wilson, R. E. 1994, PASP, 106, 921
Zahn, J.-P. 1975, A&A, 41, 329
