G3.3

Interacting feedbacks play an important role in governing the stability of the Greenland Ice Sheet under global warming. Here we study the interaction between the positive melt-elevation feedback and the negative feedback from glacial isostatic adjustment (GIA), and how they affect the ice volume of the Greenland Ice Sheet on long time scales. We therefore use the Parallel Ice Sheet Model (PISM) coupled to a simple solid Earth model (Lingle-Clark) in idealized step-warming experiments. Our results suggest that for warming levels above 2°C, Greenland could become essentially ice-free on the long-term, mainly as a result of surface melting and acceleration of ice flow. The negative GIA feedback can mitigate ice losses and promote a partial recovery of the ice volume.

Exploring the full factorial parameter space which determines the relative strength of the two feedbacks reveals that four distinct dynamic regimes are possible: from stabilization, via recovery and self-sustained oscillations to the irreversible collapse of the Greenland Ice Sheet. In the recovery regime an initial ice loss is reversed and the ice volume stabilized at 61-93% of the present day volume. For certain combinations of temperature increase, atmospheric lapse rate and Earth mantle viscosity, the interaction of the GIA feedback and the melt-elevation feedback leads to self-sustained, long-term oscillations in ice-sheet volume with oscillation periods of tens to hundreds of thousands of years and oscillation amplitudes between 15-70% of present-day ice volume. This oscillatory regime reveals a possible mode of internal climatic variability in the Earth system on time scales on the order of 100,000 years that may be excited by or synchronized with orbital forcing or interact with glacial cycles and other slow modes of variability.

How to cite: Zeitz, M., Haacker, J. M., Donges, J. F., Albrecht, T., and Winkelmann, R.: Interacting melt-elevation and glacial isostatic adjustment feedbacks allow for distinct dynamic regimes of the Greenland Ice Sheet, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-9968, https://doi.org/10.5194/egusphere-egu22-9968, 2022.

Recent studies suggest the hotspot currently under Iceland was located beneath eastern Greenland at ~40 Ma BP and that the upwelling of hot material from the Iceland plume towards Greenland is ongoing. A warm upper mantle has a low viscosity, which in turn causes the solid Earth to rebound much faster to deglaciation. In the area of the Kangerlussuaq glacier, a large GPS velocities residual after removing predicted purely elastic deformations caused by present-day ice loss suggests the possibility of such fast rebound to little ice age (LIA) deglaciation. Here we investigate the lithospheric thickness and the mantle viscosity structure beneath SE-Greenland by means of model predictions of solid Earth deformation driven by a low viscosity mantle excited by the LIA deglaciation to the present day. From the comparison of such modeled deformations with the GPS residual, we conclude that 1) a rather thick lithosphere is preferred (90-100 km) 2) and the upper mantle most likely has a viscosity that changes with depth. Assuming a two layer upper mantle, it is not well constrained which part of the upper mantle has to be low, with a preference for low viscosity in the deeper upper mantle.

To understand such results we implemented forward modelling with more realistic earth models, relying on improvements in seismic models, petrology and gravity data. This yields 3D viscosity maps that can be compared to inferences based on the 1D model and forms the basis for 3D GIA models. The conclusion based on the 1D model can be explained with 3D Earth models. In the area of the Kangerlussuaq glacier the seismic derived viscosities prefer a higher viscosity layer above a lower viscosity one. This stems from the slow decrease in viscosity with depth. The layer that is characterized as shallow upper mantle still contains shallow regions with low temperatures, while the deeper upper mantle reaches low viscosities. Generally, for GIA earth models the “higher above lower” viscosity layering is unusual. However, the analysis of the 1D model clearly shows this to be one of the preferred model regions, in combination with a large lithosphere thickness of 100 km. This is a notable result that draws attention to the importance of shallow layering in GIA models. 

How to cite: Barletta, V. R., van der Wal, W., Bordoni, A., and Khan, S. A.: Effect of Icelandic hotspot on Mantle viscosity in southeast Greenland, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-10884, https://doi.org/10.5194/egusphere-egu22-10884, 2022.

