The bedrock deformation in response to a melting ice sheet provides negative feedback on ice mass loss. When modelling the future behaviour of the Antarctic Ice Sheet, the impact of bed deformation on ice dynamics varies but can reduce projections of future sea-level rise by up to 40 % in comparison with scenarios that assume a rigid Earth. The rate of the solid Earth response is mainly dependent on the viscosity of the Earth's mantle, which varies laterally and radially with several orders of magnitude across Antarctica. Because modelling the response for a varying viscosity is computationally expensive and has only recently been shown to be necessary over centennial time scales, sea-level projection ensembles often exclude the Earth's response or apply a globally constant relaxation time or viscosity. We use a coupled model to investigate the accuracy of various approaches to modelling the bedrock deformation to ice load change. Specifically, we compare the sea-level projections from an ice-sheet model coupled to (i) an elastic lithosphere, relaxed asthenosphere (ELRA) model, with either uniform and laterally varying relaxation times, (ii) a glacial isostatic adjustment (GIA) model with a radially varying Earth structure (1D GIA model), and (iii) a GIA model with laterally varying earth structures (3D GIA model). Furthermore, using the 3D GIA model we determine a relation between relaxation time and viscosity which can be used in ELRA and 1D models. We conduct 500-year projections of Antarctic Ice Sheet evolution using two different climate models and two emissions scenarios: the high emission scenario SSP5-8.5 and the low emission scenario SSP1-2.6. Using a rigid Earth model, this results in ∼3–7.5 m of barystatic sea-level rise with significant retreat in various basins due to marine ice sheet instability. The results show that using a uniform relaxation time of 300 years in an ELRA model leads to a total sea-level rise that deviates less than 40 cm (6 %) from the average of the 3D GIA models in 2500. This difference in the projected sea-level rise can be further reduced to 20 cm (4 %) by using an upper mantle viscosity of 1019 Pa s in the 1D GIA model, and to 10 cm (2 %) in 2500 by using a laterally varying relaxation time map in an ELRA model. Our results show that the Antarctic Ice Sheet contribution to sea-level rise can be approximated sufficiently accurate using ELRA or a 1D GIA model when the recommended parameters derived from the full 3D GIA model are used.

