© The Authors 2024

1 Introduction

In ultra-hot-Jupiter systems, the tidal interaction between the star and the planet results in the deformation of both bodies into ellipsoidal shapes. The deformation of the star leads to a periodic modulation of the observed flux as the stellar tide, raised by the planet, rotates with the planet’s orbital phase (Russell 1912; Morris 1985; Mislis et al. 2012; Esteves et al. 2013). This periodic signal in the light curve is referred to as stellar ellipsoidal variation, a phenomenon that has been empirically confirmed in several hot-Jupiter systems (see, e.g., Shporer et al. 2019; Wong et al. 2021; Parviainen et al. 2022). Similarly, the stellar tidal force deforms the planet along the axis pointing toward the star. The tidal deformation of the planet is challenging to measure during transit as the elongated axis is predominantly aligned with the line of sight, thereby resulting in a nearly circular sky projection of the planet. Nonetheless, slight deviations from the standard spherical planet transit light curve can be observed as the shape and cross-sectional area projected by the deformed planet varies with phase, thereby allowing the detection of deformation in high-precision light curves (Akinsanmi et al. 2019; Hellard et al. 2019; Barros et al. 2022). Since the projected area of the deformed planet varies throughout its orbit around the host star, it must also impact the observed orbital phase curve.

Although many effects are accounted for, including stellar ellipsoidal variation, when analyzing phase curve observations (see, e.g., Morris 1985; Loeb & Gaudi 2003; Esteves et al. 2013; Shporer 2017), the planet shape is usually assumed to be spherical due to the expected low amplitude of tidal deformation (Bourrier et al. 2020). However, with the precision of present and upcoming instruments, the subtle effects of deformation may now be detectable in phase curves. Since the deformation of a planet depends on its interior structure, its detection will offer valuable insight into the composition of the planet (Kellermann et al. 2018; Hellard et al. 2019; Akinsanmi et al. 2019). Furthermore, at high levels of precision, disregarding the effect of planetary deformation may bias interpretations of phase curves of ultra-hot Jupiters (UHJs).

In this paper, we drop the assumption of sphericity, but instead model the deformed planet shape as a triaxial ellipsoid to investigate the contribution of planet tidal deformation to phase curve observations of ultra-hot Jupiters. In Sect. 2 we model the shape of a deformed planet and describe a numerical and analytical approach to account for tidal deformation in phase curves of UHJs. In Sect. 3 we characterize the signature of tidal deformation across the full planetary orbit. We also compare the detectability of deformation in transit and phase curve observations for several space-based observatories. In Sect. 4 we illustrate how the measurement of tidal deformation can be used to constrain the core mass fraction of a planet, while in Sect. 5 we investigate biases in the parameters derived from phase curve analysis when sphericity is assumed for a deformed planet. We summarize our main conclusions in Sect. 6.

2 Modeling tidal deformation

2.1 Deformation and Love number

The deformation of a planet in response to stellar tidal forces depends on its interior structure and can be parameterized by the second-degree fluid Love number for radial displacement h_2f. In hydrostatic equilibrium, the giant planet deforms as a fluid body, and we have that h_2f = 1 + k_2f, where k_2f is the second-degree Love number for potential that measures the distribution of mass within the planet¹ (Love 1911; Kellermann et al. 2018). Since the Love number h_2f depends on the distribution of mass within the planet, it measures the degree of central condensation of a body, thereby allowing us to constrain the interior structure (Kramm et al. 2011, 2012). The value of h_2f ranges from 0 to 2.5, although certain nonequilibrium processes in the interior and nonlinear perturbative effects can result in values outside this range (Ragozzine & Wolf 2009; Wahl et al. 2021). A body with a fully homogeneous mass distribution (constant density) will have the maximum h_2f value of 2.5, whereas a highly differentiated body with most of its mass condensed in a core will have a lower value (h_2f ⋍ 1 for main sequence stars; Ragozzine & Wolf 2009). More massive objects tend to have smaller h_2f values since they are more compressible, and thus more centrally condensed (Leconte et al. 2011a). Measurements of h_2f can therefore be used to infer the presence and mass of a planetary core. For instance, Saturn’s lower h_2f of 1.39 (Lainey et al. 2017) is indicative of a higher core mass fraction than Jupiter’s h_2f = 1.565 (Durante et al. 2020). Similarly, limits on the Love number of HAT-P-13 b, obtained from its eccentricity due to its unique orbital configuration with a highly eccentric outer planet, allowed constraints to be placed on the maximum core mass of the planet (Batygin et al. 2009; Kramm et al. 2012; Buhler et al. 2016; Hardy et al. 2017). However, these works derive different limits on the Love number of the planet, highlighting the difficulty in accurately measuring the Love number of exoplanets. Therefore, measuring the shape of a tidally deformed giant planet from its light curve allows h_2f to be estimated, which in turn provides direct constraints on the interior structure of the planet.

2.2 Shape model

The shape of a deformed planet is described here as a triaxial ellipsoid with semi-axes [r₁, r₂, r₃] where r₁ is the planet radius oriented along the star–planet (substellar) axis, r₂ is along the orbital direction (dawn–dusk axis), and r₃ is along the polar axis. The volumetric radius of the ellipsoid is given by R_v = (r₁ r₂ r₃)^1/3. According to Correia (2014), the axes of the ellipsoid are related by r₂ = R_v(1 − 2q/3), r₁ = r₂(1 + 3q), and r₃ = r₂(1 − q), where q is an asymmetry parameter given by $$q=\frac{h_{2f}}{2Q_{\text{M}}}{\left(\frac{R_{\text{v}}}{a}\right)}^3$$(1)

and depends on the Love number h_2f, the planet-to-star mass ratio Q_M = M_p/M_⋆, and the orbital semimajor axis a.

The nonspherical shape of the ellipsoidal planet implies that it projects a varying cross-sectional area as it rotates with orbital phase. The projected area as a function of orbital phase angle (ϕ = 2π × phase) is given (e.g., by Leconte et al. 2011a) as $$A(\varphi )=\pi \sqrt{r_1^2{r}_2^2{\cos }^{2}i+r_3^2{\sin }^{2}i\left(r_1^2{\sin }^2\varphi +r_2^2{\cos }^2\varphi \right)},$$(2)

