Three-Dimensional Interior Inversion of Icy Moons from Gravity, Topography, and Tidal Modeling

Introduction: The characterization of icy moon interiors is essential for assessing the physical conditions of their subsurface oceans and their potential habitability. Icy moons represent key targets because of their potential active geology, evidence for internal oceans, and relevance to ongoing and future planetary missions. Constraining their internal structure requires estimates of the properties and geometry of the rocky core, ocean, and ice shell.

Mass and low-degree gravity coefficients derived from spacecraft tracking data provide constraints on the bulk distribution of mass and on the total thickness of the hydrosphere. However, these quantities are not sufficient to independently determine the thicknesses and properties of the ice shell and subsurface ocean. The tidal Love number k₂, which depends on the global viscoelastic response of the body, provides an additional constraint on the rheological properties of the internal layers. Lateral variations in ice-shell thickness can be investigated using the high-degree components of the topography, which are sensitive to shallow and short-wavelength structures, together with an assumed compensation mechanism at the ice-ocean interface.

This study applies a Bayesian inversion algorithm based on a Markov Chain Monte Carlo (MCMC) method [1] to combine geophysical constraints and estimate the radii, densities, rheological properties, shell thickness, and synthetic gravity field of Enceladus and Europa. The goal is to develop a generalized framework for the three-dimensional characterization of icy moon interiors.

Modeling: Interior structure inversions of icy moons are commonly based on one-dimensional radial models [2], which are well suited to investigate the global-scale properties of the internal layers, including their radii, densities, and rheological parameters. However, this approach does not exploit geophysical measurements that contain information on lateral variations, such as gravity and topography at short and intermediate wavelengths. To include these additional constraints, this study adopts an enhanced interior model that accounts for the shape of the surface and uses the gravity–topography relationship to indirectly constrain lateral variations in the ice shell and the ice–ocean interface.

This extension requires assumptions on the degree and mechanism of internal compensation. At long wavelengths, the gravity and topography of synchronously rotating icy moons are expected to be dominated by the hydrostatic response to rotation and tidal forcing, with only small deviations from hydrostatic equilibrium. This component is represented by a set of concentric triaxial ellipsoids, whose semiaxes are computed from the rotation period, reference layer radii, and densities by imposing equipotential conditions at the endpoints of the principal axes, following the minimization algorithm proposed by [3]. The remaining short-wavelength component is treated as non-hydrostatic and is associated with shallow sources within the ice shell. As a preliminary analysis, this non-hydrostatic contribution is modeled as the result of the viscoelastic response of the icy shell and of compensation at the ice–ocean interface.

The ice shell is therefore represented as the sum of a hydrostatic reference figure and a non-hydrostatic contribution associated with surface and basal topography. The latter represents lateral variations in shell thickness and is linked to the gravity field through a degree-dependent admittance model. This formulation allows the inversion to retain the efficiency of a layered radial description for the global interior while incorporating lateral information from gravity and topography to constrain the three-dimensional structure of the outer shell.

 

Results: We present a Bayesian inversion framework for the interior structure of icy moons based on a MCMC approach. The inversion uses the planetary mass and low-degree gravity field as primary geophysical constraints, together with interior parameters consistent with available observational bounds and accuracies based on legacy exploration missions. Since the tidal Love numbers of icy moons are not yet sufficiently constrained by observations, k₂ and h₂ are not imposed as target quantities. Instead, they are treated as model predictions derived from the ensemble of accepted interior structures.

The interior models are generated by sampling the mean radii, densities, and rheological properties of the core, ocean, and ice shell. For each model, the hydrostatic shape of the internal interfaces is computed using a fourth-order formulation, while the ice-shell bottom topography is inferred from the observed surface topography and model parameters by imposing a compensation mechanism at the ice-ocean interface. The resulting three-dimensional shell-thickness distributions are used to compute the associated gravity field, which is included as a target observable in the inversion.

For each accepted model, the degree-2 tidal Love numbers k₂ and h₂ are computed with PyALMA [4], which provides the viscoelastic tidal response of multilayer planetary interiors. Additional calculations with LOV3D [5] are used to evaluate the effect of lateral heterogeneity and compressibility on the tidal response, including whether shell-thickness variations inferred from gravity and topography produce measurable perturbations in k₂ and h₂. This comparison quantifies whether the tidal response can be interpreted using one-dimensional interior models or whether three-dimensional structure introduces non-negligible corrections.

The inversion provides posterior distributions for the mean radii, densities, and rheological properties of the internal layers, together with maps of ice-shell thickness and their associated uncertainties. Rather than fitting observed tidal parameters, the study predicts the expected ranges of k₂ and h₂ and evaluates their sensitivity to interior structure, shell-thickness variations, and compressibility. This approach supports the definition of geophysical observables that future missions could use to distinguish between competing interior models of icy moons.

 

Conclusion: This study presents a Bayesian framework for the three-dimensional characterization of icy moon interiors by combining geophysical constraints, gravity-topography relationships, and tidal modeling. The approach retrieves interior models consistent with the mass and gravity field while predicting the corresponding tidal Love numbers k₂ and h₂ from the posterior distributions and layers’ rheology and compressibility. This framework can support the identification of the measurements required to distinguish among competing interior models and to improve the characterization of the habitability-relevant properties of icy moons.

 

[1] Mosegaard et al. (1995) JGR Solid Earth, 100.B7: 12431-12447

[2] Genova et al. (2019) GRL, 46.7: 3625-3633

[3] Tricarico et al. (2014) The Astrophysical Journal, 782.2: 99-111

[4] Petricca et al. (2024) Icarus 417

[5] Navarro et al. (2024) The Planetary Science Journal 5: 129-152