Io’s widespread volcanic activity is driven by strong tidal heating within its interior [1]. The distribution of tidal heating within Io's interior is closely linked to the moon's thermal state and its evolution [2,3]. Yet, despite years of observations [4], where inside Io heat is being dissipated remains largely unknown. For example, it remains unclear whether an asthenosphere (a low-viscosity layer beneath the crust) exists [4]. Ground-based and space-borne observations provide a window into Io's interior. With the recent measurement of Io’s tidal Love number k2 [5], there are now four independent gravity observations: C20, C22, and both parts of k2 [5,6]. This makes Io, together with Titan [7], the outer solar system moon for which we have the most gravity observations. Furthermore, there is a large database of volcanic activity observations [e.g., 4], assumed to link more directly to Io's tidal heating pattern. In turn, the heating pattern is sensitive to Io’s interior radial rheology structure [2], complementing the tidal response observations. However, no comprehensive Bayesian inversion exists that incorporates both gravity and volcanic activity observations. We investigate how combining all gravity observations and the observed volcanic heat flux constrains Io’s interior. As a measure of how Io’s heat output is distributed, we use the ratio of polar to equatorial heat flux. Our interior model consists of four layers: core, mantle, asthenosphere, and lithosphere, and we compute the tidal response using Andrade rheology [8]. We use a Markov Chain Monte Carlo approach to integrate all gravity observations and the observed heat flux ratio to constrain Io’s interior. Different observables are sensitive to different interior parameters. The static gravity coefficients, C20 and C22, constrain the core density and radius, k2-measurements the effective rigidity and viscosity, and the heat flux ratio the radial viscosity profile. The latter requires the majority of the heating (>90%) to occur in a low viscosity asthenosphere (1013 Pa s), even with conservative assumptions on the distribution of Io’s unmeasured background heat flux. We thus show that, within current observational limits, Io must have an asthenosphere where a significant amount of its heat is generated. The existence of an asthenosphere is in line with previous results [9] and also expected from melt advection models [e.g., 10]. While we used Andrade rheology, we find viscosities that imply melt fractions at the edge of what is considered solid [e.g., 11], raising questions about the validity of rheology models typically used to model tidal deformation in planetary objects. However, experimental work on anelastic deformation under Io-like conditions is likely needed to resolve this problem [8]. Finally, our study provides statistically consistent ranges for Io’s interior parameters, which could be used in future studies. Figure 1. The posterior PDFs of the ratio between energy dissipation in the asthenosphere (Ea) and mantle (Em) (subplot a), the polar-to-equatorial heat flux ratio (subplot c), and the ratio between asthenosphere and mantle viscosity ηx (subplot f). The 2-dimensional PDFs, indicating the correlations between two parameters, are given in the off-diagonal subplots. The three cases, without the observed polar flux ratio as one of the observables (blue), with the observed ratio (red), and with a uniform flux ratio (green), are plotted on top of each other. The dotted lines in the histograms and the contours in the 2D PDFs contain 68.3% of the cases. References [1] Peale, S. J., Cassen, P., & Reynolds, R. T. 1979, Science, 203, 892 [2] Segatz, M., Spohn, T., Ross, M. N., & Schubert, G. 1988, Icarus, 75, 187 [3] Hussmann, H. & Spohn, T. 2004, Icarus, 171, 391 [4] Davies, A. G., Perry, J. E., Williams, D. A., Veeder, G. J., & Nelson, D. M. 2024, The Planetary Science Journal, 5, 121 [5] Park, R. S., Jacobson, R. A., Gomez Casajus, L., et al. 2025, Nature, 638, 69 [6] Lainey, V., Arlot, J.-E., Karatekin, ¨O., & van Hoolst, T. 2009, Nature, 459, 957 [7] Petricca, F., Vance, S. D., Parisi, M., et al. 2025, Nature, 648, 556 [8] Bierson, C. J. 2024, Icarus, 414, 116026 [9] Bierson, C. J. & Nimmo, F. 2016, Journal of Geophysical Research (Planets), 121, 2211 [10] Spencer, D. C., Katz, R. F., Hewitt, I. J., May, D. A., & Keszthelyi, L. P. 2020, Journal of Geophysical Research (Planets), 125, e06604 [11] Scott, T. & Kohlstedt, D. L. 2006, Earth and Planetary Science Letters, 246, 177

Northern Europe experiences vertical land motion and sea level changes as a consequence of past changes in ice sheet cover in Fennoscandia and the British Isles. The process, called glacial isostatic adjustment (GIA), is controlled by the subsurface structure. Numerical models of GIA can be compared to observations of uplift or past sea level changes to constrain the subsurface structure, and such models can also be used to correct present-day sea level observations to reveal sea level changes due to climate change. GIA models for northern Europe usually adopt a homogeneous upper mantle viscosity even though seismic studies indicate contrasting elastic lithosphere thickness and upper mantle structure between Northwestern Europe and Eastern Europe. This raises the question whether the effect of lateral variations in structure (3-D viscosity) can be detected in observations of GIA and whether including such variations can improve GIA model predictions. In this study, we compare model output from a finite element GIA model with 3-D viscosity to observations of paleo sea level and current vertical land motion. We use two different methods to derive 3-D viscosities, based on seismic velocity anomalies and upper mantle temperature estimates. We use three different reconstructions of the Eurasian ice sheet, one based on an inversion using a 1-D viscosity model, and two others based on glacial geology and modelling. When we use these two reconstructions, we find that the data are fit better using 3-D viscosity models. Models with two separate 1-D viscosities for Fennoscandia and for the British Isles cannot replicate a 3-D model because a 3-D model redistributes GIA-induced stresses differently from a combination of models with 1-D viscosities. The fit to data across Fennoscandia is improved when, as indicated by seismic models, the upper mantle viscosity is higher than for the rest of Northern Europe. The best fit is obtained with a model with dry olivine rheology, in agreement with other evidence from Fennoscandia.

The innermost Galilean moon, Io, exhibits widespread tidally-driven volcanism. Monitoring of its volcanoes has revealed that they are not homogeneously distributed across its surface: volcanic activity is higher at low latitudes and peaks east of the sub- and anti-Jovian points. Dissipation in a radially symmetric solid body cannot explain the observed longitudinal shift but dissipation in a magma ocean can. However, recent observations show that Io does not have one. Here, we demonstrate that a longitudinal shift in the heating pattern naturally arises from the feedback between tidal heating and melt production. The feedback between tidal dissipation and interior properties that results in interiors that deviate from radial symmetry is expected to drive the interior evolution of other tidally-active worlds, including icy moons such as Europa and Enceladus and exo-planets/moons with high eccentricity or obliquity.