where i is the orbital inclination of the planet. The projected area of the ellipsoid is plotted in Fig. 1a for a planet with the size, mass, and h_2f of Jupiter orbiting a Sun-like star at a distance of 3 stellar radii. The area is normalized to the projected area at mid-transit. We see that the ellipsoid projects its minimum area at mid-transit and mid-eclipse (phases 0 and 0.5), while the maximum is at quadrature where the projected area is 9% larger than at mid-transit. As shown in previous studies (e.g., Correia 2014; Akinsanmi et al. 2019; Hellard et al. 2019), the varying projected shape and area of the ellipsoidal planet modifies the transit light curve. This allows the measurement of the planet’s shape and its Love number, but requires high-precision transit data. Barros et al. (2022) were only able to detect the deformation of WASP-103 b by combining transit observations from several space-based instruments for improved precision. They obtained a 3σ significant measurement of h_2f. Hellard et al. (2020) obtained a less significant measurement for WASP-121 b using only a few transit observations from HST. Therefore, significantly detecting tidal deformation from transit observations remains quite challenging. We can further extend the previous works by probing for the varying area of the deformed planet in full-orbit phase curve observations, opening the potential for a more robust measurement of h_2f that facilitates improved modeling of planetary interior structures.

Fig. 1 Deformed planet phase curve components. (a) Phase-dependent projected area of an ellipsoid normalized to the projected area at phase zero (mid-transit). (b) Schematic of the projected ellipsoidal planet shape at different phases and the illuminated regions. (c) Atmospheric phase variation of a spherical (solid) and ellipsoidal (dashed) planet calculated analytically using Eq. (11). The cases with (green) and without (black) phase offset are shown. (d) Difference between the ellipsoidal and spherical planet atmospheric phase variations for the cases with and without phase offset. This represents the contribution of deformation to the planetary signal. For comparison, the stellar ellipsoidal variation FEV with a peak-to-peak amplitude of 100 ppm is also shown as the brown dotted curve. The cyan shaded regions in all panels represent the in-transit and in-eclipse phases.

2.3 Phase curve model

The phase curve model is a combination of planetary and stellar signals, which we summarize in this section with a modification introduced to account for planet deformation. The phase curve model is composed of transit F_tra, occultation F_occ, planet atmospheric phase variation F_p, and stellar phase variation F_⋆ signals, and is given as a function of ϕ as $$F(\varphi )=F_{\text{tra}}\times F_{\star }(\varphi )+F_{\text{occ}}\times F_{\text{p}}(\varphi ),$$(3)

where the stellar phase variation is composed of the stellar ellipsoidal variation F_EV (Morris 1985) and Doppler boosting F_DB (Loeb & Gaudi 2003) signals as $$F_{\star }(\varphi )=1+F_{\text{EV}}+F_{\text{DB}}.$$(4)

The components of Eqs. (3) and (4) are described in the following sections.

2.3.1 Transit and occultation

The planet’s transit in front of the star and the occultation (secondary eclipse) behind the star can be modeled by several available tools. Here we use the publicly available ellc transit tool (Maxted 2018), where our transit model F_tra adopts the power-2 limb-darkening law (Hestroffer 1997) with c_LD and α_LD as limb-darkening coefficients (LDCs). The power-2 law has been shown to outperform other two-parameter laws in modeling intensity profiles generated by stellar atmosphere models (Morello et al. 2017; Claret & Southworth 2022). The occultation model, F_occ, is the same as the transit, but without limb darkening and normalized to a depth of unity.

To model the transit signal of a deformed planet, we followed the implementation of Akinsanmi et al. (2019) that incorporates the Correia (2014) ellipsoidal shape model (Sect. 2.2) into ellc. Therefore, our generated transit light curve takes into account the varying shape and cross-sectional area of the ellipsoid during transit. In addition to the usual transit parameters, the ellipsoidal transit model takes h_2f and Q_M as inputs, while replacing the typical spherical planet radius R_p by the volumetric radius R_v (Eq. (1)). Spherical planet transit and occultation signals can also be generated by setting h_2f = 0, which makes R_v=R_p.

2.3.2 Ellipsoidal variation

As previously mentioned, an ellipsoidal variation (EV) signal can be observed in phase curves of UHJs as a result of stellar deformation by a massive planet (Morris 1985; Mislis et al. 2012). Detailed expressions for the flux variations due to EV can be found in Esteves et al. (2013) and Csizmadia et al. (2023). However, the EV signal can be approximated by a cosine function at the first harmonic of the orbital period (Shporer 2017), so it has two peaks (at phases 0.25 and 0.75) during an orbit. It is given as $$F_{\text{EV}}=A_{\text{EV}}(1-\cos 2\varphi ),$$(5)

where the semi-amplitude A_EV depends on the planet-to-star mass ratio Q_M, the scaled semimajor axis a/R_*, and the inclination i as $$A_{\text{EV}}={\alpha }_{\text{EV}}{Q}_{\text{M}}{\left(\frac{R_{*}}{a}\right)}^3{\sin }^{2}i$$(6)

where α_EV is a factor that depends on the linear limb-darkening coefficient and gravity-darkening coefficient.

2.3.3 Doppler beaming

This is an orbital photometric modulation caused by the motion of the star around the system’s center of mass. The motion results in a Doppler shift of stellar light leading to a variation in the amount of photons observed in a particular passband (Loeb & Gaudi 2003; Shporer 2017). The signal is described by a sine curve at the orbital period, and so it has one peak (at phase 0.25) during an orbit. It is given by $$F_{\text{DB}}=A_{\text{DB}}\sin (\varphi ),$$(7)

where the semi-amplitude A_DB is given by $$A_{\text{DB}}={\alpha }_{\text{DB}}4\frac{K}{c}$$(8)

Here K is the RV semi-amplitude, c is the speed of light, and α_DB is a coefficient of order unity that depends on the target’s spectrum in the observed passband. For most planets with sub−km s⁻¹ RV amplitudes, A_DB is very small (<5 ppm), and so Doppler beaming has little effect on the phase curve signal.

2.3.4 Planetary atmospheric phase variation

As a planet orbits its star, the planetary flux emanating from different fractions of the projected hemisphere contributes to the total observed flux from the system, resulting in a variation in the observed flux with phase. The planetary flux is composed of reflected light and thermal emission from the atmosphere, but these components are difficult to disentangle in observations as they produce variations that are similar in shape. However, since planets susceptible to tidal deformation are strongly irradiated, most of the observed flux is expected to be from thermal emission rather than reflected light. Here, we calculate the total planetary phase variation using numerical and analytical methods in order to investigate the contribution of planetary tidal deformation.

Numerical calculation

To numerically simulate the planetary phase variation signal, we adopted a pixelation approach where the planet (ellipsoidal or spherical) is projected on a two-dimensional Cartesian grid (XY). The projected shape of the ellipsoidal planet on the grid is an ellipse with semi-axes, r_x and r_y, that varies with orbital phase and can be derived from Eq. (2). On the other hand, the projected shape of the spherical planet is a circle with radius r_x=r_y, which remains constant with orbital phase.