The Earth’s mass redistribution due to deglaciation and recent ice sheet melting causes changes in the Earth’s gravity field and vertical land motion in Greenland. The changes are because of ongoing mass redistribution and related elastic (on a short time scale) and viscoelastic (on time scales of a few thousands of years) responses. These signatures can be used to determine the mantle viscosity. In this study, we infer the mantle viscosity associated with the glacial isostatic adjustment (GIA) and long-wavelength geoid beneath the Greenland lithosphere. The viscosity is determined based on a spatio-spectral analysis of the Earth’s gravity field and the land uplift rate in order to find the GIA-related gravity field. We used and evaluated different land uplift data, i.e. the vertical land motions obtained by the Greenland Global Positioning System (GPS) Network (GNET), GRACE and Glacial Isostatic Adjustment (GIA) data. In addition, a  combined land uplift rate using the Kalman filtering technique is presented in this study. We extract the GIA-related gravity signals by filtering the other effects due to the deeper masses i.e. core-mantle (related to long-wavelengths) and topography (related to short-wavelengths). To do this, we applied correlation analysis to detect the best harmonic window. Finally, the mantle viscosity using the obtained GIA-related gravity field is estimated. Using different land uplift rates, one can obtain different GIA-related gravity fields. For example, different harmonic windows were obtained by employing different land uplift datasets, e.g. the truncated geoid model with a harmonic window between degrees 10 to 39 and 10 to 25 showed a maximum correlation with the GIA model ICE-6G (VM5a) and the combined land uplift rates, respectively. As shown in this study, the mantle viscosities of 1.6×10²² Pa s and 0.9×10²² Pa s for a depth of 200  to 650  km are obtained using ICE-6G (VM5a) model and the combined land uplift model, respectively, and the GIA-related gravity potential signal.

How to cite: Bagherbandi, M., Amin, H., Wang, L., and Shirazian, M.: Mantle viscosity derived from geoid and different land uplift data in Greenland, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-13323, https://doi.org/10.5194/egusphere-egu22-13323, 2022.

The western Russian Arctic was partially covered by the Eurasian ice sheet complex during the Last Glacial Maximum (~26 ka BP) and is a focus area for Glacial Isostatic Adjustment (GIA) studies. However, there have been few GIA studies conducted in the Russian Arctic due to the lack of high quality deglacial relative sea-level (RSL) data. Recently, Baranskaya et al. (2018) released a quality-controlled deglacial RSL database for the Russian Arctic that consists of ~400 sea-level index points and ~250 marine and terrestrial limiting data that constrain RSL since 20 ka BP. Here, we use the RSL database to constrain the 3D Earth structure beneath the Russian Arctic, with consideration of the uncertainty in ice model ICE-7G_NA, which is assessed via iteratively refining the ice model with fixed 1D Earth model to achieve a best fit with the RSL data. Also, the uncertainties in 3D Earth parameters and RSL predictions are investigated.

We find an optimal 3D Earth model (Vis3D) improves the fit with the deglacial RSL data compared with the VM7 1D model when fixed with the ICE-7G_NA ice model. Similarly, we show improved fit in the White Sea area, where 1D model shows notable misfits, with the refined ice model ICE-7G_WSR when fixed with VM7 Earth model. The comparable fits of ICE-7G_NA (Vis3D) and ICE-7G_WSR (VM7) implies that the uncertainty in the ice model might be improperly mapped into 3D viscosity structure when a fixed ice model is employed. Furthermore, fixed with refined ice model ICE-7G_WSR, we find an optimal 3D Earth model (Vis3D_R), which fits better than ICE-7G_WSR (VM7), and the magnitude of lateral heterogeneity decreases significantly from Vis3D to Vis3D_R.  We conclude that uncertainty in the ice model needs to be considered in 3D GIA studies.

How to cite: Li, T., Peltier, W. R., Stuhne, G., Khan, N., Baranskaya, A., Shaw, T., Wu, P., and Horton, B.: The inclusion of ice model uncertainty in 3D Glacial Isostatic Adjustment modelling: a case study from the Russian Arctic, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-852, https://doi.org/10.5194/egusphere-egu22-852, 2022.