Rates of relative sea-level rise during the final stage of the last deglaciation, the early Holocene, are key to understanding future ice melt and sea-level change under a warming climate1. Data about these rates are scarce2, and this limits insight into the relative contributions of the North American and Antarctic ice sheets to global sea-level rise during the early Holocene. Here we present an early Holocene sea-level curve based on 88 sea-level data points (13.7-6.2 thousand years ago (ka)) from the North Sea (Doggerland3,4). After removing the pattern of regional glacial isostatic adjustment caused by the melting of the Eurasian Ice Sheet, the residual sea-level signal highlights two phases of accelerated sea-level rise. Meltwater sourced from the North American and Antarctic ice sheets drove these two phases, peaking around 10.3 ka and 8.3 ka with rates between 8 mm yr-1 and 9 mm yr-1. Our results also show that global mean sea-level rise between 11 ka and 3 ka amounted to 37.7 m (2σ range, 29.3-42.2 m), reconciling the mismatch that existed between estimates of global mean sea-level rise based on ice-sheet reconstructions and previously limited early Holocene sea-level data. With its broad spatiotemporal coverage, the North Sea dataset provides critical constraints on the patterns and rates of the late-stage deglaciation of the North American and Antarctic ice sheets, improving our understanding of the Earth-system response to climate change.

This study is concerned with the total rate of mass change in the oceans as seen by the GRACE and GRACE-FO missions between Mar-2002 extending to beyond 2023 into the GRACE-FO observation window. We use the same monthly source data in the form of level-2 spherical harmonic coefficients and find that the selected glacial isostatic adjustment models, ocean mask definitions, and other de-aliasing corrections have a significant effect on the outcome of the barystatic ocean mass rates. In addition to the ICE6G models, and GIA models by Caron and Ivins we also implement GIA models derived from a 3D model forced by ICE-7G ice history. In the end we find barystatic ocean mass rates of change varying between 1.5 to 2.4 mm/yr. The relevance of this result is that the mass component of the sea level change has a larger uncertainty than that we earlier expected. Our conclusion is that it is more difficult to close the sea level budget which involves ocean volume estimates from satellite altimetry, density estimates from Argo profiling float data and the mass component of sea level changes observed by satellite gravimetry.

The contribution of the Antarctic Ice Sheet to barystatic sea-level rise could be as high as eight metres around 2300 but remains deeply uncertain. Ice sheet retreat causes bedrock uplift, which can exert a stabilising effect on the grounding line. Yet, sea-level projections exclude bedrock adjustment, use simplified Earth structures or omit the uncertainty in climate response and Earth structure. We show that the grounding line retreat is delayed by 50 to 130 years and the barystatic sea-level contribution reduced by 9-23% when the heterogeneity of the solid Earth is included in a coupled ice - bedrock model under different emission scenarios till 2500. The effect of the solid Earth feedback in ice sheet projections can be twice as large as the uncertainty due to differences between climate models. We emphasise that realistic Earth structures should be considered when projecting the Antarctic contribution to barystatic sea-level rise on centennial time scales.

Intensive tidal heating makes the Galilean satellite Io to an outstanding example of a volcanically active world. Most of the heat generated in the interior is lost through a large number of active volcanoes. The distribution of Io's volcanoes on the surface could help us to constrain properties below Io’s crust, regulating the heat transport mechanism. For this study, we assume that (1) the presence of global volcanism is linked to the presence of melt in the upper mantle; that (2) the large-scale variation in volcanic density is inherited from non-uniform tidal heating and smoothed by vigorous convection; and that (3) the total number of hot-spots is controlled by the spatial frequency of thermal instabilities in the convecting layer. Three unknown parameters are explored: the fraction of convective heat transport compared to magmatic heat transport, the mantle viscosity, and the thickness of the heated layer. In order to evaluate which combinations of interior properties can explain Io's present volcanic distribution, we develop a model based on parameterised heat flow scalings to approximate different spatial characteristics of Io's interior convection pattern. Our model combines internal heating, and convective and magmatic heat flow. Parts of the parameter space that are not in agreement with the observation-derived conditions are ruled out.Our results show that the observed small- and large-scale characteristics of Io's volcanic pattern can be explained by sub-lithospheric anomalies caused by convection. Solutions that allow for active volcanism and Io's specific large-scale variations in volcanic activity range from a thick mantle of high viscosity (1017 Pa s) to a thin asthenosphere of low viscosity (1012 Pa s). If Io's volcanos are correlated to the spatial frequency of thermal instabilities, the range of Io's total volcanic features between 250 and 3030 can further constrain the parameter space. This favours a mantle with a low melt fraction, a low mantle viscosity, and a magmatic heat transport of >80%.