At each phase a flux map is assigned to the projected planet. Different flux maps can be adopted (see, e.g., Louden & Kreidberg 2018, and references therein). Here we define a flux map given as a function of co-latitude (θ) and longitude (λ) of the visible hemisphere as which simulates the simple sinusoidal phase function typically used to model the total brightness variation of planets (e.g., Cowan et al. 2012; Kreidberg et al. 2018; Shporer et al. 2019; Bell et al. 2019; Wong et al. 2021; Zieba et al. 2022). Here θ varies in the range 0–180° between the poles, while λ varies from λ_c − 90° to λ_c + 90° with λ_c denoting the central longitude of the planet at a particular orbital phase (see Fig. 2). The longitudes are defined such that the substellar point is on λ = 0° and the antistellar point is on λ = 180°. The hotspot offset is denoted by δ. The coefficients A₀ and A₁ define the flux value of pixels at the dawn–dusk limbs of the planet, and the flux difference between these limbs and the center at eclipse, respectively. For illustration, we set A₀ = 0 to simulate a total nightside flux of zero, while we set A₁ = 1 so as to have a flux value of 1 at the substellar point. The planet’s phase variation F_p(ϕ) can then be obtained by integrating the flux map across the projected surface at each phase. $$ℱ(\theta ,\lambda )=A_0+A_1\sin \theta \cos (\lambda +\delta ),$$(9)

Assuming a planet with the size, mass, and h_2f of Jupiter with a/R_⋆ = 3, Fig. 2 shows the flux map for the ellipsoidal and spherical planet cases at transit, quadrature, and eclipse phases. We note the larger projected area of the ellipsoidal planet at quadrature where there are more brightly emitting pixels compared to the spherical planet. The generated phase variation signals are also shown, where we see that the varying projected area of the ellipsoidal planet leads to additional flux between transit and eclipse since the emitting surface is larger at these phases. The maximum difference between the ellipsoidal and spherical planet occurs slightly after quadrature where the observer sees the largest illuminated or emitting surface of the ellipsoidal planet. Meanwhile, by construction, the spherical planet signal is equivalent to a cosine function. Numerically computing the planetary contribution to the phase curve can be computationally intensive, especially when fitting to observations. Therefore, we consider an analytical solution.

Fig. 2 Schematic of the numerical calculation of planetary phase variation for an ellipsoidal and a spherical planet. Top: planet flux map for an ellipsoidal planet at transit (left), quadrature (middle), and eclipse (right) phases. The projected area of the planet changes with the number of pixels at each phase given as npix. Middle: same as top panels, but for a spherical planet. Bottom: planet flux variation obtained from integrating the flux map at the different phases. The black and red circles indicate the integrated flux value at each illustrated phase for the ellipsoidal and spherical planets, respectively.

Analytical calculation

Following previous phase curve studies (e.g., Cowan et al. 2012; Shporer et al. 2019; Wong et al. 2021), we analytically modeled the total atmospheric phase variation with a cosine function between the minimum flux, F_min, and maximum flux F_max, as $$F_{\text{p}}(\varphi )=\left(F_{\max }-F_{\min }\right)\frac{1-\cos (\varphi +\delta )}{2}+F_{\min }.$$(10)

The dayside flux F_d (i.e., occultation depth) and nightside flux F_n of the planet are derived as the value of F_p(ϕ) at ϕ = π and ϕ = 0, respectively, while the semi-amplitude of the atmospheric variation A_atm is (F_max − F_min)/2.

As seen in the numerical simulation, the main effect introduced by tidal deformation in phase curves is geometric owing to the changing projected area of the ellipsoid. Therefore, we account for tidal deformation by multiplying the planet’s phase variation (Eq. (10)) by the normalized phase-dependent projected area of the ellipsoid² (Eq. (2)): $$F_{\text{p}}{(\varphi )}^{\text{def}}=F_{\text{p}}(\varphi )\frac{A(\varphi )}{A(0)}$$(11)

We note that other phase variation models (e.g., the kinematic temperature model of Zhang & Showman 2017, or other analytical functions, as given in Csizmadia et al. 2023) can be adopted for F_p(ϕ) and the same operation can be performed to incorporate tidal deformation. The analytical model produces an equivalent planetary phase variation signal to that obtained numerically (see Fig. A.1), thereby validating the simplified model that we use hereafter.

Figure 1c shows the atmospheric phase variations with zero phase offset for a spherical planet (green solid curve) and an ellipsoidal planet (green dashed curve) both with F_d = 3600 ppm³. The phase variation signals with eastward hotspot offset of 20° (black curves) are also shown. Figure 1d shows the difference between the phase variations of the ellipsoidal and spherical planets (i.e., F_p(ϕ)^deſ − F_p(ϕ)), highlighting the contribution of tidal deformation. In the case without hotspot offset (green curve), the contribution is symmetric about the mid-eclipse phase, peaking between quadrature and eclipse on both sides with an amplitude of ~200 ppm. In the case of a hotspot offset (black curve), the deformation contribution is asymmetric with some of its power shifted to the peak before eclipse giving a higher amplitude of ~250 ppm before eclipse and 150 ppm after eclipse. The plot also shows the stellar ellipsoidal variation signal F_EV with an arbitrarily chosen peak-to-peak amplitude of 100 ppm. The EV signal peaks at quadrature (i.e., spaced by 0.5 in phase) making it distinguishable from the planetary deformation signal, which has peaks at different phases with shorter spacing between them.

Our total phase curve implementation accounts for deformation by replacing F_tra, F_occ, and F_p(ϕ) in Eq. (3) by their ellipsoidal planet equivalents⁴. Figure 3 shows the ellipsoidal planet phase curve for different values of h_2f with h_2f = 0 corresponding to a spherical planet. We see that the transit depth varies for different values of h_2f since larger h_2f implies more deformation, which causes the ellipsoidal planet to project a smaller cross section of its shape during transit. However, larger values of h_2f increase the projected area of the ellipsoid toward quadrature, thereby causing an increase in the flux received from the planet between transit and occultation.

Fig. 3 Ellipsoidal planet phase curves with different values of h2f, where h2f = 0 is equivalent to a spherical planet. The inset zooms in on the phase curve region after transit.

3 Signature and detectability of tidal deformation

Tidal deformation is only significant for UHJs orbiting close to the Roche limit of their stars, which limits the detectability to only a handful of confirmed planets. Previous studies (e.g., Correia 2014; Akinsanmi et al. 2019; Hellard et al. 2019; Berardo & de Wit 2022a) have identified WASP-12b, WASP-103 b, and WASP-121 b as some of the most promising targets to measure deformation. For the analyses presented below, we adopt the parameters of WASP-12b as representative of such UHJs.