New research aims to improve relative sea-level (RSL) projections for the Norwegian coast. The main objectives are to: i) collect observations of past RSL changes, ranging from the end of the last ice age to the last century, ii) develop a high-quality database of post-glacial sea-level index points (SLIPs) for the Norwegian coast, and to iii) improve our understanding of past and future vertical land motion using glacial isostatic adjustment (GIA) modelling. To now, our collection of new empirical data has focussed on three significant, but enigmatic RSL histories that are not adequately reproduced in existing GIA models: very recent stillstands and transgressions documented by historical tide gauge records, rapid transgressions during the early- to mid-Holocene Tapes period, and abrupt transgressions during the latest Pleistocene Younger Dryas chronozone. Ongoing field sampling is focussed on developing high-resolution RSL trends from salt marshes, isolation basins, and raised beaches, using multiple biostratigraphic and geochemical proxies (i.e. micropaleontology, macrofossils, x-ray fluorescence, C/N) and dating techniques (i.e. Pb-210, Cs-137, C-14, tephrochronology, geochemical markers). Results from various localities spanning the Norwegian coast provide robust constraints for the timing and rate of RSL change during the Younger Dryas and Tapes chronozones. Additional results providing new estimates of very recent RSL trends in southwest Norway are presented by Holthuis et al. (Late Holocene sea-level change and storms in southwestern Norway based on new data from intertidal basins and salt marshes; Session CL5.2.2). These new and emerging constraints are being integrated into a post-glacial RSL database that incorporates high-quality data from the entire Norwegian coastline. Over 1000 SLIPs have been assembled from published studies. These existing data were updated using current radiocarbon calibration curves, high-resolution digital elevation models, new field observations, and new quantitative estimates of relevant uncertainties. Ongoing GIA modelling is utilizing the new RSL database, a glaciological model that freely simulates ice sheet changes, as well as geodetic and ice margin chronology constraints, to develop rigorous uncertainty estimates for present and future GIA along the Norwegian coast.

How to cite: Lakeman, T. R., Nixon, F. C., Romundset, A., Simpson, M. J. R., Svendsen, J. I., Vasskog, K., Bondevik, S., Milne, G., and Tarasov, L.: Improving past and future relative sea-level constraints for the Norwegian coast, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-12689, https://doi.org/10.5194/egusphere-egu22-12689, 2022.

How to cite: Pallisgaard-Olesen, G., Pedersen, V. K., Gomez, N., and Mitrovica, J. X.: Sea level response to Quaternary erosion and deposition in Scandinavia, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-4969, https://doi.org/10.5194/egusphere-egu22-4969, 2022.

Subsidence is a land use problem in the western and northern Netherlands, especially where both shallow soft soil subsidence and deeper subsidence components, including glacio-isostatic adjustment (GIA), add up. The aim of this study is to improve the estimation of the GIA component within the total subsidence signal across the Netherlands during the Holocene, using coastal plain paleo-water level markers. Throughout the Holocene, the GIA induced subsidence in the Netherlands has been spatially and temporally variant, as shown by previous studies that used GIA modelling and geological relative sea-level rise reconstructions. From the latter work, many field data points are available based on radiocarbon dated coastal basal peats of different age and vertical position. These reveal Holocene relative sea-level rise to have been strongest in the Wadden Sea in the Northern Netherlands. This matches post-glacial GIA subsidence (forebulge collapse) as modelled for the Southern North Sea, being located in the near-field of Scandinavian and British former ice masses.

In this study, geological data analysis of RSL and other paleo-water level data available from the Dutch coastal plain for the Holocene period is considered in addition. The analysis takes the form of designing and executing a 3D interpolation (kriging with a trend: KT), where paleo-water level Z_(x,y,age) is predicted and the field-data points are the observations (Age, X, Y and Z as knowns). We use a spatio-temporal 3D grid that covers the Dutch coastal plain, and reproduces and unifies earlier constructed sea level curves and high-resolution sampled individual sites (e.g. Rotterdam). The function describing the trend part of the interpolation separates linear and non-linear components of relative water level rise, i.e.: long-term background subsidence and shorter-term GIA subsidence signal and postglacial water level rise. The kriging part then processes remaining subregional patterns. The combined reconstruction thus yields a spatially continuous parameterization of regional trends that (i) allows to separate subsidence from water level rise terms, and (ii) is produced independently of GIA modelling to enable cross-comparison. Results are presented for the coastal plain of the Netherlands ([SW] Zeeland – Rotterdam – Holland – Wadden Sea – Groningen [NE]). The percentage of the total coastal-prism accommodation space that appears due to subsidence, from the south to the north of the study area increases by 20%. Holocene-averaged subsidence rates from the first analysis ranged from ca. 0.1 m/kyr (Zeeland) to 0.4 m/kyr (Groningen), which is 5-10 times larger than present-day GPS/GNSS-measured rates.