Models of glacial isostatic adjustment (GIA) play a central role in the interpretation of various geologic and geodetic data to understand and simulate past and future changes in ice sheets and sea level, and infer rheological properties of the deep Earth. A relatively recent advance has been the development of models that include 3D Earth structure, as opposed to 1D, spherically symmetric structure. However, a major limitation in employing 3D GIA models is their high computational expense. As such, we have developed a method using artificial neural networks (ANNs) and the Tensorflow library to emulate the influence of 3D Earth models with the goal of more affordably constraining the parameter space of these models: specifically the radial (1D) viscosity profile upon which the lateral variations are added. This study provides an initial “proof of concept” assessment of using ANNs to emulate the influence of lateral Earth structure on GIA model output. Our goal is to test whether the fast surrogate model can accurately predict the difference in these outputs (i.e., RSL and uplift rates) for the 3D case relative to the SS case. If so, the surrogate model can be used with a computationally efficient SS (Earth) GIA model to generate output that reproduces output from a 3D (Earth) GIA model. Evaluation of the surrogate model performance for deglacial RSL indicates that it is able to provide useful estimates of this field throughout the parameter space when trained on only ≈ 15 % (≈ 50) of the parameter vectors considered (330 in total). Our results indicate that the ANN:model misfits, while not negligible, are of a scale such that useful predictions of deglacial RSL changes can be made. We applied the surrogate model in a model:data comparison exercise using RSL data distributed along the North American coasts from the Canadian Arctic to the US Gulf coast. We find that the surrogate model is able to successfully reproduce the data:model misfit values such that the region of minimum misfit either overlaps the 3D GIA model results, or is within two increments in the parameter space. The surrogate model can, therefore, be used to accurately explore this aspect of the 3D Earth model parameter space. While the 3D Earth models can outperform the SS Earth models for some regional subsets of the RSL data set, the SS Earth models still produce better fits overall. In summary, this work demonstrates the utility of machine learning in 3D Earth GIA modelling and so future work to expand on this analysis is warranted.

Abstract. The bedrock response to a melting ice sheet provides a negative feedback on ice mass loss. When modelling the future behaviour of the Antarctic Ice Sheet, accounting for the impact of bed deformation on ice dynamics can reduce predictions of future sea level rise by up to 40 % in comparison with scenarios that assume a rigid Earth. The rate of the solid Earth response is mainly dependent on the viscosity of the Earth’s mantle, which varies laterally and radially with several orders of magnitude across Antarctica. Because modelling the response for a varying viscosity is complex, sea level projections often exclude the Earth’s response, or apply a globally constant relaxation time or viscosity. We investigate how accurate such approximations are using an ice sheet model coupled with a glacial isostatic adjustment (GIA) model that simulates the bedrock response to changes in ice loading including lateral and radial variations in viscosity (3D GIA model). Using the 3D model we determined a relation between relaxation time and viscosity which can be used in simpler models. We compare output from an elastic lithosphere relaxed asthenosphere (ELRA) with uniform and laterally varying relaxation times, and a GIA model with a radially varying Earth structure (1D GIA model) and a 3D GIA model. We conducted 500 year projections of the Antarctic ice sheet evolution using two different climate models and two emissions scenarios: the high emission scenario SSP5-8.5 and the low emission scenario SSP1-2.6. The results show that using a uniform relaxation time of 300 years in the ELRA model or an upper mantle viscosity of 1019 Pa∙s in the 1D GIA model leads to a total sea level rise that deviates less than 40 cm from the average of the 3D GIA models. The difference in the sea level rise predicted with 1D and 3D GIA models can be further reduced to 10 cm by using laterally varying relaxation time maps in an ELRA model. Our results show that the effect of 3D viscosity variations on the AIS contribution to sea level rise can be approximated using the ELRA model or a 1D GIA model when the recommended parameters derived from the full 3D GIA model are used.