3.1 The signature of deformation

The detectability of deformation requires that the observed light curve of a deformed planet is distinct such that it cannot be accurately reproduced by a spherical planet model. From visual inspection of Fig. 3, one can imagine that most of the differences between the spherical and deformed planet phase curves can perhaps be explained by a modification of some of the spherical planet model parameters (e.g., radius, hotspot offset). However, in actual observations, the value of these parameters would be initially unknown, but determined by performing a fit to the data. Therefore, the detectable signature of deformation is the residuals obtained from fitting a deformed planet observation (transit or phase curve) with a spherical planet model. The amplitude and timescale of the residual anomalies indicate the prominence of deformation in the data and hence its detectability.

To estimate and compare the signature of deformation in transit and phase curve observations, we assume broadband observations of WASP-12 b in the optical (e.g., with TESS, CHEOPS, or PLATO) and band-integrated observations in the near-infrared (NIR; e.g., with JWST or HST). We compared the deformation signatures by simulating ideal transit and phase curve observations (but no uncertainties) of a deformed planet and performing a least-squares fit⁵ using a spherical planet model. We selected TESS for the representative simulation of the optical observations based on the available measurements and JWST/NIRSpec in PRISM mode for the NIR. Although we used TESS and NIRSpec-PRISM as examples of observing instruments, the approach is applicable to other instruments observing in the optical and NIR. The adopted parameters of the simulation are given in Table 1; they are all allowed to vary freely except for the LDCs and A_DB , for which we used Gaussian priors. In the transit case, the out-of-transit baseline is kept short to minimize the effect of orbital phase variation on the fit, and we put strong priors on the parameters (F_d, F_n, δ, A_EV) since they cannot be accurately estimated from the transit fit. For real observations, the phase variation in transit-only observations can be accounted for either by applying priors on these parameters from previous phase curve analyses and observations of the target or by detrending using a time-dependent function.

The results of the fits are shown in Fig. 4. The residuals of the transit fits display anomalies that represent the signature of deformation in transit light curves as described in Correia (2014) and Akinsanmi et al. (2019): transit ingress and egress anomalies due to the longer transit of the ellipsoidal planet, and in-transit anomalies due to varying projected shape and eclipsed area by the rotating ellipsoidal planet. The residuals of the phase curve fits show an additional anomaly across the orbit also due to the varying projected area of the ellipsoidal planet, which contributes additional flux at different phases that cannot be easily accounted for by the spherical planet model. There are also ingress and egress anomalies in the occultation due to the longer occultation duration of the ellipsoidal planet.

To quantify the deviation of the spherical model fits from the deformed planet simulated observations, we computed the root sum of squares (RSS) of the residuals. We find a greater deviation in the PRISM transit residuals than for TESS indicating a larger deformation signal in the NIR than the optical due to the reduced effect of limb darkening in the NIR that amplifies the ingress and egress anomalies. The phase curve residuals also exhibit larger deviations than the corresponding transit residuals (>2× for PRISM) due to the wider phase range of the phase curve deformation signature. This implies that the effect of deformation is more prominent in phase curves than in transit-only observations. Comparing the phase curve deformation signature between the instruments, we see that the out-of-transit amplitude is just as high as the in-transit amplitude for PRISM, whereas it is much smaller for TESS. The RSS of the TESS phase curve signature does not provide much improvement compared to the RSS of the transit due to the small contribution from the out-of-transit signature.

In summary, phase curves provide a longer phase range than transits, which allows the shape of the planet to be probed better as it rotates with phase. The longer phase range can also be more easily sampled to attain better precision, which facilitates the detection of the anomaly. We quantify the detectability of deformation from transits and phase curves in the next section.

Fig. 4 Fitting the phase curve and transit of a deformed planet with a spherical planet model. Left: simulated deformed planet phase curve of WASP-12 b in the TESS (blue) and NIRSpec-PRISM (red) bands. The spherical planet model fit is overplotted for each light curve. The bottom panel shows the residuals of the fit for each light curve. Right: same as the left panel, but for observation of the transit only. The residuals are shown in the bottom panel for each light curve with the root sum of squares (RSS) quoted.

Fig. 5 Detectability of deformation in the optical (TESS; left) and NIR (JWST/NIRSpec-PRISM; right). The horizontal dashed line indicates the simulated value of h2f =1.565. The red points show the retrieved value of h2f from phase curves and the 1σ uncertainty while the blue points show the retrieved h2f from transits. The gray dotted error bars show the 3σ uncertainty of the measurements. A similar plot for HST/WFC3 is shown in Fig. A.2.

3.2 Detectability of tidal deformation in transit and phase curve observations

Detecting tidal deformation implies measuring a statistically significant value of h_2f from fitting transit and/or phase curve observations with the ellipsoidal planet model described in Sect. 2. We investigated the detectability of tidal deformation at optical and NIR passbands by performing injection-retrieval simulations. We simulated the transit and phase curve of deformed WASP- 12b at a cadence of two minutes using the parameters listed in Table 1. For ease of fitting, we did not include any phase variation in the simulated transit signals, assuming that this will be accounted for as previously mentioned.

In order to investigate the photometric precision required to detect deformation and assess the detectability by different instruments, we added random Gaussian noise of different levels to the simulated transit and phase curve, which we then fit with ellipsoidal and spherical planet models. The dynesty Python package (Speagle 2020) was used to sample the parameter space and derive parameter posteriors. All parameters for each model were allowed to vary freely except for the power-2 LDCs, for which we used Gaussian priors with arbitrary uncertainties of 0.01, Q_M whose prior was based on derived values in Collins et al. (2017), and A_DB whose prior was derived from Eq. (8).

For each noise level we generated 20 simulated observations by adding random realizations of the noise. We then performed the model fit to each observation. The posterior distributions from the fits with the different noise realizations were then merged. This method marginalizes over several possible noise realizations, and thus mitigates the impact of obtaining biased results based on a single randomly generated favorable or unfavorable noise realization.

The detectability of tidal deformation in the optical and NIR is summarized in Fig. 5. It shows the retrieved value of h_2f from the fits to the simulated optical and NIR observations at different noise levels. For the optical simulations, the transit fits allow h_2f to be recovered with 3σ significance for noise levels up to 70 ppm/2 min, while the phase curve fits extends this to a higher noise level of 85 ppm/2 min. At noise levels ≤50 ppm/2 min the phase curve fits provide up to 15% improvement in the precision of retrieved h_2f compared to the transit fits. However, the improvement decreases below 5% at higher noise levels since the relatively small amplitude of the deformation signature outside of transit phases becomes indiscernible from the noise. For the NIR simulations, h_2f is recovered at 3σ in the transit fits for noise levels up to 120 ppm/2 min. Across all added noise levels, the phase curve fits provide a significant improvement (76% at 50 ppm and >30% at 150 ppm) in the precision of retrieved h_2f compared to the transit fits due to the long phase range and large amplitude deformation signature outside transit.