The research presented in this abstract is part of the project Living on soft soils: subsidence and society (grantnr.: NWA.1160.18.259). 

How to cite: de Wit, K., van de Wal, R. S. W., and Cohen, K. M.: Reconstructing large scale differential subsidence in the Netherlands using a spatio-temporal 3D paleo-groundwater level interpolation, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-8350, https://doi.org/10.5194/egusphere-egu22-8350, 2022.

A challenge to understanding Late Quaternary glaciation history is the mechanism(s) responsible for the asymmetry in an individual glaciation cycle between the slow pace of glaciation and the more rapid pace of deglaciation (e.g., Broecker and Van Donk, 1970). It is increasingly clear that a major contributor to the rate of global deglaciation is the instability of marine terminating ice streams. Recent analyses by Velay-Vitow et al. (2020) suggest that these instabilities were often triggered by ocean tides of anomalously high amplitude. Examples include the Hudson Strait Ice Stream responsible for Heinrich Event 1 (H1) and the Amundsen Gulf Ice Stream. Here, we analyse the instability of the Laurentian Channel and St Lawrence River Channel ice stream system. Our analysis begins with the recognition of highly significant misfits of up to 60 m at ~9,000 calendar years ago between deglacial relative sea-level histories inferred by Vacchi et al. (2018) at sites along the St Lawrence River Channel and those predicted by the ICE-6G_C (VM5a) and ICE-7G_NA (VM7) models of the Glacial Isostatic Adjustment process. We suggest that these disagreements between models and data may be due to the St Lawrence River Channel ice stream becoming unstable during the deglaciation of the Laurentide Ice Sheet (LIS) due to the hypothesized tidal mechanism for ice stream destabilization. We investigate a sequence of scenarios designed to provide a best estimate of the timing of this event. Since this ice stream penetrated deeply into the interior of the LIS and was connected to the Laurentian Channel ice stream, the instability of the latter was required in order for destabilization of the St Lawrence River channel ice stream to be possible. We explore the consistency of the implied sequence of events with the observational constraints.

How to cite: Peltier, R., Li, T., Stuhnne, G., Velay-Vitow, J., Vacchi, M., Englehart, S., and Horton, B.: Resolving the Influence of Ice Stream Instability on Postglacial Relative Sea-Level Histories: the case of the St Lawrence River Channel Ice Stream, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-1343, https://doi.org/10.5194/egusphere-egu22-1343, 2022.