Introduction Io exhibits widespread volcanism powered by tides raised by Jupiter. The distribution of volcanoes offers a window into the interior of the moon. The distribution shows more volcanism at the equator as well as peak volcanic output which is shifted by roughly 30-60 degrees to the east of the subjovian point [1]. Models of tidal dissipation that assume a spherically symmetric, solid Io cannot reproduce this shift [2]. More recently, it has been proposed that tidal dissipation in a magma ocean [3] or in a non-spherically symmetric, solid Io [4] can induce this lag. In this study, we explore the second option and show that solid-body dissipation can induce an eastward shift of the tidal dissipation pattern. Method The amount of tidal dissipation experienced by a planet or moon is a function of its interior properties (e.g. the shear modulus and viscosity). These properties depend on the body's thermal state, which for Io is mainly dictated by tidal heating. This results in feedback between interior properties and tidal response, modelling this requires coupling a thermal and tidal model. For the former, we parameterize the dependence of viscosity, η(θ, φ), and shear modulus, μ(θ, φ), on melt-fraction using laws derived from laboratory experiments and we assume that tidal dissipation and melt anomalies are linearly related by a factor c as δΦ(θ, φ) = c(θ, φ) δQ(θ, φ) [2]. The proportionality factor c(θ, φ) ranges between 0 and 0.06 [2] and its dependency on latitude and longitude mimics convection and other forms of lateral heat transport. In this work, we assume c=0.015 and consider that heating patterns of spherical harmonic degrees higher than 4 are blurred by such processes [5]. We compute tidal dissipation using the spectral code LOV3D [6]. We assume there is only tidal heating in the asthenosphere [1] and set the rheological properties such that the dissipation results in an average surface heat flux consistent with observations, at 2.4 Wm-2 [7]. We further assume that the asthenosphere has a globally constant, average melt-fraction of 10% [8]. We aim to obtain a steady state. To do so we first convert the tidal dissipation pattern of a uniform Io into a map of viscosity and shear modulus anomalies. The resulting maps are used to recompute tidal heating and iterate. To ensure a steady state and prevent runaway heating, we normalize the tidal dissipation pattern such that the average does not vary between iterations. Figure 1: The surface heat flux pattern after one iteration (so with the introduction of lateral variations in the rheology) minus the initial surface heat flux pattern. Results Figure 1 shows the difference in surface heat flux, which is the radially integrated tidal dissipation, after one iteration compared to the initial pattern. The pattern shows a slight shift to the east of the anti-subjovian point, located at zero degrees longitude. Additionally, the local increase in melt-fraction induces higher peak-to-peak variations of tidal heating. After 89 iterations, a steady solution is found and we retrieve a difference in surface heat flux as shown in Fig. 2. The eastward shift has grown to roughly 35 degrees at the equator, and the peak-to-peak variations have also grown. Tests with different values for c and other coupling methods, not shown here, also show this eastward shift but the resulting pattern can be different. In general, a stronger coupling between dissipation and melt-fraction, represented by a higher c, causes a larger shift and bigger peak-to-peak values which is consistent with previous results that used finite element modelling [4]. Given that we started with a spherically symmetric interior it might be surprising that an asymmetry arises. We suspect the cause for the eastward shift is an inherent asymmetry in the forcing since the tidal potential can be represented as having an eastward and westward propagating component, of which the eastward component has a larger amplitude [9]. To confirm this, we set the amplitude of the westward and eastward components to the same value and no longer observe a shift. Overall, our results show that the feedback between tidal heating and rheology can induce an eastward lateral shift in tidal heating and hence melt-fraction. This demonstrates that solid-body tides can account for the observed eastward shift in the volcanic heat flow pattern which is especially relevant in anticipation of Juno’s measurement of Io’s k2 Love number, which could confirm or rule out the presence of a magma ocean. Our results are not particular to Io, hence suggesting that the same mechanism might lead to 3D structure in tidally heated exoplanets and exomoons. Figure 2: Surface heat flux pattern at convergence after 89 iterations minus the initial pattern. Bibliography [1] Davies, A. G., Perry, J. E., Williams, D. A., & Nelson, D. M. 2024, Nature Astronomy, 8, 94; [2] Steinke, T., Hu, H., H ̈oning, D., van der Wal, W., & Vermeersen, B. 2020a, Icarus, 335, 113299; [3] Tyler, R. H., Henning, W. G., & Hamilton, C. W. 2015, The Astrophysical Journal Supplement Series, 218, 22; [4] Steinke, T. 2021, PhD thesis, Delft University of Technology; [5] Beuthe, M. 2013, Icarus, 223, 308; [6] Rovira-Navarro, M., Matsuyama, I., & Berne, A. 2024, Planetary Science Journal, 5; [7] Lainey, V., Arlot, J.-E., Karatekin, ̈O., & van Hoolst, T. 2009, Nature, 459, 957; [8] Spencer, D. C., Katz, R. F., & Hewitt, I. J. 2020, Journal of Geophysical Research (Planets), 125, e06443; [9] Chen, E. M. A., Nimmo, F., & Glatzmaier, G. A. 2014, Icarus, 229, 11