Our injection and retrieval simulations thus imply that detecting tidal deformation is more favorable at NIR passbands compared to the optical since we are able to measure h_2f significantly at higher noise levels in the NIR, particularly for phase curve observations. As phase curves require more observing time than transits, we compared h_2f retrieval from phase curves to the same observing time in transit observations. This allowed us to determine the most effective observing strategy for the precise measurement of h_2f . For the same observing time, we find that phase curves provide improved h_2f precision over the combined transit by up to 10% in the optical and up to 50% for JWST.

3.2.1 Impact of limb darkening on detectability

The effect of limb darkening is most significant at ingress and egress, where the signature of tidal deformation also manifests strongly. limb darkening is therefore capable of compensating for part of the deformation-induced signal, thus reducing detectability, particularly for transit observations (Hellard et al. 2019; Akinsanmi et al. 2019). Therefore, detecting deformation requires proper treatment of limb darkening.

We find in our injection–retrieval simulations that increasing the width of our limb-darkening priors results in reduced detection significance. For instance, when the widths of the LDC priors are slightly increased from 0.01 to 0.015 in the optical simulations, the noise level required for a 3σ detection reduces from 70 to 50 ppm/2 min in the transit fits and from 85 to 55 ppm/2 min in the phase curve fits. The reduction in both cases occurs because most of the detection power in the optical comes from the transit phases where limb darkening acts. However, the increased LDC prior width in the NIR causes only a slight reduction in noise level from 120 to 110 ppm/2 min for a 3σ transit detection due to the weaker limb-darkening effect in the infrared. The NIR phase curve fits remain largely unchanged, even with a further increase in the LDC prior width, since the deformation signature spans a wide phase range with a large out-of-transit amplitude that allows a precise constraint of the deformation in spite of limb darkening.

3.2.2 Disentangling deformation from atmospheric phenomena

It is possible for other atmospheric effects to confound the deformation signal or make its detection challenging. For example, highly irradiated planets may exhibit strong day-to-night temperature gradients that can result in different chemical compositions and scale heights between the two hemispheres (see, e.g., Pluriel et al. 2020; Helling et al. 2020). This leads to limb asymmetries (Espinoza & Jones 2021) that can affect the detectability of deformation in transit.

To investigate whether the effects of day-to-night inhomogeneities can mimic tidal deformation, we simulated the transit of a planet with limb asymmetries using catwoman (Jones & Espinoza 2020) and fit it with a deformed planet model. We find that tidal deformation is unable to account for limb asymmetries since tidal deformation does not cause asymmetric transits. Conversely, we also found that catwoman is unable to accurately fit a deformed planet transit indicating that the two effects are distinguishable.

Furthermore, the effect of dayside–nightside inhomogeneities is chromatic and is dependent on the compositions of the two regions. In contrast, the Love number is an intrinsic property of the planet that remains constant across wavelengths. Therefore, we can differentiate between them by ensuring that consistent h_2f is derived across different spectral bins or passbands. Therefore, spectroscopic observations with HST and JWST will be very useful in this regard.

3.2.3 Detectability of deformation by different instruments

Figure 5 can be used to estimate the detectability of deformation by different instruments observing phase curves or transits of close-in giant planets in the optical and NIR. Here, we estimate the detectability for observations with JWST, PLATO, TESS, CHEOPS, and HST.

JWST. Using the JWST Exposure Time Calculator,⁶ we estimated a precision of 40 – 50 ppm/2 min for observations of WASP-12 b with different JWST instruments, for example NIRSpec-PRISM. At the similar noise level of 50 ppm/2 min, our NIR simulations recovered a highly significant detection of $h_{2f}={1.52}_{-0.087}^{+0.094}$ (17σ) from the phase curve fit and $h_{2f}={1.55}_{-0.37}^{+0.32}$ (4.2σ) from the transit fit. Figure 6 shows the ellipsoidal and spherical planet model fits to the simulated PRISM observations with one realization of the 50 ppm noise. The residual of the spherical planet fit shows a clear modulation induced by tidal deformation that cannot be accounted for by the spherical planet. Figure 6 further shows the contribution of planet deformation (cf. Fig. 1d) to the total observed phase curve in comparison to the best-fit stellar ellipsoidal variation signals F_EV . We observe that the amplitude of F_EV in the spherical planet model fit increases (from the simulated value of 60 ppm to 210 ppm) in an effort to account for the deformation signal of the planet. However, it is unable to do so since the deformation contribution peaks at different phases with a shorter inter-peak interval compared to F_EV. This is still the case when a tidal phase lag parameter is added to F_EV. This implies that the assumption of sphericity for a deformed planet can lead to biases in the estimation of other system parameters. In Sect. 5 we further investigate the bias on different system parameters when assuming sphericity for a deformed planet.

PLATO. For stars brighter than V ≃ 11, PLATO will be capable of achieving a photometric precision of 34 ppm h⁻¹ (Rauer et al. 2014; Rauer & Heras 2018) equivalent to 186ppm/2 min. The observation strategy of PLATO consists of a long-pointing phase where two different fields are observed for two to three years each, and a short observation phase (known as step and stare), where each stare would last two to three months (Rauer & Heras 2018). Assuming a conservative scenario where WASP-12 is observed in the short observation mode for two months, this will amount to ~50 orbits of the planet leading to a precision of ~30 ppm/2 min. At this noise level in our optical simulations, we recovered a 4.3σ detection with $h_{2f}={1.59}_{-0.37}^{+0.40}$ from the phase curve fit and a 3.6σ detection with h_2f = 1.54 ± 0.42 from the transit fit. Figure A.4 shows the best-fit models to the simulated PLATO observations. PLATO long-pointing observations of 2–3 yr of UHJs susceptible to deformation will provide exquisite photometric precision that will only be limited by stellar noise and limb-darkening modeling.

TESS. We estimated the precision of TESS from its actual observations of WASP-12 b for which it obtained a precision of ~ 1800 ppm/2 min (Wong et al. 2022). A total of 78 orbits of the planet were observed (across four TESS sectors), which improved the precision to ~200 ppm/2 min. At this noise level in our optical simulations, we obtained a ~2.3σ detection of h_2f from the phase curve fit and a slightly lower ~2σ detection from the transit fit. Similar to Fig. 6, Fig. A.3 shows the best-fit models to the simulated TESS observations with 200 ppm/2 min noise. In order to detect deformation with higher significance, more TESS observations will be required to lower the noise level.