Many coastal cities are an early casualty in climate-related coastal flooding because of processes resulting in land subsidence and thus enhanced relative sea-level (RSL) rise. Much of the Atlantic coast of North America has been sinking for thousands of years, at a maximum rate of ~20 cm per century as a consequence of solid Earth deformation in response to deglaciation of the Laurentide ice sheet (between ~18,000 and ~7,000 years ago) [e.g. Love et al., Earth's Future, 4(10), 2016]. Karegar et al. [Geophysical Research Letters, 43(7), 2016] have shown that vertical land motion along the Atlantic coast of the USA is an important control on nuisance flooding. A key finding in this study is that while glacial isostatic adjustment (GIA) is the dominant process driving land subsidence in most areas, there can be large deviations from this signal due to the influence of anthropogenic activity impacting hydrological processes. For example, between Maine (45°N) and New Hampshire (43°N), the GPS data show uplift while geological data show long-term subsidence. The cause of this discrepancy is not clear, but one hypothesis is increasing water mass associated with the James Bay Hydroelectric Project in Quebec [Karegar et al., Scientific Reports, 7, 2017].

The primary aim of this study is to better constrain and understand the processes that contribute to contemporary and future vertical land motion in this region to produce improved projections of mean sea-level change and nuisance flooding. The first step towards achieving these aims is to determine a GIA model parameter set that is compatible with observations of past sea-level change for this region. We make use of two regional RSL data compilations: Engelhart and Horton [Quaternary Science Reviews, 54, 2012] for northern USA and Vacchi et al. [Quaternary Science Reviews, 201, 2018] for Eastern Canada, comprising a total of 1013 data points (i.e., sea level index points and limiting data points) over 38 regions distributed throughout our study region. These data are well suited to determine optimal GIA model parameters due to the magnitude of other signals being much smaller, particularly in near-field regions such as Eastern Canada. We consider a suite of 32 ice history models that is comprised mainly of a subset from Tarasov et al. [Earth and Planetary Science Letters, 315–316, 2012] as well as the ICE-6G and ANU models. We have computed RSL for these ice histories using a state-of-the-art sea-level calculator and 440 1-D Earth viscosity models per each ice history model to identify a set of Earth model parameters that is compatible with the observations.

How to cite: Parang, S., Milne, G. A., Karegar, M. A., and Tarasov, L.: Towards an improved understanding of vertical land motion and sea-level change in eastern North America, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-50, https://doi.org/10.5194/egusphere-egu22-50, 2022.

Earth’s topography and bathymetry is shaped by a complex interplay between solid-Earth processes that deform the Earth from within and the surface processes that modify the outer shape of the Earth. At the surface, an ultimate baselevel set by global sea level marks the defining transition from erosion to deposition. Over geological time scales, this baselevel has resulted in a distinct hypsometric distribution (distribution of surface area with elevation), with the highest concentration of surface area focused in a narrow elevation range near present-day sea level.

This particular feature in Earth’s hypsometry makes the global land fraction very sensitive to changes in sea level. Indeed, a sea-level change will result in a significant change in the land fraction as dictated by the hypsometric distribution, thereby modulating the very same sea-level change. However, it remains unexplored exactly how sea-level changes have modified the global land fraction over past glacial cycles and into the future.

Here we analyse how Earth’s hypsometry has changed over the last glacial cycle as large ice sheets waxed and waned particularly in Scandinavia and North America. These changes in global ice volume resulted in a significant global excursion in sea level, modulated regionally by solid-Earth deformation, gravitational effects, and effects from Earth’s rotation. These changes modified Earth’s hypsometry, and therefore the global land fraction at any given time. Consequently, we can map out how Earth’s hypsometry has influenced global mean sea level (GMSL) over time. To examine this relationship between Earth’s hypsometry and sea level further, we look to the deep future, to a scenario where both the Greenland Ice Sheet and the Antarctic Ice Sheets will melt away completely over multi-millennial timescales. This scenario is not meant to represent a realistic future scenario per se, but it allows us to define the hypsometric GMSL correction needed for any GMSL that the Earth has experienced recently or will experience in the future.

How to cite: Pedersen, V. K., Gomez, N., Pallisgaard-Olesen, G., Garbe, J., Aschwanden, A., Winkelmann, R., and Mitrovica, J.: The influence of Earth’s hypsometry on global sea level through a glacial cycle and into the future, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-11569, https://doi.org/10.5194/egusphere-egu22-11569, 2022.

GIA is a global process, because of gravitational effects, its interplay with earth rotation, and the large spatial extent of ice-sheet and ocean loading. However, mainly due to the presence of heterogeneities in the structure of crust and upper mantle, modelling of GIA signals often requires a regional approach. This is particularly true in the light of continuous advances in earth observation techniques, that allow increasingly accurate determination of land deformation, coastal sea level change, and mass balance of glaciers and ice sheets.

This talk will address a number of open issues related to regional GIA models, such as the effect of transient and non-linear rheologies, and the complementary role of forward and semi-empirical approaches, with an eye on the needs of the geodetic, sea level and cryosphere communities.

How to cite: Riva, R.: Regional GIA: modelling choices and community needs, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-918, https://doi.org/10.5194/egusphere-egu22-918, 2022.

A new finite element method called FEMIBSF is presented that is capable of modelling Glacial Isostatic Adjustment (GIA) on compressible earth models with three-dimensional (3D) structures. This method takes advantage of the classical finite element techniques to calculate the deformational and gravitational responses to the driving forces of GIA (including body forces and pressures on Earth’s surface and core-mantle boundary, namely CMB). Following Wu (2004) and Wong & Wu (2019), we implement the GIA driving forces in the commercial finite element software Abaqus and solve the equation of motion in an iterative manner. Different from those two studies, all formulations and calculations in this study are not associated with spherical harmonics but are performed in the spatial domain. Due to this, FEMIBSF is free from expanding the load, displacement, and potential into spherical harmonics with the short-wavelength components (of high degree and order) neglected. We compare the loading Love numbers (LLNs) generated by FEMIBSF with their analytical solutions for homogeneous models and numerical solutions for layered models calculated by the normal-mode approach/code, ICEAGE (Kaufmann, 2004), the iterative body force approach/code, IBF (Wong & Wu, 2019) and the spectral-finite element approach/code, VILMA-C (Martinec, 2000; Tanaka et al., 2011). We find that FEMIBSF agrees well with analytical and numerical LLN results of these codes. In addition, we show how to compute the degree-1 deformation directly in the spatial domain with the finite element approach and how to implement it in a GIA model using Abaqus. Finally, we demonstrate that the CMB pressure related to the gravitational potential change in the fluid core only influences the long-wavelength surface displacement and potential such as the degree-2 component.

References

Kaufmann, G. (2004). Program Package ICEAGE, Version 2004. Manuscript. Institut für Geophysik der Universität Göttingen.

Martinec, Z. (2000). Spectral–finite element approach to three-dimensional viscoelastic relaxation in a spherical earth. Geophysical Journal International, 142(1), 117-141.

Tanaka, Y., Klemann, V., Martinec, Z. & Riva, R. E. M. (2011). Spectral-finite element approach to viscoelastic relaxation in a spherical compressible Earth: application to GIA modelling. Geophysical Journal International, 184(1), 220-234.

Wong, M. C. & Wu, P. (2019). Using commercial finite-element packages for the study of Glacial Isostatic Adjustment on a compressible self-gravitating spherical earth–1: harmonic loads. Geophysical Journal International, 217(3), 1798-1820.

Wu, P. (2004). Using commercial finite element packages for the study of earth deformations, sea levels and the state of stress. Geophysical Journal International, 158(2), 401-408.

How to cite: Huang, P., Steffen, R., Steffen, H., Klemann, V., van der Wal, W., Reusen, J., Tanaka, Y., Martinec, Z., and Thomas, M.: A finite element approach to modelling Glacial Isostatic Adjustment on three-dimensional compressible earth models, EGU General Assembly 2022, Vienna, Austria, 23–27 May 2022, EGU22-6013, https://doi.org/10.5194/egusphere-egu22-6013, 2022.

Uncertainty in present-day glacial isostatic adjustment (GIA) rates represent at least 44% of the total gravity-based ice mass balance signal over Antarctica. Meanwhile, physical couplings between solid Earth, sea level and ice dynamics enhance the dependency of the spatiotemporally varying GIA signal on 3D rheology. For example, the presence of low-viscosity mantle beneath melting marine-based ice sheet sectors such as the Amundsen Sea Embayment may delay or even prevent unstable grounding line retreat. Improved knowledge of upper mantle thermomechanical structure is therefore required to refine estimates of current and projected ice mass balance.

Here, we present a Bayesian inverse method for mapping shear wave velocities from high-resolution adjoint tomography into thermomechanical structure using a calibrated parameterisation of anelasticity at seismic frequency. We constrain the model using regional geophysical data sets containing information on upper mantle temperature, attenuation and viscosity structure. The Globally Adaptive Scaling Within Adaptive Metropolis (GASWAM) modification of the Metropolis-Hastings algorithm is utilised to allow efficient exploration of the multi-dimensional parameter space. Our treatment allows formal quantification of parameter covariances, and naturally permits us to propagate uncertainties in material parameters into uncertainty in thermomechanical structure.

We find that it is possible to improve agreement on steady state viscosity structure between tomographic models by approximately 30%, and reduce its uncertainty by an order of magnitude as compared to a forward-modelling approach. Direct access to temperature structure allows us to estimate lateral variations in lithospheric thickness, geothermal heat flow, and their associated uncertainties.