Polar wander is the reorientation of a body following a perturbation in the moment of inertia due to a mass anomaly. It has occurred on Mars as a result of the formation of the Tharsis region (e.g. Matsuyama and Manga 2010), on Pluto due to volatile ice deposition (Keane et al. 2016) and has also been proposed for icy moons such as Europa (Schenk et al. 2020). Pluto and icy moons have a tidal bulge in addition to being flattened due to rotation. In that case the polar wander path is more complicated and involves both a reorientation around the tidal axis (pointing to the central body) and one around the rotational axis. The path and speed of the anomaly depend on the size and location of the load and the mantle viscosity, with the end location determined by the strength of the outer shell with very high viscosity. Thus, if there are observational clues about polar wander, different scenarios in terms of size and timing of load, and internal structure can be tested. Most studies assume approximate solutions or focus on the end points of polar wander. Sometimes small angle polar wander is assumed, such as for polar wander due to Pleistonene deglaciation for Earth. In that case the linearized form of the Liouville equation is used where it is assumed that the rotational axis is the z-axis and the equatorial bulge is perpendicular to that. Here we present sensitivity studies of the complete path for large angle polar wander to changes in loading.We use a semi-analytical solution for reorientation of a rotating tidally deformed body with a high-viscosity visco-elastic shell (Hu et al. 2019) valid for large polar wander (>10 degrees). The method does not assume that the body is fully relaxed at any moment of the reorientation (the so-called fluid limit solution). Therefore, it works for faster loads (e.g. due to impacts) and can provide the complete reorientation path. Results of this method have recently been reproduced by Patočka and Kihoulou (2023) using a simpler approach.We use Triton as case study. Polar wander on Triton can occur because of volatile migration (Rubincam 2003). Here we assume a point load. For Triton we use a 5-layer model as given in Hu et al. (2019) with a mantle viscosity of 1019 Pa s, overlain by ice shells with a viscosity of 1021 Pa s.Figure 1: The movement of a mass anomaly emplaced on the Triton model. The rotational axis is at the origin pointing out of plane, the tidal axis pointing to the central body is at 0° longitude, pointing to the right. The starting position of the mass anomaly is at 15 ° colatitude and 15 ° longitude. Circles: Heaviside load for a point mass of 1.5×1017 kg, open squares: Heaviside load for a point mass of 3×1017 kg which roughly corresponds to a 100m thick disc of ice of 100 km radius. Crosses: load increasing linearly in time over 15 Ma to a magnitude of 3×1017 kg.We apply a point mass with a magnitude of 3×1017 kg which roughly corresponds to a 100 m thick disc of ice of 100 km radius. The load is applied as instantaneous (Heaviside) forcing or linearly increasing (ramp). The path of a mass anomaly in the so-called bulge-fixed reference frame is shown in figure 1. For a larger mass the anomaly ends up closer to the tidal axis, as is expected. The ramp load leads to a path that is not very different from that of the Heaviside load. However, if the time over which the load increases is increased, the path is closer to that of the smaller Heaviside load (not shown). Thus, the path of the mass anomaly is an important constraint on the magnitude and timing of the load.In future work we aim to study Pluto. The Sputnik Planitia region is the remnant of an impact. This by itself would cause a negative mass anomaly which would induce polar wander such that it moves towards the north pole. However, if the impact basin that was compensated by uplift in the subsurface ocean (Keane et al. 2016) possible compounded with the buried remnant of the impactor (Ballantyne et al 2024) it would shift towards the sub-Charon point. We will investigate various initial positions and thickness of a nitrogen ice layer and study possible end locations.ReferencesBallantyne, H.A., Asphaug, E., Denton, C.A., Emsenhuber, A. and Jutzi, M., 2024. Sputnik Planitia as an impactor remnant indicative of an ancient rocky mascon in an oceanless Pluto. Nature Astronomy, pp.1-8.Hu, H.X.S., van der Wal, W. and Vermeersen, L.L.A., 2019. Rotational dynamics of tidally deformed planetary bodies and validity of fluid limit and quasi-fluid approximation. Icarus, 321, pp.583-592.Keane, J.T., Matsuyama, I., Kamata, S. and Steckloff, J.K., 2016. Reorientation and faulting of Pluto due to volatile loading within Sputnik Planitia. Nature, 540(7631), pp.90-93.Matsuyama, I. and Manga, M., 2010. Mars without the equilibrium rotational figure, Tharsis, and the remnant rotational figure. Journal of Geophysical Research: Planets, 115(E12).Patočka, V. and Kihoulou, M., 2023. Dynamic reorientation of tidally locked bodies: Application to Pluto. Earth and Planetary Science Letters, 617, p.118270.Rubincam, D.P., 2003. Polar wander on Triton and Pluto due to volatile migration. Icarus, 163(2), pp.469-478.Schenk, P., Matsuyama, I. and Nimmo, F., 2020. A very young age for true polar wander on Europa from related fracturing. Geophysical Research Letters, 47(17), p.e2020GL088