CHEOPS. Using the CHEOPS Exposure Time Calculator⁷, we estimated a photometric precision of 450 ppm/2 min⁸ for observations of WASP-12 b. For a month-long observation of WASP-12 b (or 28 accumulated phase curves), a precision of 85 ppm/2 min can be attained. At this noise level in our optical simulations, we recovered a 3σ detection of h_2f from the phase curve fit and 2.5σ from the transit fit. Figure A.5 shows the best-fit models to the simulated CHEOPS observations.

HST. Observations of WASP-12 have been performed with different HST instruments. Using the G141 grism of HST/WFC3 Kreidberg et al. (2015) obtained a precision equivalent to ~110ppm/2min on three transit visits. Phase curve and transit simulations at this noise level result in a 4.3σ and 2.6σ detection of h_2f, respectively. We see that a noise level of 110 ppm for WFC3 constrains h_2f better than 85 ppm from CHEOPS due to the larger NIR phase curve amplitude (Fig. A.6). Doubling the number of visits to six improves the detection significance in the phase curve and transit to 5.6σ and 3σ, respectively.

We note that gaps within the data (particularly for HST and CHEOPS) will reduce the deformation detection significance estimated for the different instruments or increase the number of visits required to reach a certain detection significance. For example, even after combining three HST/WFC3-G141 transit visits of WASP-12 b in Kreidberg et al. (2015), the ingress and egress phases still have poor coverage, making it challenging to measure tidal deformation from these datasets within a reasonable number of HST visits. However, observations can be scheduled so that the combined visits attain better phase coverage. Instrumental and/or astrophysical noise (e.g., spacecraft jitter or stellar activity) can also impact the detection of deformation from the light curves. These will have to be mitigated while analyzing the data using methods such as polynomial decorrelations, Gaussian processes (Barros et al. 2012), or wavelets (Csizmadia et al. 2023).

Fig. 6 Spherical planet (SP) and ellipsoidal planet (EP) model fits to the simulated deformed planet signals with a PRISM noise level of 50 ppm/2 min. The top panel shows the simulated phase curve (left) and transit (right) observations with the fitted models overplotted. The phase curve fit also shows the contribution of deformation to the light curve (blue) and the FEV component of the EP and SP model fits. The middle and bottom panels show the residuals from fitting a spherical planet model (cyan) and an ellipsoidal planet model to the simulations.

4 Constraining the core mass from h2f

As mentioned previously, measuring the Love number of a planet sheds insight into its interior structure. For instance, Baumeister et al. (2020) found that adding the Love number as input to interior structure models significantly decreases the degeneracy of the possible interior structure. Furthermore, simulations of the interior structure of low-mass planets by Baumeister & Tosi (2023) found that a Love number precision of 10σ can constrain the core and mantle size of an Earth analog to ~13% of its true value. Here we show how the retrieved value of h_2f from the transit and phase curve fits of the different instruments can be used to place constraints on the core mass fraction M_c /M of the planet. To do this, we calculate the theoretical value of h_2f as a function of M_c /M using Eqs. (1) and (2) of Buhler et al. (2016) and adopting the radial density profile in Eq. (7) of Helled et al. (2014), which consists of a constant density core and a stratified density envelope.

Figure 7a shows the dependence of h_2f on M_c /M. As M_c /M increases the planet is more centrally condensed and less homogeneous, which causes h_2f to decrease. We overplot the retrieved values of h_2f from the transit and phase curve fits of the JWSTNIRSpec-PRISM simulations (Sect. 3.2.3). We see from the 1σ lower limits of h_2f that the phase curve fit constrains the maximum core mass fraction of the planet as 0.19, whereas the larger h_2f uncertainties from the transit fit give a maximum M_c /M of 0.56. Figure 7b compares the maximum core mass fraction derived from the retrieved h_2f from the phase curve and transit fit of the different instruments. We see that the JWST phase curve fit provides the tightest constraint on the core mass of the planet followed by the HST and PLATO phase curve fits, and then the JWST and PLATO transit fits. We see that the <3σ h_2f measurements from CHEOPS and TESS provides no constraint on the core mass with a maximum core mass being equal to the planetary mass. We note that although the measurements from the different instruments are still consistent with a zero core mass, more precise h_2f measurements (e.g., from more than one JWST phase curve or PLATO long-pointing) will also allow the lower limit of M_c/M to be determined, definitively constraining the planetary core and improving interior structure models.

5 Bias on derived system parameters

Parameters of the model adopted in fitting observations are often correlated with one another, making it possible for other phase curve parameters to dampen or entirely mask the deformation signal depending on the noise level of the light curve. Assuming sphericity when fitting the phase curve of a deformed planet could introduce biases in the estimation of various other parameters. The deformation then manifests itself as an astrophysical source of systematic error in the measurement of the parameters.

We investigate the parameter biases by fitting a spherical planet model to the simulated TESS, PLATO, CHEOPS, HST, and PRISM deformed planet phase curves with increasing values of h_2f up to 1.8. The noise level of each instrument estimated in Sect. 3.2.3 was added to the simulated light curves. For each parameter we compare the obtained posterior to the true simulated value by calculating the number of standard deviations by which the median posterior value differs from the true value: (fit – truth)/σ_fit. This normalization allows us to easily compare the parameter deviation for each instrument.

The result of the fits is shown in Fig. 8. As expected, for h_2f = 0 (spherical planet) the fitted parameters are close to (<1σ) the true simulated values, leading to deviations close to zero. A zero deviation is also obtained for cases where the deformed planet phase curve is indistinguishable from that of a spherical planet, due to high noise level in the light curve compared to the deformation signal.

For all instruments we observe significant biases in the estimation of R_p/R_★, a/R_★, and b, which all show clear trends with h_2f. Since a tidally deformed planet projects only a small cross section of its shape during transit, the best-fit spherical planet model has a smaller R_p than the volumetric radius R_v of the ellipsoid. Therefore, as h_2f increases and the planet is more deformed, it projects an even smaller cross-sectional area, which requires a smaller spherical R_p to match the transit depth. We find that the derived R_p/R_★ increasingly deviates from (i.e., underestimates) the simulated R_v/R_★ as a function of h_2f. At the highest simulated h_2f value of 1.8, R_p/R_★ deviates by up to 25σ for the TESS fit; 50σ for the CHEOPS fit; and >80σ for the HST, PLATO, and PRISM phase curve fits, which have the most precise radius measurements. When combined with the mass estimate from radial velocity (RV) measurements, the radius deviation can thus lead to an overestimation of the density of the planet by more than 10% if sphericity is assumed (Burton et al. 2014; Correia 2014; Barros et al. 2022). As the precision of RV mass estimates increases, the bias in derived densities due to tidal deformation will constitute a bottleneck in constraining the composition of exoplanets susceptible to deformation (Berardo & De Wit 2022b).

