Changes in sea level in response to ice-mass loss are spatially variable because of glacial isostatic adjustment (GIA), the deformational, gravitational and rotational response of the Earth to changes in surface ice and water distribution. The response of the bedrock is viscoelastic, beginning with an instantaneous elastic response of the solid planet's lithosphere and mantle and transitioning to a longer timescale viscous relaxation of the mantle towards isostatic equilibrium. GIA models produce predictions of changes in the height of the sea surface equipotential and solid Earth surface (i.e. sea-level changes) in response to surface ice and water loading changes, which are in turn used to interpret satellite, geophysical and geological records and serve as input to models of ice-sheet dynamics.
Accurate GIA modelling is required to constrain the sea-level and solid Earth feedbacks on ice dynamics in the coming centuries, especially along unstable marine-grounded ice fronts in Antarctica where bedrock uplift and gravitational drawdown of the sea surface due to ice loss act as a negative feedback to stabilise the retreat of marine-grounded ice-sheet grounding lines (e.g. Gomez et al., 2010, 2013, 2015; De Boer et al., 2014; Konrad et al., 2015; Larour et al., 2019). Furthermore, the Earth's response to past and modern ice cover changes makes a significant contribution to satellite records of modern mass changes in marine sectors of West Antarctica that are actively experiencing ice loss (e.g. King et al., 2012; the IMBIE team, 2018).
GIA modelling in Antarctica is complicated by the fact that the continent is characterised by strong lateral variability in lithospheric thickness and upper-mantle viscosity, with low viscosities in the west and high viscosities in the east (e.g. An et al., 2015a; Heeszel et al., 2016; Lloyd et al., 2020). In particular, the low-viscosity mantle and thin lithosphere observed under the West Antarctic Ice Sheet (WAIS) identified from increasingly resolved seismic tomography and geodetic and geologic constraints (Ritzwoller et al., 2001; Morelli and Danesi, 2004; Kaufmann et al., 2005; Nield et al., 2014; Heeszel et al., 2016; Barletta et al., 2018; Shen et al., 2018; Lloyd et al., 2020) lead to a more spatially localised (short wavelength) and faster viscoelastic response to surface loading than cratonic regions covered by Late Pleistocene ice sheets (e.g. Hay et al., 2017; Powell et al., 2020). Over West Antarctica, upper-mantle viscosities are thought to vary by several orders of magnitude over short spatial scales (∼ hundreds of kilometres) reaching as low as 1018 Pa s in the Amundsen Sea Embayment (ASE) beneath areas of active marine ice loss (e.g. Nield et al., 2014; Barletta et al., 2018). This implies that viscous effects due to 20th century and more recent ice loss will become significant on annual to decadal timescales and accelerate during the time frame of instrumental records (Barletta et al., 2018; Powell et al., 2020). Viscous deformation due to ongoing ice loss also has the potential to influence ice-sheet grounding lines in the coming centuries (Gomez et al., 2015; Kachuck et al., 2020; Coulson et al., 2021) but has not been included in recent high-resolution coupled ice-sheet-sea-level model projections (Larour et al., 2019).
To accurately capture the timing and wavelength of GIA effects across Antarctica, models must be capable of both accounting for 3-D Earth structure (i.e. 3-D GIA models such as Latychev et al., 2005, or van der Wal et al., 2015) and be of sufficient spatiotemporal resolution to capture the geometry of grounded ice cover (e.g. Han et al., 2022). Commonly, GIA, ice-sheet and coupled ice-sheet-GIA modelling (e.g. Gomez et al., 2015; Konrad et al. 2015) studies that consider modern and future ice-sheet changes have been performed with only 1-D (radially varying) Earth structure models (e.g. Kendall et al., 2005; Spada and Stocchi, 2007; Adhikhari et al., 2016) or with coarse spatial resolutions of > 20 km due to the computational expense. GIA models capable of kilometre to sub-kilometre resolution have also been developed (e.g. the 1-D GIA model by Adhikhari et al., 2016; the 3-D GIA model by Latychev et al., 2005, with updates described in the Supplement of Gomez et al., 2018). For computational efficiency, some of these models implement regional grid refinement techniques which allow for a higher resolution along ice retreat margins. Alongside this, improvements in ice-sheet models (e.g. Nowicki et al., 2016; Seroussi et al., 2020; DeConto et al., 2021) and observational products (e.g. Studinger, 2014; Bamber et al., 2020; Smith et al., 2020; Morlinghem et al., 2020) have been made to provide similarly high-resolution (kilometre to sub-kilometre) ice thickness and bedrock topography datasets that serve as input to GIA models. These advancements together allow for GIA models to capture short-wavelength bedrock deformation and input ice loading changes at unprecedented detail but at still heavy computational expense, particularly for global 3-D GIA models.
It is well-established that dynamic ice-sheet models are sensitive to the chosen grid resolution (e.g. Durand et al., 2009; Van den Berg and Van de Wal, 2006), requiring kilometre to sub-kilometre resolution to accurately represent ice dynamics and grounding-line migration in some applications (e.g. Gladstone et al., 2012; Pattyn et al., 2013; Cornford et al., 2016). It has also been suggested that bedrock topography as fine as 1 km resolution may be needed to capture the influence of fine-scale topographic features on the ice-sheet evolution (Durand et al., 2009), and high resolution may also be needed to represent some embayment walls and pinning points that act to slow down retreat (e.g. Favier et al., 2012; Joughin et al., 2014; Berger et al., 2016).
While topographic features themselves can be very fine scale, the changes in bedrock elevation and sea level in response to ice cover changes tend to be longer wavelength, and the corresponding spatial resolution required to accurately resolve these changes in GIA models and their influence on ice dynamics remains poorly understood. Larour et al. (2019) suggested that kilometre-scale resolution may be required to capture the elastic component of deformation in response to ice loss. However, the idealised tests they performed considered an isolated and increasingly localised load, and their conclusion may not hold for more realistic, spatially coherent ice loss scenarios. Furthermore, their model did not include viscous deformation in response to ongoing ice loss during the simulation or account for lateral variations in Earth structure. There have been limited studies assessing the length scale of realistic viscoelastic bedrock response in the rheologically complex region beneath the WAIS, though a recent high-resolution bedrock deformation modelling study by Zwinger et al. (2020) suggests a convergence in modelled viscoelastic deformation at resolutions finer than 5 km. The broad spatial nature of bedrock deformation and the spatially coherent nature of ice-sheet retreat, which becomes less localised over longer timescales, suggest that sub-kilometre to kilometre grid resolution, which comes at great computational cost, may not be necessary for accurate GIA model calculations.
The aim of this study is to assess the sensitivity of GIA model predictions of Earth deformation and sea-level change in response to modern and future ice loss to spatial resolution in the rheologically complex marine sectors of the WAIS. We build a 3-D viscosity model based on the most recent Antarctic-wide seismic tomography model (Lloyd et al., 2020) to serve as input to a 3-D finite volume, global GIA model (Latychev et al., 2005) to assess the performance of 3-D GIA model predictions across surface grid resolutions of 1.9-15 km. We repeat calculations with a range of Earth models, considering the contribution from elastic and viscous deformation separately. We focus on the response to observed modern ice loss over the last two decades (Shepherd et al., 2019) and projected future ice-sheet retreat in the coming century (Golledge et al., 2019; DeConto et al., 2021) in the Amundsen Sea Embayment of West Antarctica. Our study is motivated by the following questions. What 3-D grid resolution is necessary to adequately capture the elastic and viscous deformation and sea-level changes in response to ice loading changes? How significant is the effect of grid resolution compared to sources of uncertainty and simplifications made in some previous modelling, in particular the neglect of viscous deformation?