The Cassini mission has provided measurements of the gravity of several moons of Saturn as well as an estimate of the tidal response, which is expressed as the degree 2 Love number k2 of its largest moon, Titan. The first estimates of Titan's Love number were larger than pre-Cassini expectations. Interior modelling suggested it may be explained with a dense ocean, but the interpretation remains unclear. We analysed Cassini tracking data to determine Titan's gravity field and its Love number. Our gravity results are consistent with earlier studies, but we find a lower Love number for Titan of k2 = 0.375 +/- 0.06. This lower value follows from an elaborate investigation of the tidal effects. We show that a dense ocean is not implied by the obtained Love number; instead, a water or ammonia ocean is more probable. A lower density ocean can increase the likeliness of contact between the silicate core and ocean, which can leach minerals into the ocean and could promote its habitability.

ANET-POLENET (Antarctic Network of the Polar Earth Observing Network) bedrock GNSS sites in the Ross Sea region of Antarctica surround an LGM load center in the Siple region of the Ross Embayment and record crustal motion due to GIA. Rather than a radial pattern of horizontal motion away from the former load, we instead observe three primary patterns of deformation; 1) motions are reversed towards the load in the southern region of the Transantarctic Mountains (TAM), 2) motions are radially away from the load in the Marie Byrd Land (MBL) region, and 3) an overall gradient in motion is present, with magnitudes progressively increasing from East to West Antarctica. We investigate the effects of alternative Earth model and ice loading scenarios, with the goal of understanding these distinct patterns of horizontal bedrock motion and their drivers. Using GIA models with a range of 1D Earth models, alternative ice loading scenarios for the Wilkes Subglacial Basin (LGM time scale) and the Siple Coast (centennial and millennial time scales) are explored. We find that no 1D model, regardless of the Earth model and ice loading scenario used, reproduces all three distinct patterns of observed motion at the same time. For select ice loading scenarios we also examine the influence of more complex rheology by invoking a boundary in Earth properties beneath the Transantarctic Mountains. This approach accounts for the strong lateral gradient in Earth properties across the continent by effectively separating East and West Antarctica into two different Earth model profiles. Some of our GIA models utilizing 3D Earth structure reproduce predicted motions that match all three observed patterns of deformation, and we find that a multiple order magnitude of change in upper mantle viscosity between East and West Antarctica is required to fit the observations.