Both a/R_★ and b deviate at h_2f =1.8 by >2.5σ for TESS, CHEOPS, and HST, while the higher precision phase curves of PLATO and PRISM deviate by > 12σ. These parameters deviate from the true inputs in order to compensate for the longer transit duration of the deformed planet (Arkhypov et al. 2021). The limb-darkening parameters c_LD and α_LD also show biases as these parameters attempt to absorb the deformation signal concentrated around ingress and egress.

Biases are also observed in the dayside and nightside fluxes, F_d and F_n, which show opposing trends as a function of h_2f. Both parameters act together to increase the amplitude A_atm of the atmospheric phase variation in order to account for the additional flux from the deformed planet at phases with larger projected areas compared to the sphere (see Fig. 6). Biases on these parameters affect the temperatures derived for the hemispheres of the planet. The temperature of the dayside is overestimated and that of the nightside is underestimated, which results in a larger day-to-night contrast. The deviations in F_d and F_n are up to 20σ at h_2f =1.8 for the PRISM fit where the amplitudes of these parameters are much higher than the light curve scatter. Converting the dayside and nightside fluxes to temperatures⁹ (T_d and T_n), we find that their deviations lead to a 5.4σ over-estimation of T_d and a 13σ underestimation of T_n for PRISM, and ≲3σ deviations for HST. For the optical phase curves with lower simulated values of F_d = 466 ppm and F_n = 0 ppm, deviations in the derived F_d decrease to ~8σ (T_d deviation of 2.2σ) for PLATO and <2σ for TESS and CHEOPS (T_d deviations of < 1σ), whereas deviations in F_n are close to zero.

The hotspot offset δ also exhibits a positive deviation of up to 11σ at h_2f =1.8 for the PRISM phase curves with the largest simulated offset of 20°. The larger eastward hotspot of the spherical planet attempts to compensate for the deformation-induced flux before eclipse. However, the optical simulations with small offsets of 6.2° showed only <1.5σ deviations. As mentioned in Sect. 3.2.3, the amplitude of ellipsoidal variation amplitude increases in an attempt to compensate for the deformation contribution. This leads to a positive deviation from the true value by up to 50σ at h_2f = 1.8 for the PRISM fit, ~l0σ for the HST and PLATO fits, 5σ for the CHEOPS fit, and 2σ for the TESS fit. There are no significant deviations in the derived transit time and Doppler beaming amplitude.

In general, the PRISM phase curve shows the largest parameter deviations due to the high precision of JWST data and the large amplitude deformation signal that the spherical planet model parameters try to compensate for. Conversely, the TESS phase curve shows the smallest parameter deviations due to the high noise level of the simulated observations, which makes the phase curve similar to that of a spherical planet and makes detecting deformation challenging. Nevertheless, we see that even for such low or nondetection cases, deformation can still lead to significant deviations (ranging from 3 to 25σ) in some of the derived parameters. Therefore, we recommend accounting for planetary deformation when fitting the phase curves of UHJs by multiplying the planet’s atmospheric phase variation by the phase-dependent projected area of the ellipsoid (Eq. (2)). This introduces a single new parameter h_2f to fit, if the mass ratio is kept fixed to RV values, which allows marginalizing over the possible shapes of the planet. In analyzing the HST and Spitzer phase curves of WASP-103 b, Kreidberg (2018) accounted for the planet’s deformation by similarly including the phase-dependent projected ellipsoid area in their model. However, instead of fitting for shape parameters of the ellipsoid, the planet shape was fixed using tabulated predictions in Leconte et al. (2011a) assuming a planetary radius and age.

Fig. 7 Constraints on the core mass fraction of WASP-12 b derived from the retrieved h2f. Left: Constraints from JWST phase curve (red) and transit (blue). The black curve shows the theoretical Love number as a function of core mass fraction for a two-layer model using the formulations of Buhler et al. (2016) and Helled et al. (2014). Right: Upper limits on core mass fraction based on the precision of retrieved h2f from different instruments. The number of phase curve and transits vists used are JWST −1 visit, HST −3 visits, PLATO −50 visits, CHEOPS −28 visits, TESS −78 visits.

Fig. 8 Parameter biases from fitting a spherical planet model to simulated phase curve of deformed planets of different h2f. The parameter names are given in Table 1. For each parameter, we compare the deviation of the fitted spherical planet parameters (fit) from the true input parameters of the deformed planet (truth) and scale this by the standard deviation of the fit (σfit). The colors show the biases for the different instruments.

6 Conclusions

In this work we investigated the effects of tidal deformation on the phase curves of UHJs orbiting close to the Roche limit of their stars. We expanded on previous works (e.g., Leconte et al. 2011a; Burton et al. 2014; Correia 2014; Akinsanmi et al. 2019; Hellard et al. 2019) studying tidal deformation in transit observations by introducing a simple approach to additionally account for tidal deformation when modeling phase curve observations of UHJs. Similar to Correia (2014) and Akinsanmi et al. (2019), our model describes the planet as a triaxial ellipsoid parameterized by the Love number h_2f, thereby allowing us to estimate its value from fitting observations. In addition to the signature of deformation in transit caused by the varying projected shape and area of the planet, we showed, using numerical and analytical calculations, that the planet’s deformation also manifests itself outside of transit for the same reason. Therefore, the extended phase coverage of the deformation signal in phase curves allows the planet’s shape to be probed better.

Using injection–retrieval simulations, we find that detecting tidal deformation is more favorable at NIR passbands compared to the optical as we are able to measure h_2f significantly at higher noise levels in the NIR. Furthermore, NIR phase curve observations significantly improve the deformation detection compared to transit-only observations due to larger planetary phase curve amplitudes in the NIR. For the precision expected from several JWST instruments and modes, we find that the measurement of h_2f for WASP-12 b improves from a precision of ~4σ from one transit-only analysis to 17σ from a single phase curve. Combining several transits equivalent to the phase curve time investment still does not attain the h_2f measurement significance of a phase curve.

The high h_2f precision obtainable from JWST phase curves of relevant targets will provide unprecedented constraints on the planet’s interior structure (e.g., core mass fraction Buhler et al. 2016; Kramm et al. 2012). Additionally, spectrophotometric light curves from HST and JWST will allow us to verify that the derived h_2f is not wavelength dependent, effectively disentangling tidal deformation from atmospheric phenomena that are chromatic.

