Issue 
A&A
Volume 675, July 2023



Article Number  A57  
Number of page(s)  11  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/202346082  
Published online  06 July 2023 
Tidally heated exomoons around ϵ Eridani b: Observability and prospects for characterization
^{1}
Leiden Observatory, Leiden University,
PO Box 9513,
2300 RA
Leiden, The Netherlands
email: kleisioti@mail.strw.leidenuniv.nl
^{2}
Faculty of Aerospace Engineering, TU Delft,
Building 62 Kluyverweg 1,
2629 HS
Delft, The Netherlands
^{3}
Lunar and Planetary Laboratory, University of Arizona,
Tucson, AZ
85721, USA
Received:
6
February
2023
Accepted:
2
May
2023
Context. Exomoons are expected to orbit gas giant exoplanets just as moons orbit Solar System planets. Tidal heating is present in Solar System satellites, and it can heat up their interior, depending on their orbital and interior properties.
Aims. We aim to identify a tidally heated exomoon’s (THEM) orbital parameter space that would make it observable in infrared wavelengths with MIRI/JWST around ϵ Eridani b. We study the possible constraints on orbital eccentricity and interior properties that a successful THEM detection in infrared wavelengths can bring. We also investigate what exomoon properties need to be independently known in order to place these constraints.
Methods. We used a coupled thermaltidal model to find stable equilibrium points between the tidally produced heat and the heat transported within a moon. For the latter, we considered a spherical and radially symmetric satellite with heat being transported via magma advection in a sublayer of melt (asthenosphere) and convection in the lower mantle. We incorporated uncertainties in the interior and tidal model parameters to assess the fraction of simulated moons that would be observable with MIRI.
Results. We find that a 2R_{Io} THEM orbiting ϵ Eridani b with an eccentricity of 0.02 would need to have a semimajor axis of 4 planetary Roche radii for 100% of the simulations to produce an observable moon. These values are comparable with the orbital properties of the satellites of the Solar System gas giants. We placed similar constraints for eccentricities up to 0.1. We conclude that if the semimajor axis and radius of the moon are known (e.g., with exomoon transits), tidal dissipation can constrain the orbital eccentricity and interior properties of the satellite, such as the presence of melt and the thickness of the meltcontaining sublayer.
Key words: planets and satellites: interiors / planets and satellites: detection / infrared: planetary systems / planets and satellites: individual: ϵ Eridani b
© The Authors 2023
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1 Introduction
Since all giant planets in the Solar System host moon systems, it is plausible that exomoons orbit exoplanets as well. Even though detection techniques have shown tremendous success in identifying over 5000 exoplanets, exomoon detection is more challenging because the deviation on the exoplanet’s signal caused by their presence is relatively small. However, progress in observational astronomy has made the detection of exomoons a possibility for the near future (Lazzoni et al. 2022; Teachey & Kipping 2018; Heller 2017; Heller et al. 2016).
Several indirect and direct methods have been proposed for the detection of exomoons. Indirect detection methods study the effect that a moon has on a planet’s signal. These methods include Transit Timing Variations (TTV; Kipping 2009; Simon et al. 2007; Sartoretti & Schneider 1999), Transit Duration Variations (TDV; Kipping 2009), and other photometric effects on the stellar signal (Heller 2014), centroid shifts (Agol et al. 2015), and Doppler monitoring of directly imaged exoplanets (Ruffio et al. 2023; Vanderburg et al. 2018) caused by the gravitational interaction between the planet and its companion. Other methods for exomoon detection include dynamical sculpting of the circumplanetary disk (Kenworthy & Mamajek 2015) and the detection of planetsatellite mutual eclipses (Cabrera & Schneider 2007).
Though there have been several tentative detections (see e.g., Kipping et al. 2022; Oza et al. 2019; Teachey & Kipping 2018; Kenworthy & Mamajek 2015) and despite the numerous known exoplanets and the surveys searching for companions around them (Kipping et al. 2022, 2012), no exomoon detection has been confirmed yet. Apart from being due to their small signal, this may also be due to the biases of detection methods to detect exoplanets closer to their host star, meaning that any potential satellite would be stripped off during the planet’s migration (see e.g., Dobos et al. 2021; Trani et al. 2020).
Direct imaging, which is sensitive to young exoplanets further out from their host star, could offer an observational window to far out systems. Only recently, Benisty et al. (2021) directly imaged the thermal emission from a circumplanetary disk at submillimeter wavelengths in the PDS 70 system, a strong prerequisite for exomoon formation. Direct imaging of tidally heated exomoons (THEMs) is challenging because, except for the need to suppress the star’s light, the final signal contains both the flux from the moon and the planet. However, THEMs (Limbach & Turner 2013) can be brighter than their host planet’s signal in some infrared bands, due to tidal interactions between the exomoon and the host planet. A direct detection via thermal spectral energy distribution (SED), in contrast to other proposed methods such as transits, does not rely on a specific satellite orbital phase nor require multiple observations.
This method offers a detection window for exomoons around gas giant planets at distances of several tens of astronomical units from the parent star. The method can also directly measure the flux emitted from the surface, which under the assumption of thermal equilibrium can be equal to the tidal heat flux produced in the interior. Even though the mass of the satellite is not measurable via this method, the advantage of constraining the tidally generated flux can offer possibilities for characterization because the latter depends on the interior structure and orbital properties of the moon.
Tidal models with a broad range of complexity have been used to model Solar System bodies and moons (RoviraNavarro et al. 2022; Steinke et al. 2020; Bierson & Nimmo 2016; CastilloRogez et al. 2011; Hussmann & Spohn 2004; Fischer & Spohn 1990; Reynolds et al. 1987) and have already been applied to exoplanets (Barr et al. 2018; Shoji & Kurita 2014; Henning et al. 2009). For the first time, Dobos & Turner (2015) modeled THEMs while taking into account viscoelasticity, a characteristic that describes the tidal response of the material making up the interior of planetary bodies. RoviraNavarro et al. (2021) used a coupled thermalorbital evolution model and compared the effect that different rheological models have on THEM surface temperatures while studying their longevity. Finally, Jäger & Szabó (2021) took into account the presence of a hotspot on the surface of homogeneous THEMs without melt and modeled the increased radiation in thermal wavelengths to infer detectability.
In this work, we study the effect of tidal heating on the direct imaging observations of exoIos. We define an exoIo as an exomoon with the same density and interior structure as the moon Io, that is, with a thin, meltcontaining asthenosphere beneath the lithosphere (Moore 2001). We use models developed to describe Io’s tidal heat flux (RoviraNavarro et al. 2021; Moore 2003) and apply them to THEMs with the assumption of thermal equilibrium. We use ϵ Eridani b as a test case for a potential host planet. The exoplanet ϵ Eridani b is a closeby gas giant with a distance of 3.2 parsec (MacGregor et al. 2015) and a minimum mass of 0.65–0.78 M_{J} (Rosenthal et al. 2021; LlopSayson et al. 2021; Mawet et al. 2019), making it a plausible host to relatively large satellites (Canup & Ward 2006). In Sect. 2, we explain the tidalthermal model and the interior structure used in this work as well as the algorithm that solves for the surface temperature of moons in thermal equilibrium. In Sect. 3, we describe how we define an observable moon around ϵ Eridani b. In Sect. 4, we present the orbitalinterior configurations that would lead a THEM around ϵ Eridani b to be detectable with MIRI under the assumptions of our model. We then extend this analysis through a sensitivity analysis of the various free parameters of our interior model, which allows us to quantify the interior constraints a detection may provide. Finally, in Sect. 5, we discuss possible ways to validate our analysis; factors that could affect the interpretation of our results, such as the presence of moon atmospheres and different moon sizes; and exomoon constraints that our method can yield in combination with other exomoon detection methods.
2 Models
To evaluate a THEM’s detectability and infer its surface heat flux, we used an interior model that describes the moon’s layers and heat transfer mechanisms and a tidal model that calculates the tidal dissipation given the interior structure and orbital parameters. A schematic overview of how these two models interact is given in Fig. 1.
We modeled THEMs using models initially developed to describe Io’s interior, which have been shown to be compatible with Io heat flux observations (Spencer et al. 2000; Veeder et al. 1994). In our model, we assumed a metallic core, a solid silicate mantle, and a partially molten sublayer of melt beneath a stiff lithosphere. The moons are assumed to be spherically symmetric, and the properties of each layer are uniform. We assumed that tidal heat is generated in the deep mantle and asthenosphere (viscoelastic layers) and that the lithosphere behaves elastically. Our model implementation was adapted from RoviraNavarro et al. (2021). In the following subsections, we introduce both the interior and tidal models we used.
Fig. 1 Feedback between the heat transfer and the tidal models. The interior structure affects the Im(k_{2}) value, which the tidal dissipation depends on. High tidal dissipation can cause the formation of melt, leading to a new sublayer (new interior structure; see the black arrows in the upper part of the plot). Thermodynamic equilibrium is reached when the moon’s interior transfers the same amount of heat as what is generated via tidal interactions (see the grey arrows in the lower part of the plot). 
2.1 Interior model
We assumed that different heat transfer mechanisms are present in the interior layers of the exoIos (see Fig. 2). The heat that is generated through tidal friction in the viscoelastic layers (mantle, asthenosphere) is ultimately emitted via radiation from the surface. We assumed that the heat is transported via convection in the lower mantle and via conduction in the stiffer lithosphere, compatible with bodies in a stagnant lid regime (e.g., Reese et al. 1999; Schubert et al. 1979).
For high tidal heating rates, the interior temperature rises and a partially molten sublayer, the asthenosphere, forms beneath the moon’s lithosphere. As introduced by Moore (2001), we modeled the heat transport in the asthenosphere with melt segregation. This model successfully describes Io’s current state. If no melt segregation is assumed, convection alone cannot explain the observed heat flux of Io when assuming thermal equilibrium (Moore 2003). As such, the exoIos discussed in our work are modeled with the same mechanism.
Overall, the parameters of our interior model are shown in Fig. 3. We assumed a core of radius R_{c} equal to 0.52 × R (Anderson et al. 1996), where R is the moon radius, with density ρ_{c}; a mantle with density ρ_{m}; a shear modulus µ_{m}; solidus viscosity η_{s,m}; and a temperature T_{m}. The dependence of the mantle viscosity (η_{m}) on the temperature was parameterized with the activation energy E_{a} (Eq. (1a)). The asthenosphere’s viscosity (η_{α}) also depends on the presence of melt through the B parameter, as shown in Eq. ((1b); Mei et al. 2002). Finally, the melt segregation is a function of the permeability exponent n and the scale velocity γ. (1a) (1b)
where R_{g} is the ideal gas constant, ϕ_{α} is the melt fraction of the asthenosphere, and T_{s} is the average of the mantle solidus temperature. The input values used can be found in Table 1. Table 2 shows the parameters included in our sensitivity analysis (see Sect. 4.2).
The interior model output parameters, which were obtained once thermal equilibrium was reached (also listed in Table 3), consist of the asthenosphere thickness (h_{α}) and the melt fraction needed to transfer heat via melt segregation for a value of T_{m}. The material properties thatmake up the moon’s interiorare temperature dependent. Thus, all the output parameters listed in Table 3 and Fig. 3 are a function of T_{m}. RoviraNavarro et al. (2021) explain the interior properties’ dependency on the melt fraction and temperature in detail, and they include a detailed description of the heat transfer model (including convection and conduction). We highlight that when the mantle temperature crosses a threshold, which corresponds to the fraction of melt in the asthenosphere rising above 0.45 (Moore 2003), the rheology of the mantle is no longer adequately described using viscoelasticity and is more accurately portrayed as a magma ocean, a regime that is not modeled in our work. We only discuss the relevance of a magma ocean regime qualitatively in terms of detectability.
Fig. 2 Interior layers and the corresponding heat transfer mechanisms assumed in our model. 
Fig. 3 Input and output parameters used in the thermal model. Descriptions of these values can be found in Tables 1, 3 and 2. 
Interior parameters – inputs.
2.2 Tidal model
The amount of tidal heat dissipated in the interior of a moon with zero obliquity is given through (e.g., Segatz et al. 1988; Makarov & Efroimsky 2014): (2)
where n is the mean motion, e the orbital eccentricity of the satellite, and G the universal gravitational constant. We assumed synchronous rotating moons and performed simulations for eccentricity values up to 0.1. Forlargereccentricities, highorder terms (𝒪(e^{4})) must be taken into account since Eq. (2) is accurate until the second eccentricity order (Renaud et al. 2021). We assumed that the moon has zero obliquity, close to Io’s (Baland et al. 2012), and thus we did not take it into account in the calculation of the total amount of tidal heat generated.
The imaginary part of the k_{2} Love number, Im(k_{2}) describes the deformation of the moon’s gravity field due to the tidal interactions with the planet. It depends on the interior structure, properties, rheology of the moon, and the moon’s orbital period (Segatz et al. 1988). Efforts to obtain Im(k_{2}) for Io (Lainey et al. 2009) with astrometric observations constrained it at 0.015 ± 0.003.
Equation (2) demonstrates how the tidally generated heat flux of a moon depends on the orbital properties of the satellite. Selfluminous exomoons closer to the host planet and with higher orbital eccentricity would produce higher tidal heat fluxes, which is favorable to their observability.
2.2.1 Different approaches to Im(k_{2}) calculation
A simplistic approach to calculate the imaginary k_{2} Love number assumes that Im(k_{2}) = const., where Q (tidal quality factor) quantifies the fraction of the orbital energy that is dissipated per orbit due to friction. This parameter is not well constrained, however rocky bodies are expected to have values from 10 to 500 (Goldreich & Soter 1966). Such an approach, does not take into account the feedback between the thermal state of a planet and its response to tidal forces. In addition, it does not consider the tidal quality factor’s dependency on the orbital period of the satellite (Renaud & Henning 2018). Ithas also been found to underestimate tidal heat production (see e.g. Meyer & Wisdom 2007; Dobos & Turner 2015).
Another approach, which has been widely used for Solar System moons (Fischer & Spohn 1990; Hussmann & Spohn 2004), is to model tidal dissipation using the viscoelastic theory for selfgravitating bodies (Peltier 1974; Sabadini R. 2016), where a rheological law is needed to couple stress and strain. This method allows for the incorporation of the material properties’ temperature dependence (e.g. Karato & Wu 1993). The most simple viscoelastic model is the Maxwell model. The value of Im(k_{2}) is closely related to the Maxwell time (). When the tidal period is close to this value, tidal dissipation reaches its maximum. For periods longer than what corresponds to the Maxwell time, the moon responds as a viscous fluid and for shorter tidal periods, as an elastic body. This approach has also been implemented for the study of exoplanet tidal responses (Shoji & Kurita 2014; Henning et al. 2009; Barr et al. 2018).
A viscoelastic approach to describe tidal heating in exomoons was used for the first time by Dobos & Turner (2015). They implemented the Maxwell viscoelastic model, which, however, does not accurately describe the laboratory behavior of olivine (Jackson & Faul 2010) and does not take into account the anelastic transient creep response over timescales shorter than the Maxwell time. The transient creep is described as a regime where the strain rate is a function of time. On the contrary, the (more) advanced Andrade rheological model (Andrade 1910) adopts this behavior and is more realistic. The Andrade model has been used in studies of Solar System bodies (CastilloRogez et al. 2011; Bierson & Nimmo 2016) and exoplanets (Renaud & Henning 2018; Walterová & Běhounková 2020). Renaud & Henning (2018) found significantly higher tidal heating when assuming Andrade rheology for Iolike moons, thus the incorporation of it – instead of the Maxwell one – has implications on THEM detectability.
Parameters and their ranges used for the sensitivity analysis.
Interior parameters – outputs.
2.2.2 Calculation of Im(k_{2})
We used the viscoelastic theory for selfgravitating bodies and incorporated the Andade model to calculate Im(k_{2}). The tidal response of the moon can be obtained by solving the equations of motion for the deformation of each layer in the Fourier domain, via the correspondence principle (Peltier 1974). The latter results in a set of differential equations that describe the deformation of each layer.
To solve the differential equations and obtain the gravitational potential at the surface, we used the propagator matrix technique (JaraOrué & Vermeersen 2011; Sabadini R. 2016) and the Andrade rheological law. For more details on how we solved the propagator matrix technique, see Appendix A of RoviraNavarro et al. (2021). The viscoelastic response of the material in each layer depends on the shear modulus µ and viscosity η. The Fourier transformed shear modulus is related to the creep function , which for the Andrade law is (Efroimsky 2012): (3) (4)
The transient creep response is modeled in the last term of Eq. (4) and is described by two parameters ζ and α. We assumed that ζ = 1 for the rest of our work. For assumptions on α, see Sect. (4.2).
2.3 Thermal equilibrium
The tidal model describes the heat that is generated through tidal interactions, and the heat transfer model simulates how this energy is transferred through the different layers of the interior. The two models are not independent but interact with each other through feedback (Fig. 1, upper part). An Im(k_{2}) value corresponds to a particular interior structure but also affects the amount of tidal dissipation that is produced (Eq. (2)). For high values of tidal dissipation, the interior temperature increases, and consequently the interior properties, layers, and melt fraction change. This means that the interior structure is modified, affecting the Im(k_{2}).
Therefore, equilibrium is reached once the amount of heat emitted from the moon’s surface into space is equal to the amount of heat that is generated (Fig. 1, lower part). Henning et al. (2009) showed that significantly tidally active planetary bodies achieve thermal equilibrium in a few million years. The overall thermal equilibrium equation that we solve for is: (5)
where is the surface heat flow and is the internal heat production. The value of was computed via the equations for melt segregation presented in Moore (2003) and via Eq. (2).
In thermal equilibrium, the average surface temperature T_{surf} can be calculated by StefanBoltzmann law: (6)
where σ is the StefanBoltzmann constant, q_{s} is the equilibrium surface heat flux (where , A is the moon’s bond albedo (assumed equal to Io’s), and S is the stellar irradiation. For ϵ Eridani, this corresponds to 0.34 × L_{Sun} (Saumon et al. 1996), where L_{Sun} is the solar luminosity. The second term of Eq. (6) describes the heat flux that is not reflected by the moon’s surface and is thus taken into account in the thermal budget.
Thermal equilibrium points were defined as the intersections of the tidal heat and the transferred heat (i.e., the heat generated and the heat advected to the surface of the moon). These points can be stable or nonstable (Moore 2003), as shown in Fig. 4. At a stable equilibrium point, a temperature increase would lead to a higher heat advection rate compared to the heat generation one. The moon, thus, tends to restore its equilibrium. The reversed scenario is described by a . In this case, a slight increase in the interior mantle temperature would drive the moon out of equilibrium since the heat production would exceed the advection capabilities of the moon, warming up the interior (Fig. 4). High temperatures coming from formation would potentially make the moon approach Fig. 4 from the right side of the graph. This is a favorable condition for the moon to reach the higher stable equilibrium states first, such as the warmer stable equilibrium point.
In this work, we present results for the stable equilibrium points that correspond to higher temperatures. Scenarios where the moon would reach the lower stable equilibrium point exist (e.g., going out and in a mean motion resonance, see Fuller et al. 2016); however, studying these cases was beyond the scope of this work, which is mainly focused on assessing whether and for what interior and orbital properties putative THEMs could reach temperatures that allow for observations and under what assumptions we can place constraints on their properties.
Fig. 4 Generated tidal heat flux as a function of the mantle temperature (T_{m}) and the log of the asthenosphere’s (η_{α}) and mantle viscosity (η_{m}). Stable and nonstable equilibria points are shown as black and gray dots, correspondingly. At a nonstable equilibrium point, a deviation in the mantle temperature drives the moon out of equilibrium, whereas a body at a stable thermal equilibrium state tends to restore its equilibrium. The tidal heat flux was calculated for a moon with a semimajor axis of a = 6.59 Roche radii, e = 0.02, α = 0.5, and B = 30 with the rest of the sensitivity parameters set to their nominal values (Table 2). 
3 Observability
In this section, we define what we consider to be an observable moon around ϵ Eridani b. Figure 5 shows the expected star, planet, and moon fluxes in this system. The exoplanet ϵ Eridani b was indirectly detected through radial velocity with a mass of ≈0.65M_{J} (Rosenthal et al. 2021; LlopSayson et al. 2021). The planet’s spectrum was approximated via the models for young gas giants of Spiegel & Burrows (2012), assuming no clouds, 500 Myr age, 1 M_{J}, and solar metallicity (ϵ Eridani b: −0.04 dex; Rosenthal et al. 2021). The star’s SED is shown as a black body of temperature 5084 K (Kovtyukh et al. 2003). In some wavelength regions, the flux of the THEM surpasses that of the planet – for a 270 K exomoon, these regions are at 3.6 and 6 µm, and widen significantly at higher exomoon temperatures.
The instrument MIRI is equipped with four coronagraphs: one Lyot and three fourquadrant phase masks (4QPMs; Boccaletti et al. 2022). Of these, the 4QPMs have the smallest achievable inner working angles (IWAs), between 0.33 and 0.49arcsec (Boccaletti et al. 2015). For a semimajor axis of 3.5 AU (LlopSayson et al. 2021), the planet has a separation of ≈1.09 arcsec from its host star, which is outside of the inner working angle of the 4QPM coronagraphs. Boccaletti et al. (2022) measured the onsky performance of MIRI’s coronagraphs and concluded that postprocessing techniques can bring the final contrast down to the background and detector limited noise floor at separations larger than one arcsecond. Future instruments, like METIS (Brandl et al. 2021), are expected to reach even smaller contrast ratios.
We therefore define a “detectable exomoon” as the coldest exomoon that reaches the 10σ point source detection limit for MIRI in a 10000 s integration at ϵ Eridani’s distance of 3.2 parsecs (MacGregor et al. 2015). The average moon surface temperature in our model is calculated through Eq. (6) and represents the theoretical total thermal output of the moon, assuming no localized variations due to volcanic activity or hotspots on the surface. Figure 5 shows that this limit is reached for an average surface temperature of 270 K for a 2R_{Io} exomoon orbiting ϵ Eridani b. Moons of the same size with a higher effective temperature are also considered detectable.
The moon’s thermal flux could be constrained from observed excess of heat flux at different wavelengths where the planet is not expected to be bright, for example, at water or methane absorption bands. It is beyond the scope of this work to disentangle the moon’s signal from the planet’s in specific observational bands. Rather, our goal is to assess the detectability and orbitalinterior configurations that can lead to observations. The use of different planetary atmosphere models could also potentially affect the combined planetmoon flux measured. We ignored any effect the latter would have in disentangling the putative detected flux in its planet and moon components.
Fig. 5 Fluxes of hypothetical 2R_{Io} moons around e Eridani b. The black lines are the 10 σ and 10 000 s detection limits for MIRI (Glasse et al. 2015). The planet model (Spiegel & Burrows 2012) is shown in gray. The red and purple lines correspond to the black bodies of a 270 and 470 K 2 × R_{Io} exomoon, and the yellow dashed line corresponds to the black body of the star. 
Fig. 6 Surface equilibrium temperature of a 2R_{Io} THEM as a function of its orbital eccentricity (e) and semimajor axis (a). The black dot indicates Io’s current orbital parameters. The nominal values of the interior parameters that we used are presented in Table 2. The 270 K isotherm divides the parameter space into detectable and nondetectable exomoons. 
4 Results
In this section, we explore how interior and orbital parameters could be constrained if the thermal heat flux of a THEM was known. We first explore how this could be known if our model was perfect and all the free parameters were well determined. Thus, the first subsection refers to results obtained assuming the nominal values of Table 2. We then take into account modeling parameter uncertainties to explore how constraints can be placed, as shown in Sect. 4.2.
4.1 Interior structure dependence on orbital properties
Figure 6 shows the resulting surface temperatures using the models from Sect. 2 for a 2R_{Io} exoIo in a stable thermal equilibrium state as a function of the moon’s semimajor axis a and orbital eccentricity e. For moons with the same density as Io, this corresponds to 8 M_{Io}. Satellite formation theories indicate an upper limit of M_{moon}/M_{planet} ≈ 10^{−4} (Canup & Ward 2006). According to this limit, exomoons of such mass could be found around 1.1 M_{J} planets or higher. This is compatible with the ϵ Eridani b minimum mass between 0.65 and 1.55 M_{J} (Rosenthal et al. 2021; Hatzes et al. 2000). We focused our analysis on exomoons of the same size because the moon’s radius needs to be constrained in order to derive any interior and orbital conclusions with our method. The way to constrain the moon’s radius, as well as the implications it has on the interpretation of our results, is further discussed in Sects. 5.2 and 5.3.
Tidal interactions between a host exoplanet and the modeled exomoons could heat their surface temperature up to the order of 600 K for the shown ranges of orbital parameters (Fig. 6). Moons experiencing such an amount of tidal heating would be in a magma ocean regime (see below). The values of a and e that would make a SuperIo reach MIRI’s detection limit (Sect. 3) are demonstrated with the 270 K isotherm, which divides Fig. 6 into two areas: “detectable” and “nondetectable” moons. For example, a THEM orbiting epsilon Eridani b at 5.5 Roche radii would need to maintain an eccentricity higher than 0.009 to be detectable with MIRI. This value is of comparable magnitude to Solar System moons in mean motion resonance (MMR).
Figure 6 can thus be used to predict surface equilibrium temperatures of THEMs under the assumptions of our model (see Tables 1,2). The dark red area defines the area of the parameter space for which the melt fraction of the mantle exceeds 0.45 (Moore 2003). When this excess occurs, the asthenosphere behavior resembles a magma ocean. We note that tidal dissipation in such a regime is largely unstudied and depends on poorly constrained parameters. SuperIos with orbital characteristics that fall within this area produce higher tidal heat rates than what their interiors can advect, either through convection or heat piping. As a result, they are in a “magma ocean” state where they reach the limits of our models since liquid tides Tyler (2011) or poroviscoelasticity RoviraNavarro et al. (2022) would need to be considered. Nevertheless, such moons would be detectable (under our assumptions) since their surface would reach higher temperatures than those with lower asthenosphere melt. Because all moons with asthenosphere properties on the boundary of a magma ocean are detectable, moons with a magma ocean would be as well. Such moons could reach warmer thermal equilibrium states with a more efficient heat transport mechanism (convection in low viscous layer; Solomatov 2007). They could also still be approaching thermal equilibrium at the moment of detection, radiating excess heat (e.g., from formation).
In Figs. 11a and b, the surface temperature and Im(k_{2}), respectively, of a SuperIo orbiting its host star at 4.65 Roche radii are plotted as a function of its orbital eccentricity. An assumption of a constant Im(k_{2}) would lead to significantly different results compared to the correspondence principle technique. Calculating Im(k_{2}) via the correspondence principle technique, as is done in this work, can affect our interpretation of whether an exomoon is detectable or not. For example, for the studied moon and a semimajor axis of 4.65 Roche radii the equilibrium, Im(k_{2}) can change by a factor of ten, depending on the eccentricity (see Fig. 11b, blue line). From Eq. (6), this difference would correspond to a calculated temperature that is off by 500 K, a temperature difference that, as seen in Sect. 3, can define whether a moon is observable.
The value of Im(k_{2}) provides a direct link between an exomoon observation with its interior structure and could be used to draw conclusions on the interior properties shown in Table 3. Figure 7 shows two different equilibrium interior structures that correspond to eccentricities of 0.004 and 0.02 and a semimajor axis of 4.65 Roche radii. The two different eccentricities have equilibrium surface temperatures of 270 K and 470 K, respectively. The Im(k_{2}) of the 470 K exoIo is approximately twice as high as the 270 K exomoon, for a semimajor axis of 4.65 Roche radii.
Given that moons would reach stable thermal equilibrium at different tidal heat rates for each pair of a and e, figures similar to Fig. 6 can be produced for every interior property in Table 3. Figure 8 shows how the melt fraction of the asthenosphere changes with varying semimajor axis and eccentricity. Jäger & Szabó (2021) modeled THEMs with one hotspot, assuming a nonhomogeneous heat flux distributed among their surface. Such hotspots enhance detectability in thermal wavelengths. While our model does not take this into account, we highlight that the various hotspots of a THEM could be used to place constraints on its orbital period, based on the variability of the signal that the hotspots would induce in direct imaging observations.
If the surface temperature and semimajor axis are known, the moon’s eccentricity (Fig. 6) and interior properties (Fig. 8) can be inferred because a surface temperature value corresponds to one stable equilibrium Im(k_{2}) value for a given eccentricity. In other words, the surface temperature of a putative exomoon with a given semimajor axis and interior structure would be a function of the orbital eccentricity. The blue lines in Fig. 11 show how for a semimajor axis of 4.65 Roche radii the inferred surface temperature (Fig. 11a) and the interior properties (Figs. 11b–f) change with the orbital eccentricity. Since some interior properties exhibit uncertainties, we performed a sensitivity analysis to assess whether and to what extent constraints on the orbital eccentricity and interior can be placed with a putative heat flux observation.
Fig. 7 Two different equilibrium interior structures for a = 4.65 Roche radii. The one on the left (e = 0.004) has a thinner asthenosphere compared to the one on the right (e = 0.02), resulting in a smaller Im(k_{2}). The parameter symbols are explained in Table 3. 
Fig. 8 Asthenosphere melt fraction (ϕ_{α}) of stable thermal equilibrium states of exoIos as a function of their orbital eccentricity and semimajor axis. 
4.2 Orbitai and interior constraints from observations
To summarize Sect. 4.1, our coupled tidalthermal analysis can relate the surface heat flux of putative exoIos with the corresponding equilibrium Im(k_{2}) values for a given moon’s semimajor axis. This also means that constraints on the interior equilibrium states could possibly be deduced. Given a semimajor axis, a single solution for the orbital eccentricity exists when assuming the nominal values of Table 2. However, the values of the parameters in Table 2 are not known exactly, and thus their uncertainty needs to be incorporated in the conclusions that can be drawn with our model.
In this section, we present the results of this analysis when varying the parameters presented in Table 2 for the shown number of points and ranges. We obtained surface equilibrium temperatures as well as the corresponding interior parameters in Table 3 with the goal to assess (i) how a surface temperature alone can constrain the moon semimajor axis and eccentricity; (ii) the fraction of simulations for a given eccentricity and semimajor axis that lead to an observable exoIo; and (iii) the constraints in the interior properties and the orbital eccentricity if the surface temperature and the semimajor axis are known.
For all three purposes, we varied the mantle shear modulus from 2 × 10^{9} (Steinke et al. 2020) to 6.5 × 10^{10} Pa, which corresponds to the end case of a mantle shear modulus equal to the lithospheric one, assuming no melt (Segatz et al. 1988; Peale et al. 1979). The Andrade parameter α can take values from 0.1 to 0.5 for olivine (Gribb & Cooper 1998), while the melt fraction coefficient can take values from 10 to 40 (Henning et al. 2009). Finally, we varied the parameters that control the magma advection between typical values (η : 2–3 Moore 2003, γ : 10^{−6}−10^{−5} Katz 2008). Because the mantle solidus viscosity is highly uncertain, we varied it by a factor of ten from the commonly used value of 10^{16} Pas (Renaud & Henning 2018; Moore 2003; Fischer & Spohn 1990). We note that higher values of η_{s,m} would need to be considered for internal structures with an asthenosphere melt fraction larger than the threshold melt fraction defined in Sect. 2 (Kervazo et al. 2022). Each set of the sensitivity parameters led to different stable thermal equilibrium states and thus surface temperatures. The sensitivity parameters were varied altogether along the discussed ranges. Ultimately, each a and e pair could result in 720 different equilibrium states.
The influence of the sensitivity analysis on the final equilibrium temperature is demonstrated in Fig. 9, where the a and e values that can lead to 270 K (observable with MIRI; Sect. 3) are shown. This already puts constraints on the orbital properties that a particular heat flux observation would correspond to. For example, an exoIo with an eccentricity of 0.01 and a surface temperature observation of 270 K could be within ≈4−7.2 Roche radii away from the host planet in thermal equilibrium.
To assess the fraction of simulations that result in an observable moon, we calculated the equilibrium surface temperature for the sensitivity parameters’ ranges shown in Table 2. Each set of parameters led to a different equilibrium state for the same pair of a, e. If the resulting surface temperature was higher than 270 K, we considered the moon observable (Sect. 3). We then calculated the fraction of the 720 simulated exoIos that are observable for each pair of semimajor axis and eccentricity, and the results are presented in Fig. 10. For an eccentricity of 0.02, an exoIo around ϵ Eridani b would need to have a semimajor axis of ≈4 Roche radii for 100% of the simulations to produce an observable moon. We note that this eccentricity is about five times higher than that of Io, and the semimajor axis is around 1.6 times smaller than Io’s semimajor axis.
Incorporating the uncertainties of the sensitivity parameters places more realistic constraints on the orbital eccentricity and interior properties for a theoretical T_{surf} and a detection. Figure 11 shows the possible solutions for the interior properties of an exoIo orbiting ϵ Eridani b at 4.65 Roche radii as a function of the orbital eccentricity (light blue shading). The blue lines in the figure panels correspond to the nominal values analysis of Sect. 4.1. For an observation of 270 K (black dotted lines) only exoIos with e smaller than 0.022 would be in thermal equilibrium. This constraint is looser for higher temperatures. For 470 K (red dotted lines) the corresponding range would be 0.01−0.1. Moreover, the constraints on the interior properties become stricter for higher temperatures. For example, only bodies with melt fractions larger than 0.1 can explain an observation of 470 K for the given semimajor axis.
As the eccentricity increases the equilibrium surface temperature increases since the moon experiences more tidal heating. As we approached larger eccentricities, the range of possible surface temperatures (and the other interior parameters) tightened. This is because the equilibrium Im(k_{2}) experiences bigger variations close to the Maxwell time (see Fig. 11b) for the considered range of parameters. The peak of the curve, which corresponds to the Maxwell time, can be seen for lower eccentricities of the shown semimajor axis. In the nominal values analysis, the equilibrium structure with viscosity and shear modulus tuned close to the Maxwell time for a = 4.65 Roche radii forms for an eccentricity of 0.015.
Fig. 9 Range of exoIo orbital properties that are compatible with a surface temperature of 270 K and a radius of 2R_{Io} as a result of the sensitivity analysis. 
Fig. 10 Fraction of simulations that lead to an observable moon (warmer than 270 K; Sect. 3) as a function of the moon's orbital semimajor axis and eccentricity. 
5 Discussion
The direct detection of THEMs around exoplanets relies on disentangling the different fluxes (star, planet, moon) in a system at a given separation. The flux of the moon needs to be comparable to (if not higher than) the planet’s flux at a selected wavelength band, and contrast limitations with the star must also be taken into account. In this work, we focus on modeling the moon’s expected thermal flux, which depends on the amount of tidal heating that is produced in the moon interior, which is in turn a function of the interior structure and orbit of the moon.
The expected tidal heat of THEMs was modeled following a more simplistic approach regarding tidal interactions and interior modeling by Limbach & Turner (2013). For the first time, Dobos & Turner (2015) took into account the viscoelastic behavior that describes the tidal response of the material making up the interior of planetary bodies. In their work, Dobos & Turner (2015) included the Maxwell viscoelastic approach commonly used to model the response of Solar System planetary bodies (Shoji et al. 2013) as well as in exoplanet science (Barr et al. 2018). They assumed a twolayered body and convection as the main heat transfer mechanism in the mantle. In our analysis, we used a coupled thermaltidal model combined with viscoelastic Andrade rheology, and we assumed magma advection as the main heat transfer mechanism in a partially molten sublayer.
The effect of changing the classical Maxwell approach to the Andrade one has been studied for planetary bodies (Renaud & Henning 2018; Walterová & Běhounková 2017). In general, Andrade rheology produces higher tidal heating rates, due to the transient creep mechanism. This is aligned with ourresults when compared with Dobos & Turner (2015), who for a 2R_{Io} moon at Io’s orbital period and an eccentricity of 0.1 found 273 K. We found an equilibrium temperature above 400 K for the same orbital properties.
Tides provide a direct link between the interior and the orbital properties of a planetary body. Observations that can constrain the amount of tidal heat in the interior (e.g., thermal flux) combined with accurate modeling of the interior and potentially the orbit can be helpful for drawing inferences regarding the system properties of exoplanets and/or exomoons. Direct observations of THEMs are (the only) direct method of probing the tidally heated surface. Our analysis shows the effect of tidal heating on the surface temperature of THEMs. We found regimes of orbital properties where the moon has to be in a magma ocean state to be able to dissipate the tidally generated heat. As seen in Fig. 11, we also found that we can place constraints on interior properties and eccentricity depending on the (observed) temperature when given a semimajor axis and a moon radius. Several approaches to exomoon detection have been proposed (see Sect. 1) in which the exomoon mass (Kenworthy & Mamajek 2015; Sartoretti & Schneider 1999; Vanderburg et al. 2018; Teachey & Kipping 2018), radius and semimajor axis (Teachey & Kipping 2018), and atmospheric properties (Agol et al. 2015; Oza et al. 2019) could be constrained. Our work provides a way to infer system properties, such as the orbital eccentricity, and compatible interior structures. Finally a THEM detection would indicate the likely existence of at least a second moon in resonance since tides would circularize an isolated moon’s orbit.
In the following sections, we explore some observational and modeling challenges. These include discussions on longevity, semimajor axis and radius constraints, and more massive moons that could potentially host an atmosphere. Finally, we explore how our results could be potentially validated via other observational methods.
Fig. 11 Computed ranges for the model’s output parameters (Table 3) that lead to exoIos with a = 4.65 Roche radii and a radius of 2 × R_{Io} in a stable thermal equilibrium. The blue line corresponds to the nominal values analysis (Sect. 5). The gray and red dotted lines respectively show the limits of a 270 K and 470 K surface temperature observation. 
5.1 Longevity and mean motion resonance
The ϵ Eridani system is a relatively young system (400–800 Myr; Mamajek & Hillenbrand 2008; Janson, Markus et al. 2015), and such systems have a potential disadvantage for the direct detection of THEMs, as any planet in them is expected to be warm from formation. Since the planetmoon contrast also plays a significant role in THEM observability, this could (possibly) hinder a detection.
On the other hand, the system’s young age could offer an observational advantage compared to older systems. In an isolated moonplanet system, tides quickly circularize the moon’s orbit, and close moons can migrate outward quite fast (RoviraNavarro et al. 2021), which limits the observability window. A second moon in resonance can boost the eccentricity of the first moon and maintain the surface temperature high enough for detection via tidal heating. In both Jupiter’s and Saturn’s satellite systems, there are currently MMRs that likely originated with the outward migration of inner satellites, capturing the outer ones in MMRs (Yoder & Peale 1981; Dermott et al. 1988). Even though the MMRs in the Solar System indicate that they are likely prevalent in older systems as well, it is reasonable to assume that MMRs are common in relatively young systems since moons can be captured in MMRs during their formation (Peale & Lee 2002; Ogihara & Ida 2012).
The lifetime of an MMR depends on the mass ratio of the moons, among other parameters, and a shortlived resonance could last <200 Myr (Tokadjian & Piro 2022). RoviraNavarro et al. (2021) calculated that an exoIo with 2R_{Io} could have a temperature higher than 400 K for 10 Myr (during the early stages of an MMR). This would make such a moon detectable for 10 Myr. Thus, if one considers the expected timescale of MMRs, young systems could potentially serve as good targets for THEM observations.
5.2 Combination with other exomoon detection methods: Constraining the moon radius and semimajor axis
As discussed in Fig. 11, which describes the possible constraints we can place on the orbital eccentricity and the interior properties of the moon, we assumed a fixed semimajor axis and satellite radius. This means that tidal heating can provide constraints on the interior properties and eccentricity if the orbital period and the size of the satellite are known. Therefore, the question of how can we constrain these two parameters (and to what accuracy) in a direct imaging observation appears.
Even though both LlopSayson et al. (2021) and Mawet et al. (2019) find inclinations closer to an edgeon orbit of 78.8 and 89.0, respectively (with large uncertainties), there is no consensus yet of ϵ Eridani b’s inclination. Nevertheless, we explored possible ways of constraining the properties in scenarios where a system is closer to being edgeon and faceon. Since a moon in such a system would not be spatially resolved from the planet it orbits, its semimajor axis cannot be constrained in the same way the semimajor axis of a planet is measured in exoplanet direct imaging observations (e.g., Bohn et al. 2019). However, with time series observations of the same target, plausible constraints on the moon’s orbital period could be derived. If the system is edgeon (assuming the planet’s and moon’s orbits are coplanar), the moon would transit the planet, so the task of deriving the moon’s orbital period would be straightforward. For systems close to a faceon orbit, one could constrain a moon’s orbital period via the periodicity of the signal coming from patterns on its surface, as hotspots are expected to be present on the surfaces of THEMs. Time series observations of an edgeon system can also be used to infer the satellite’s radius (Cabrera & Schneider 2007).
5.3 Larger moons
Though exomoons can be of various sizes and some are even more massive than the ones seen in our Solar System (Canup & Ward 2006), the discussed results in this work refer only to 2R_{Io} exomoons. Since tidal dissipation depends on the moon’s radius, we would expect larger satellites to reach warmer thermal equilibrium states, especially in interior regimes with more efficient heat transport than heat piping. However, larger satellites could be massive enough to sustain a substantial atmosphere (Lammer et al. 2014). In this case, direct detection of the thermal emission would be hindered, since the outer atmospheric layer would be relatively colder and detection via spectroscopic signatures would be more promising (Oza et al. 2019). Nonetheless, we point out that a substantial atmosphere would redistribute the tidal energy, making the SED resemble a black body when hotspots on the surface shift the SED to bluer wavelengths (Limbach & Turner 2013).
5.4 Future model validation
In this work we applied models previously used to reproduce Io’s averaged heat flux (Moore 2003; RoviraNavarro et al. 2021) in order to investigate the parameter space, the resulting mean surface temperatures, and the possible detection of THEMs via direct imaging. These models have not been yet validated for bodies larger than Io; thus, we explore possible observations and complementary models that could validate the applicability of our results to the analysis of observations in the future.
The observables that we propose are mostly for quantities that are related to tidally induced active volcanism. One notable option is the spectral signatures from volcanic outgassing or hotspots. Oza et al. (2019) introduced a proxy between the amount of outgassed material, the size, and the tidal heating rate of the body. This means that we can constrain the tidal heating rate by quantifying volcanically related gases (e.g., sodium, potassium) with spectroscopy, when given the size of the body. This measurement could serve as an independent way to deduce the tidal heating rate and, thus, the Im(k_{2}) of the moon via Eq. (2).
Hotspots induced by tidal heating in THEMs could also serve as a way to observe surface inhomogeneities. This has already been applied in exoplanets via single band photometry (Cowan & Fujii 2018; Majeau et al. 2012). Since different regions of the moon’s surface would be visible at different orbital phases, spatial inhomogeneities on the surface could be inferred. These variations depend on interior properties in tidally active bodies (Steinke et al. 2020). Thus, one could determine whether possible inferred inhomogeneities are compatible with the amount of tidal heating induced and the interior properties assumed. In addition, the presence of hotspots would be a possible way to constrain the moon’s orbit period (see Sect. 5.2).
Another possible consideration would be validating the calculated Im(k_{2}) through orbital constraints, in other words, quantifying whether the Im(k_{2}) is compatible with the orbital properties of the putative detection. The orbital circularization timescale depends on the Im(k_{2}) of the moon (RoviraNavarro et al. 2021). The eccentricity constraints that our analysis provides and the calculated Im(k_{2}) should be compatible with the system’s age. Such an analysis would require different scenarios to be evaluated, including the possibility for a resonance lock (Fuller et al. 2016) with the planet’s interior, as is likely in the case for Saturn’s satellite system (Lainey et al. 2020).
Finally, these considerations are not limited to moons but can also be explored in exoplanets, serving as a way to better understand the mechanism of tidal heating in rocky bodies (Henning et al. 2009; Renaud & Henning 2018). Applying the discussed model to exoplanets and validating the results would offer a bigger observational sample to test the model’s limits.
6 Conclusions
We used tidal and thermal models that describe Io’s observed surface heat flux (RoviraNavarro et al. 2021) to investigate the observability of tidally heated exomoons (THEMs) around ϵ Eridani b. We conclude that for orbital properties of the same order, such as solar system satellites (semimajor axis of 5.5 Roche radii and eccentricity of 0.009), a 2R_{Io} THEM is detectable around ϵ Eridani b with MIRI. When we took into account uncertainties in the interior properties, approximately 40% of such simulated THEMs reached a surface temperature high enough to be detectable. However, 100% of our simulations led to an observable SuperIo for a moon eccentricity of 0.02 and a semimajor axis of approximately 4 Roche radii. If no significant atmosphere is present, a larger moon would experience more intense tidal heating, and thus it would also be detectable.
With a successful thermal flux detection, the models we present could place constraints on orbital parameters, with the assumption that the satellites are in thermal equilibrium. If the satellite’s semimajor axis and radius are inferred (e.g., via exomoon transits), we can place stricter constraints on the orbital eccentricity and infer interior properties, such as melt fraction and thickness of the asthenosphere (see also Table 3) or whether the moon is in a magma ocean state. For instance, for detection of a THEM with a surface temperature of 270 K and a semimajor axis of 4.65 Roche radii, the exoIo would need to have an eccentricity smaller than 0.022 to be in thermal equilibrium. This constraint is looser for higher surface temperatures. In addition, only bodies with melt fractions larger than 0.1 can explain an observation of 470 K for this semimajor axis. A number of complementary observations of exomoons and their environments, such as hotspot properties, volcanic outgassing, and orbital configuration, could be used to test the applicability of our models after the detection of a THEM.
Acknowledgements
This research has been supported by the PEPSci Programme (Planetary and ExoPlanetary Science Programme), NWO, the Netherlands.
References
 Agol, E., Jansen, T., Lacy, B., Robinson, T. D., Meadows, V. 2015, ApJ, 812, 5 [NASA ADS] [CrossRef] [Google Scholar]
 Anderson, J. D., Sjogren, W. L., Schubert, G. 1996, Science, 272, 709 [NASA ADS] [CrossRef] [Google Scholar]
 Andrade, E. N. D. C. 1910, Proc. Roy. Soc. Lond. Ser. A, 84, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Baland, R.M., Yseboodt, M., van Hoolst, T. 2012, Icarus, 220, 435 [NASA ADS] [CrossRef] [Google Scholar]
 Barr, A. C., Dobos, V., Kiss, L. L. 2018, A&A, 613, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Benisty, M., Bae, J., Facchini, S., et al. 2021, ApJ, 916, L2 [NASA ADS] [CrossRef] [Google Scholar]
 Bierson, C. J., Nimmo, F. 2016, J. Geophys. Res.: Planets, 121, 2211 [NASA ADS] [CrossRef] [Google Scholar]
 Boccaletti, A., Lagage, P.O., Baudoz, P., et al. 2015, PASP, 127, 633 [CrossRef] [Google Scholar]
 Boccaletti, A., Cossou, C., Baudoz, P., et al. 2022, A&A, 667, A165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bohn, A. J., Kenworthy, M. A., Ginski, C., et al. 2019, MNRAS, 492, 431 [Google Scholar]
 Brandl, B., Bettonvil, F., Van Boekel, R., et al. 2021, The Messenger, 182, 22 [NASA ADS] [Google Scholar]
 Cabrera, J., Schneider, J. 2007, A&A, 464, 1133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Canup, R. M., Ward, W. R. 2006, Nature, 441, 834 [NASA ADS] [CrossRef] [Google Scholar]
 CastilloRogez, J. C., Efroimsky, M., Lainey, V. 2011, J. Geophys. Res. (Planets), 116, E09008 [CrossRef] [Google Scholar]
 Cowan, N. B., Fujii, Y. 2018, Mapping Exoplanets, eds. H. J. Deeg, J. A. Belmonte (Cham: Springer International Publishing), 1469 [CrossRef] [Google Scholar]
 Dermott, S. F., Malhotra, R., Murray, C. D. 1988, Icarus, 76, 295 [NASA ADS] [CrossRef] [Google Scholar]
 Dobos, V., Turner, E. L. 2015, ApJ, 804, 41 [NASA ADS] [CrossRef] [Google Scholar]
 Dobos, V., Charnoz, S., Pál, A., RoqueBernard, A., Szabó, G. M. 2021, PASP, 133, 094401 [NASA ADS] [CrossRef] [Google Scholar]
 Efroimsky, M. 2012, ApJ, 746, 150 [Google Scholar]
 Fischer, H.J., Spohn, T. 1990, Icarus, 83, 39 [NASA ADS] [CrossRef] [Google Scholar]
 Fuller, J., Luan, J., Quataert, E. 2016, MNRAS, 458, 3867 [NASA ADS] [CrossRef] [Google Scholar]
 Glasse, A., Rieke, G. H., Bauwens, E., et al. 2015, PASP, 127, 686 [NASA ADS] [CrossRef] [Google Scholar]
 Goldreich, P., Soter, S. 1966, Icarus, 5, 375 [NASA ADS] [CrossRef] [Google Scholar]
 Gribb, T. T., Cooper, R. F. 1998, J. Geophys. Res.: Solid Earth, 103, 27267 [NASA ADS] [CrossRef] [Google Scholar]
 Hatzes, A. P., Cochran, W. D., McArthur, B., et al. 2000, ApJ, 544, L145 [Google Scholar]
 Heller, R. 2014, ApJ, 787, 14 [CrossRef] [Google Scholar]
 Heller, R. 2017, Handbook of Exoplanets, eds. H. Deeg, J. Belmonte (Switzerland: Springer), 1 [Google Scholar]
 Heller, R., Hippke, M., Placek, B., Angerhausen, D., Agol, E. 2016, A&A, 591, A67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Henning, W. G., O’Connell, R. J., Sasselov, D. D. 2009, ApJ, 707, 1000 [NASA ADS] [CrossRef] [Google Scholar]
 Hussmann, H., Spohn, T. 2004, Icarus, 171, 391 [NASA ADS] [CrossRef] [Google Scholar]
 Jackson, I., Faul, U. H. 2010, Phys. Earth Planet. Interiors, 183, 151 [CrossRef] [Google Scholar]
 Jäger, Zoltán, J., Szabó, G. M. 2021, MNRAS, 508, 5524 [CrossRef] [Google Scholar]
 Janson, Markus, Quanz, Sascha P., Carson, Joseph C., et al. 2015, A&A, 574, A120 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 JaraOrué, H. M., Vermeersen, B. L. 2011, Icarus, 215, 417 [CrossRef] [Google Scholar]
 Karato, S.I., Wu, P. 1993, Science, 260, 771 [NASA ADS] [CrossRef] [Google Scholar]
 Katz, R. F. 2008, J. Petrol., 49, 2099 [NASA ADS] [CrossRef] [Google Scholar]
 Kenworthy, M. A., Mamajek, E. E. 2015, ApJ, 800, 126 [NASA ADS] [CrossRef] [Google Scholar]
 Kervazo, M., Tobie, G., Choblet, G., Dumoulin, C., Běhounková, M. 2022, Icarus, 373, 114737 [NASA ADS] [CrossRef] [Google Scholar]
 Kipping, D. M. 2009, MNRAS, 396, 1797 [NASA ADS] [CrossRef] [Google Scholar]
 Kipping, D. M., Bakos, G. Á., Buchhave, L., Nesvorný, D., Schmitt, A. 2012, ApJ, 750, 115 [NASA ADS] [CrossRef] [Google Scholar]
 Kipping, D., Bryson, S., Burke, C., et al. 2022, Nat. Astron., 6, 367 [NASA ADS] [CrossRef] [Google Scholar]
 Kovtyukh, V. V., Soubiran, C., Belik, S. I., Gorlova, N. I. 2003, A&A, 411, 559 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lainey, V., Arlot, J.E., Karatekin, O., Van Hoolst, T. 2009, Nature, 459, 957 [NASA ADS] [CrossRef] [Google Scholar]
 Lainey, V., Casajus, L. G., Fuller, J., et al. 2020, Nat. Astron., 4, 1053 [Google Scholar]
 Lammer, H., Schiefer, S.C., Juvan, I., et al. 2014, Origins Life Evol. Biosphere, 44, 239 [NASA ADS] [CrossRef] [Google Scholar]
 Lazzoni, C., Desidera, S., Gratton, R., et al. 2022, MNRAS, 516, 391 [NASA ADS] [CrossRef] [Google Scholar]
 Limbach, M. A., Turner, E. L. 2013, ApJ, 769, 98 [NASA ADS] [CrossRef] [Google Scholar]
 LlopSayson, J., Wang, J. J., Ruffio, J.B., et al. 2021, AJ, 162, 181 [NASA ADS] [CrossRef] [Google Scholar]
 MacGregor, M. A., Wilner, D. J., Andrews, S. M., Lestrade, J.F., Maddison, S. 2015, ApJ, 809, 47 [NASA ADS] [CrossRef] [Google Scholar]
 Majeau, C., Agol, E., Cowan, N. B. 2012, ApJ, 747, L20 [NASA ADS] [CrossRef] [Google Scholar]
 Makarov, V. V., Efroimsky, M. 2014, ApJ, 795, 7 [NASA ADS] [CrossRef] [Google Scholar]
 Mamajek, E. E., Hillenbrand, L. A. 2008, ApJ, 687, 1264 [CrossRef] [Google Scholar]
 Mawet, D., Hirsch, L., Lee, E. J., et al. 2019, AJ, 157, 33 [NASA ADS] [CrossRef] [Google Scholar]
 Mei, S., Bai, W., Hiraga, T., Kohlstedt, D. 2002, Earth Planet. Sci. Lett., 201, 491 [CrossRef] [Google Scholar]
 Meyer, J., Wisdom, J. 2007, Icarus, 188, 535 [NASA ADS] [CrossRef] [Google Scholar]
 Moore, W. B. 2001, Icarus, 154, 548 [Google Scholar]
 Moore, W. B. 2003, J. Geophys. Res.: Planets, 108 [Google Scholar]
 Ogihara, M., Ida, S. 2012, ApJ, 753, 60 [NASA ADS] [CrossRef] [Google Scholar]
 Oza, A. V., Johnson, R. E., Lellouch, E., et al. 2019, ApJ, 885, 168 [Google Scholar]
 Peale, S. J., Lee, M. H. 2002, Science, 298, 593 [NASA ADS] [CrossRef] [Google Scholar]
 Peale, S. J., Cassen, P., Reynolds, R. T. 1979, Science, 203, 892 [NASA ADS] [CrossRef] [Google Scholar]
 Peltier, W. R. 1974, Rev. Geophys., 12, 649 [NASA ADS] [CrossRef] [Google Scholar]
 Reese, C., Solomatov, V., Moresi, L.N. 1999, Icarus, 139, 67 [NASA ADS] [CrossRef] [Google Scholar]
 Renaud, J. P., Henning, W. G. 2018, ApJ, 857, 98 [NASA ADS] [CrossRef] [Google Scholar]
 Renaud, J. P., Henning, W. G., Saxena, P., et al. 2021, Planet. Sci. J., 2, 4 [NASA ADS] [CrossRef] [Google Scholar]
 Reynolds, R. T., McKay, C. P., Kasting, J. F. 1987, Adv. Space Res., 7, 125 [NASA ADS] [CrossRef] [Google Scholar]
 Rosenthal, L. J., Fulton, B. J., Hirsch, L. A., et al. 2021, ApJS, 255, 8 [NASA ADS] [CrossRef] [Google Scholar]
 RoviraNavarro, M., van der Wal, W., Steinke, T., Dirkx, D. 2021, Planet. Sci. J., 2, 119 [NASA ADS] [CrossRef] [Google Scholar]
 RoviraNavarro, M., Katz, R. F., Liao, Y., van der Wal, W., Nimmo, F. 2022, J. Geophys. Res.: Planets, 127, e07117 [NASA ADS] [CrossRef] [Google Scholar]
 Ruffio, J.B., Horstman, K., Mawet, D., et al. 2023, AJ, 165, 113 [NASA ADS] [CrossRef] [Google Scholar]
 Sabadini R., Vermeersen B., Cambiotti, G. 2016, Global Dynamics of the Earth: Applications of Viscoelastic Relaxation Theory to SolidEarth and Planetary Geophysics (2nd ed.; Berlin: Springer) [Google Scholar]
 Sartoretti, P., Schneider, J. 1999, A&AS, 134, 553 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Saumon, D., Hubbard, W. B., Burrows, A., et al. 1996, ApJ, 460, 993 [Google Scholar]
 Schubert, G., Cassen, P., Young, R. 1979, Icarus, 38, 192 [NASA ADS] [CrossRef] [Google Scholar]
 Segatz, M., Spohn, T., Ross, M. N., Schubert, G. 1988, Icarus, 75, 187 [NASA ADS] [CrossRef] [Google Scholar]
 Shoji, D., Kurita, K. 2014, ApJ, 789, 3 [NASA ADS] [CrossRef] [Google Scholar]
 Shoji, D., Hussmann, H., Kurita, K., Sohl, F. 2013, Icarus, 226, 10 [NASA ADS] [CrossRef] [Google Scholar]
 Simon, A., Szatmáry, K., Szabó, Gy. M. 2007, A&A, 470, 727 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Solomatov, V. 2007, in Treatise on Geophysics (Amsterdam: Elsevier), 91 [Google Scholar]
 Spencer, J. R., Rathbun, J. A., Travis, L. D., et al. 2000, Science, 288, 1198 [NASA ADS] [CrossRef] [Google Scholar]
 Spiegel, D. S., Burrows, A. 2012, ApJ, 745, 174 [NASA ADS] [CrossRef] [Google Scholar]
 Steinke, T., Hu, H., Höning, D., van der Wal, W., Vermeersen, B. 2020, Icarus, 335, 113299 [NASA ADS] [CrossRef] [Google Scholar]
 Teachey, A., Kipping, D. M. 2018, Sci. Adv., 4, eaav1784 [NASA ADS] [CrossRef] [Google Scholar]
 Tokadjian, A., Piro, A. L. 2022, ApJ, 929, L2 [NASA ADS] [CrossRef] [Google Scholar]
 Trani, A. A., Hamers, A. S., Geller, A., Spera, M. 2020, MNRAS, 499, 4195 [NASA ADS] [CrossRef] [Google Scholar]
 Tyler, R. 2011, Icarus, 211, 770 [NASA ADS] [CrossRef] [Google Scholar]
 Vanderburg, A., Rappaport, S. A., Mayo, A. W. 2018, AJ, 156, 184 [NASA ADS] [CrossRef] [Google Scholar]
 Veeder, G. J., Matson, D. L., Johnson, T. V., Blaney, D. L., Goguen, J. D. 1994, J. Geophys. Res.: Planets, 99, 17095 [NASA ADS] [CrossRef] [Google Scholar]
 Walterová, M., Běhounková, M. 2020, ApJ, 900, 24 [CrossRef] [Google Scholar]
 Walterová, M., Běhounková, M. 2017, Celest. Mech. Dyn. Astron., 129, 235 [CrossRef] [Google Scholar]
 Yoder, C. F., Peale, S. J. 1981, Icarus, 47, 1 [CrossRef] [Google Scholar]
All Tables
All Figures
Fig. 1 Feedback between the heat transfer and the tidal models. The interior structure affects the Im(k_{2}) value, which the tidal dissipation depends on. High tidal dissipation can cause the formation of melt, leading to a new sublayer (new interior structure; see the black arrows in the upper part of the plot). Thermodynamic equilibrium is reached when the moon’s interior transfers the same amount of heat as what is generated via tidal interactions (see the grey arrows in the lower part of the plot). 

In the text 
Fig. 2 Interior layers and the corresponding heat transfer mechanisms assumed in our model. 

In the text 
Fig. 3 Input and output parameters used in the thermal model. Descriptions of these values can be found in Tables 1, 3 and 2. 

In the text 
Fig. 4 Generated tidal heat flux as a function of the mantle temperature (T_{m}) and the log of the asthenosphere’s (η_{α}) and mantle viscosity (η_{m}). Stable and nonstable equilibria points are shown as black and gray dots, correspondingly. At a nonstable equilibrium point, a deviation in the mantle temperature drives the moon out of equilibrium, whereas a body at a stable thermal equilibrium state tends to restore its equilibrium. The tidal heat flux was calculated for a moon with a semimajor axis of a = 6.59 Roche radii, e = 0.02, α = 0.5, and B = 30 with the rest of the sensitivity parameters set to their nominal values (Table 2). 

In the text 
Fig. 5 Fluxes of hypothetical 2R_{Io} moons around e Eridani b. The black lines are the 10 σ and 10 000 s detection limits for MIRI (Glasse et al. 2015). The planet model (Spiegel & Burrows 2012) is shown in gray. The red and purple lines correspond to the black bodies of a 270 and 470 K 2 × R_{Io} exomoon, and the yellow dashed line corresponds to the black body of the star. 

In the text 
Fig. 6 Surface equilibrium temperature of a 2R_{Io} THEM as a function of its orbital eccentricity (e) and semimajor axis (a). The black dot indicates Io’s current orbital parameters. The nominal values of the interior parameters that we used are presented in Table 2. The 270 K isotherm divides the parameter space into detectable and nondetectable exomoons. 

In the text 
Fig. 7 Two different equilibrium interior structures for a = 4.65 Roche radii. The one on the left (e = 0.004) has a thinner asthenosphere compared to the one on the right (e = 0.02), resulting in a smaller Im(k_{2}). The parameter symbols are explained in Table 3. 

In the text 
Fig. 8 Asthenosphere melt fraction (ϕ_{α}) of stable thermal equilibrium states of exoIos as a function of their orbital eccentricity and semimajor axis. 

In the text 
Fig. 9 Range of exoIo orbital properties that are compatible with a surface temperature of 270 K and a radius of 2R_{Io} as a result of the sensitivity analysis. 

In the text 
Fig. 10 Fraction of simulations that lead to an observable moon (warmer than 270 K; Sect. 3) as a function of the moon's orbital semimajor axis and eccentricity. 

In the text 
Fig. 11 Computed ranges for the model’s output parameters (Table 3) that lead to exoIos with a = 4.65 Roche radii and a radius of 2 × R_{Io} in a stable thermal equilibrium. The blue line corresponds to the nominal values analysis (Sect. 5). The gray and red dotted lines respectively show the limits of a 270 K and 470 K surface temperature observation. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.