This Memoir is the first dedicated to the Antarctic mantle. It is a cross-disciplinary reference work combining geochemistry and geophysics to characterize Antarctic mantle properties. Through observations and modelling the mantle structures, compositions and dynamics are characterized at regional and continental scales by subject experts. The Memoir reviews all known occurrences of sub-continental mantle xenoliths in igneous rocks. These studies are presented by region as southern or northern Victoria Land, Marie Byrd Land, the Antarctic Peninsula, East Antarctica and the sub-Antarctic Islands. Sub-oceanic mantle in tectonically emplaced and abyssal settings is also considered where known. This is complemented by a continental-scale mantle xenolith overview, mantle characteristics from igneous rocks and a quantitative mantle fabric study. State-of-the-art, continental-scale geophysical overviews of the Antarctic mantle are presented by discipline as seismology, gravity and magnetics, magnetotellurics, rheology, glacial isostatic adjustment, mantle convection and palaeotopography. This Memoir will be the reference for all researchers interested in the Antarctic mantle and its role in dynamics that shape the Antarctic surface and ice sheets.

SUMMARY Significant land uplift and horizontal motions have been recorded with Global Navigation Satellite Systems (GNSS) in areas such as Alaska, Iceland and the Northern Antarctic Peninsula (NAP) as a result of Glacial Isostatic Adjustment (GIA) due to ice melt after the Little Ice Age. Here, analysis of horizontal displacement rates can be of extra importance, as they are more sensitive to Earth properties in shallower layers than vertical displacement rates. Proper modelling of horizontal displacement rates with dedicated GIA models requires a spherical Earth with compressible rheology. However, in these small areas, the used GIA models are often incompressible using a Cartesian geometry to ease computation and in some cases allow for lateral viscosity changes or more complex rheology. We investigate the validity of modelled horizontal displacement rates using different approximations, that is using spherical or Cartesian Earth structures, and incompressible, material compressible or compressible rheology. Although the lack of self-gravity and sphericity compensate each other in the vertical, this is less the case for the horizontal. For a disc ice sheet with a radius just over 200 km and a thickness of 1000 m, differences due to sphericity are minimal and the modelled horizontal displacement rates of compressible Cartesian models differ from those simulated by a compressible spherical model by 0.63 mm a−1. Thus, compressible Cartesian GIA models can be applied for modelling horizontal displacement rates of small ice sheets like those in Alaska, Iceland and NAP. Unfortunately, the implementation of compressibility in Abaqus that we use here cannot be extended to spherical models as gravity can not be specified for a spherical body. Other modelling approaches are recommended in such cases.

SUMMARYThis paper presents a method that modifies commercial engineering-oriented finite element packages for the modelling of Glacial Isostatic Adjustment (GIA) on a self-gravitating, compressible and spherical Earth with 3-D structures. The approach, called the iterative finite element body and surface force (FEMIBSF) approach, solves the equilibrium equation for deformation using the ABAQUS finite element package and calculates potential perturbation consistently with finite element theory, avoiding the use of spherical harmonics. The key to this approach lies in computing the mean external body forces for each finite element within the Earth and pressure on Earth's surface and core–mantle boundary (CMB). These quantities, which drive the deformation and stress perturbation of GIA but are not included in the equation of motion of commercial finite element packages, are implemented therein. The method also demonstrates how to calculate degree-1 deformation directly in the spatial domain and Earth-load system for GIA models. To validate the FEMIBSF method, loading Love numbers (LLNs) for homogeneous and layered earth models are calculated and compared with three independent GIA methodologies: the normal-mode method, the iterative body force method and the spectral-finite element method. Results show that the FEMIBSF method can accurately reproduce the unstable modes for the homogeneous compressible model and agree reasonably well with the Love number results from other methods. It is found that the accuracy of the FEMIBSF method increases with higher resolution, but a non-conformal mesh should be avoided due to creating the so-called hanging nodes. The role of a potential force at the CMB is also studied and found to only affect the long-wavelength surface potential perturbation and deformation in the viscous time regime. In conclusion, the FEMIBSF method is ready for use in realistic GIA studies, with modelled vertical and horizontal displacement rates in a disc load case showing agreement with other two GIA methods within the uncertainty level of GNSS measurements.