Our results also show that detecting deformation from phase curves is relatively unaffected by the uncertainty in the limb-darkening profile as long as there is a high-amplitude deformation signature outside the transit phases, as obtainable for NIR observations. Finally, we showed that the assumption of sphericity when analyzing the phase curve of UHJs can lead to biases in several derived parameters, thereby significantly limiting the ability to properly characterize these targets (e.g., their radii, densities, dayside and nightside temperatures). We thus recommend accounting for deformation by modifying the planet’s phase variation by the phase-dependent projected area of the ellipsoid, even for cases where the deformation is barely or not detectable. Atmospheric studies of UHJs through transmission spectroscopy could also be affected by tidal deformation. Lendl et al. (2017) suggested that planet deformation could be the cause of the enhanced features measured in the transmission spectra of WASP-103 b, but this remains to be understood and can be investigated using our ellipsoidal planet model.

Acknowledgements

This work has been carried out within the framework of the National Centre of Competence in Research PlanetS supported by the Swiss National Science Foundation under grants 51NF40_182901 and 51NF40_205606. M.L. and B.A. acknowledge the financial support of the Swiss National Science Foundation under grant number PCEFP2_194576. S.C.C.B. acknowledges support from FCT -Fundação para a Ciência e a Tecnologia through national funds and by FEDER through COMPETE2020 – Programa Operacional Competividade e Internacionalização by these grants: UIDB/04434/2020; UIDP/04434/2020, 2022.06962.PTDC.

Appendix A Additional figures

Fig. A.1 Comparison of the numeric and analytic computation of the planetary phase variation signal. The spherical planet signal is equivalent to a cosine function.

Fig. A.2 Same as Fig. 5, but for HST/WFC3/ simulations.

Fig. A.3 Same as Fig. 6, but for simulated 200 ppm/2 min noise observations for TESS.

Fig. A.4 Same as Fig. 6, but for simulated 30 ppm/2 min noise observations of PLATO.

Fig. A.5 Same as Fig. 6, but for simulated 85 ppm/2 min noise observations for CHEOPS.

Fig. A.6 Same as Fig. 6, but for simulated 110 ppm/2 min noise observations for HST/WFCЗ-G141.

References

All Tables

All Figures

Fig. 1 Deformed planet phase curve components. (a) Phase-dependent projected area of an ellipsoid normalized to the projected area at phase zero (mid-transit). (b) Schematic of the projected ellipsoidal planet shape at different phases and the illuminated regions. (c) Atmospheric phase variation of a spherical (solid) and ellipsoidal (dashed) planet calculated analytically using Eq. (11). The cases with (green) and without (black) phase offset are shown. (d) Difference between the ellipsoidal and spherical planet atmospheric phase variations for the cases with and without phase offset. This represents the contribution of deformation to the planetary signal. For comparison, the stellar ellipsoidal variation FEV with a peak-to-peak amplitude of 100 ppm is also shown as the brown dotted curve. The cyan shaded regions in all panels represent the in-transit and in-eclipse phases.

In the text

Fig. 2 Schematic of the numerical calculation of planetary phase variation for an ellipsoidal and a spherical planet. Top: planet flux map for an ellipsoidal planet at transit (left), quadrature (middle), and eclipse (right) phases. The projected area of the planet changes with the number of pixels at each phase given as npix. Middle: same as top panels, but for a spherical planet. Bottom: planet flux variation obtained from integrating the flux map at the different phases. The black and red circles indicate the integrated flux value at each illustrated phase for the ellipsoidal and spherical planets, respectively.

In the text

	Fig. 3 Ellipsoidal planet phase curves with different values of h2f, where h2f = 0 is equivalent to a spherical planet. The inset zooms in on the phase curve region after transit.
In the text

Fig. 4 Fitting the phase curve and transit of a deformed planet with a spherical planet model. Left: simulated deformed planet phase curve of WASP-12 b in the TESS (blue) and NIRSpec-PRISM (red) bands. The spherical planet model fit is overplotted for each light curve. The bottom panel shows the residuals of the fit for each light curve. Right: same as the left panel, but for observation of the transit only. The residuals are shown in the bottom panel for each light curve with the root sum of squares (RSS) quoted.

In the text

Fig. 5 Detectability of deformation in the optical (TESS; left) and NIR (JWST/NIRSpec-PRISM; right). The horizontal dashed line indicates the simulated value of h2f =1.565. The red points show the retrieved value of h2f from phase curves and the 1σ uncertainty while the blue points show the retrieved h2f from transits. The gray dotted error bars show the 3σ uncertainty of the measurements. A similar plot for HST/WFC3 is shown in Fig. A.2.

In the text

Fig. 6 Spherical planet (SP) and ellipsoidal planet (EP) model fits to the simulated deformed planet signals with a PRISM noise level of 50 ppm/2 min. The top panel shows the simulated phase curve (left) and transit (right) observations with the fitted models overplotted. The phase curve fit also shows the contribution of deformation to the light curve (blue) and the FEV component of the EP and SP model fits. The middle and bottom panels show the residuals from fitting a spherical planet model (cyan) and an ellipsoidal planet model to the simulations.

In the text

Fig. 7 Constraints on the core mass fraction of WASP-12 b derived from the retrieved h2f. Left: Constraints from JWST phase curve (red) and transit (blue). The black curve shows the theoretical Love number as a function of core mass fraction for a two-layer model using the formulations of Buhler et al. (2016) and Helled et al. (2014). Right: Upper limits on core mass fraction based on the precision of retrieved h2f from different instruments. The number of phase curve and transits vists used are JWST −1 visit, HST −3 visits, PLATO −50 visits, CHEOPS −28 visits, TESS −78 visits.

In the text

Fig. 8 Parameter biases from fitting a spherical planet model to simulated phase curve of deformed planets of different h2f. The parameter names are given in Table 1. For each parameter, we compare the deviation of the fitted spherical planet parameters (fit) from the true input parameters of the deformed planet (truth) and scale this by the standard deviation of the fit (σfit). The colors show the biases for the different instruments.

In the text

	Fig. A.1 Comparison of the numeric and analytic computation of the planetary phase variation signal. The spherical planet signal is equivalent to a cosine function.
In the text

	Fig. A.2 Same as Fig. 5, but for HST/WFC3/ simulations.
In the text

	Fig. A.3 Same as Fig. 6, but for simulated 200 ppm/2 min noise observations for TESS.
In the text

	Fig. A.4 Same as Fig. 6, but for simulated 30 ppm/2 min noise observations of PLATO.
In the text

	Fig. A.5 Same as Fig. 6, but for simulated 85 ppm/2 min noise observations for CHEOPS.
In the text

	Fig. A.6 Same as Fig. 6, but for simulated 110 ppm/2 min noise observations for HST/WFCЗ-G141.
In the text