the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 15 Sep 2025
Automated forward and adjoint modelling of viscoelastic deformation of the solid Earth
Abstract. Robust models of viscoelastic Earth deformation under evolving surface loads underscore many problems in geodynamics and are particularly critical for paleoclimate and sea-level studies through their role in Glacial Isostatic Adjustment (GIA). A long-standing challenge in GIA research is to perform computationally efficient inversions for ice-loading histories and mantle structure using a physically realistic Earth model that incorporates three-dimensional viscosity variations and/or complex rheologies. For example, recent geodetic observations from melting ice sheets appear inconsistent with long-term sea-level records and have been used to argue for transient rheologies, generating debate in the literature and leaving large uncertainties in projections of future sea-level change. Here, we extend the applicability of G-ADOPT (a Firedrake-based finite element framework for geoscientific adjoint optimisation) to these problems. Our implementation solves the equations governing viscoelastic surface loading while naturally accommodating elastic compressibility, lateral viscosity variations, and non-Maxwell rheologies (including transience). We benchmark the approach against a suite of analytical and numerical test cases, demonstrating both accuracy and computational efficiency. Crucially, G-ADOPT enables automatic derivation of adjoint sensitivity kernels, allowing gradient-based optimisation strategies that are essential for high-dimensional inverse problems. Using synthetic Earth-like experiments, we illustrate its capability to reconstruct ice histories and recover mantle viscosity variations, providing a roadmap towards data assimilation and uncertainty quantification in GIA modelling and sea-level projections.
- Preprint
(6956 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
- RC1: 'Comment on egusphere-2025-4168', Anonymous Referee #1, 25 Nov 2025
-
RC2: 'Comment on egusphere-2025-4168', Anonymous Referee #2, 08 Dec 2025
This review considers “Automated forward and adjoint modelling of viscoelastic deformation of the solid Earth” submitted to Geoscientific Model Development by Scott et al. The manuscript addresses a new computational approach for joint forward modeling of viscoelastic solid Earth deformation and inversion of Earth structure and/or loading histories using adjoint methods with the existing open-source G-ADOPT numerical framework.
Overall, this is an excellent manuscript that is well suited for the publication in Geoscientific Model Development. The authors thoroughly describe the governing equations, numerical methods, and demonstrate how the forward modeling component performs for relevant benchmark mark cases. I particularly appreciated the authors’ thoroughness in describing the effect of distinct discretizations (element types) and other numerical parameters on the benchmark results. Although it is outside my area of expertise, the demonstrated joint inversion for the ice and viscosity structure looks highly promising. The modeling framework and workflows clearly represent a notable advance that will enable the GIA community to address new questions and/or provide alternative approaches towards estimating uncertainties for both time-dependent glacial loads and present-day earth structure.
I think the manuscript is effectively ready to be published in its present form, but may benefit from first addressing the following points in the “Discussion, limitations and future work” section:
- My suspicion is that down the line many key problems of interest may require resolutions that necessitate scaling beyond 6,000 cores. What changes to the numerical framework would be required to scale to 10,000s of CPUs? A switch from an algebraic to geometric (matrix-free) multigrid preconditioner? The paper would benefit from briefly describing how this could potentially be achieved within the G-ADOPT framework.
- While the numerical framework is capable of simulating models with large viscosity variations, the underlying constitutive model is for a linear viscoelastic material. What possibility is there to also integrate nonlinear constitutive models, such as the viscoelastic-plastic models commonly used in the tectonic geodynamics community and in Shijie Zhong’s codes for modeling loading induced deformation (e.g., Zhong and Watts, 2013, https://doi.org/10.1002/2013JB010408)? There is a note on line 165 that the constitutive formulation can be extended to non-linear constitutive equations, but the manuscript would benefit from more detail on exactly what can be reasonably achieved within the current framework within the context of other relevant studies.
Aside from the points, I think the manuscript is in excellent shape and only requires at most minor revisions before publication.
Citation: https://doi.org/10.5194/egusphere-2025-4168-RC2
Model code and software
G-ADOPT code and runscripts associated with 'Automated forward and adjoint modelling of viscoelastic deformation of the solid Earth' William Scott https://zenodo.org/records/16925271
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 1,114 | 103 | 19 | 1,236 | 19 | 22 |
- HTML: 1,114
- PDF: 103
- XML: 19
- Total: 1,236
- BibTeX: 19
- EndNote: 22
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
This manuscript presents a strong and scientifically valuable contribution by extending the G-ADOPT inversion framework to transient viscoelastic surface-loading problems. The authors demonstrate how the method can be applied to ice-loading scenarios with increasing rheological complexity and provide clear validation through comparisons with existing numerical codes and (semi-)analytical solutions. Numerical stability and scalability are also thoroughly assessed, giving confidence in the robustness of the implementation.
The core advances lie in the development of adjoint-based inversions for ice thickness and lateral viscosity variations (LVVs) within a 2D annular geometry. Automatically derived discrete adjoint equations are used to compute gradients, and a synthetic dataset of surface displacements through time serves as observational input for the optimisation. The inversion for the initial ice load performs impressively, recovering the synthetic distribution with high fidelity. The subsequent inversion for LVVs is likewise promising: low-viscosity regions are well resolved, and high-viscosity anomalies are retrieved with correct magnitudes and positions, though geometric details are less sharply defined. Finally, a combined inversion for both ice thickness and LVVs highlights the challenge of constructing a cost function that retains sensitivity to multiple parameter fields that reside on different meshes. Despite this complexity, the inversion still performs convincingly and effectively demonstrates the capabilities of the framework. The discussion also provides a thoughtful roadmap for future work, including the need to move beyond the controlled synthetic setting and confront challenges associated with real-world data.
Overall, the manuscript is well written, technically sound, and highly relevant to the geoscientific modelling community. The methodology is clear, and the motivations behind modelling decisions are well articulated. A few clarifications, however, would improve accessibility, especially for readers less familiar with adjoint-based gradient computation and inverse methods:
• Line 586: Regularisation is noted as an important mechanism for stabilizing the adjoint solution (see also line 688), particularly when observations are indirect outputs rather than primary solution fields. In the LVV inversion, the gradient appears significantly steeper across the upper–lower mantle viscosity transition due to the sharp contrast. Would introducing spatial or parameter-space regularisation help improve the reconstruction of features that span this boundary? If so, could the authors comment on how such regularisation might be formulated within the objective function without excessively smoothing meaningful structures?
• Line 685: It would be helpful to clarify how observational uncertainty is intended to be handled. Will observational variance be incorporated into the cost function weighting during inversion, or is the aim instead to estimate posterior uncertainty bounds on the inferred parameters—or both? Clarifying this point would help readers understand the intended interpretational framework of the inversions.
In conclusion, this is an excellent and meaningful contribution that advances adjoint-based geodynamic inversion and provides a solid foundation for future work toward applications involving real observational datasets. I believe the manuscript is suitable for publication after minor revisions addressing the points above.