| 22 Sep 2025
Inferring subglacial topography using physics informed machine learning constrained by two conservation laws
Abstract. Subglacial topography beneath the Greenland Ice Sheet is a fundamental control on its dynamics and response to changes in the climate system. Yet, it remains challenging to measure directly, and existing representations of the subglacial topography rely on a limited number of observations. Although the use of mass conservation and the development of BedMachine Greenland substantially improved the representation of the bed topography, this approach is limited to fast-flowing sectors and is less effective in regions with complex, alpine topography. As an alternative to traditional numerical methods, recent work has explored using Physics Informed Neural Networks (PINNs), constrained by only one physical law, to solve forward and inverse problems in ice sheet modeling. Building on this work, we assess three PINN frameworks constrained by distinct conservation laws, showing that PINNs informed with a single conservation law are not sufficient for regions with sparse measurements and complex topographies. To that end, we introduce a novel approach that involves coupling two conservation laws within a PINN framework to infer the subglacial topography and test this approach for three regions with distinct environments in Greenland. This PINN is trained with both the conservation of mass and an approximation of the conservation of momentum (the Shelfy-Stream Approximation), which allows us to simultaneously infer the ice thickness and basal shear stress using observations of ice velocities, surface elevation, and the apparent mass balance in a mixed inversion problem. We compare the predicted ice thickness to ground-truth ice-penetrating radar measurements of ice thickness, showing that the PINN informed with two conservation laws is capable of inferring ice thickness in sparsely surveyed regions. Furthermore, comparisons of predicted bed topographies with BedMachine Greenland show that this approach is capable of discovering new bed features in slower-moving regions and in regions of complex topography, highlighting its potential for better constraining the bed topography of the Greenland Ice Sheet.
Competing interests: One of the authors is a member of the editorial board of journal The Cryosphere.
Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors. Views expressed in the text are those of the authors and do not necessarily reflect the views of the publisher.
=== Summary and general comments ==
This is a well-organized and well written paper describing the application of physics-informed neural networks (PINNS) for further improving the generation of sub-glacial bed topography datasets, building on previous “BedMachine” efforts that have been ongoing for the past decade or so. Two conservation equations – the continuity equation and a momentum balance equation for ice flow – are considered and introduced as additional constraints on the loss function, akin to their introduction as constraints in the “cost function” of PDE-constrained optimization (a good analogy to consider including for readers more familiar with the language of glaciology modeling / optimization?). The two constraints are considered on their own and in combination and the resulting bed topography datasets are compared and contrasted with previous BedMachine results for three glaciologically distinct regions. Overall, the authors argue convincingly that the new approach has merit and demonstrates potential for improving inferred bed topography in regions where the traditional BedMachine approach begins to break down.
Overall, the work seems very worthy of publication and readers of The Cryosphere will find it a worthy contribution. My suggestion would be to accept for publication with minor revisions, noting that these are those suggested revisions (detailed below and identified by their line number in the submitted version) are largely editorial in nature.
My one more substantial suggestion – not necessarily for this publication but possibly for a future effort – is that I think it would be very useful to redo this exercise for a single region (I realize the computational cost could be a challenge, so pick a single region, like the most challenging one discussed here) but using L1L2 (Blatter/Pattyn) for the momentum balance model as opposed to SSA (which you’ve already done and could reuse as the baseline). Because the former allows for internal deformation, a possibly very different vertical velocity profile (and hence depth-averaged velocity and flux divergence) might be implied in regions of slower moving ice (noting that the modeled 2d surface velocity field could still be used as the velocity constraint, so that there should not necessarily be any significant reformulation of the loss functions discussed herein). It would be interesting to see if the more accurate stress balance constraint helped to alleviate any of the remaining problems discussed below.
=== Detailed comments ==
21: Are these refs now considered the definitive source for defining the amount of potential SLR locked up in the ice sheets? If not, maybe consider adding one or two more from other authors for the sake of diversity?
23: “numerical ice sheet modeling”
25: It seems like a summary-level reference might also be appropriate here (?), e.g. something from one of the recent IPCC reports (that integrates results from a large number of individual publications).
28-29: It should probably be noted here that these experiments were assuming a marine ice sheet with a significant over-deepening inland (since this configuration would necessarily be more sensitive than say, an ice sheet grounded above sea level).
37-38: Is the ~2km limit proposed here coming from the Durand et al. (2011) paper? While I more-or-less agree with this idea, I don’t know that this single reference is adequate to support the precision implied by this statement. Maybe consider softening it a little bit to something less precise, e.g. “order km-scale spatial resolution”?
63: “three regions in Greenland”; maybe add a few words of clarification here that they are glaciologically distinct / different? E.g., presumably you mean regions where velocity occurs primarily via fast sliding, a region where it occurs via a mix of sliding and deformation, etc.
Figure 1 caption: “The loss function is comprised of …” or “The loss function includes data loss …”
89: “fully connected layers”, maybe use “fully connected (‘dense’) layers …” ?
101: Should it be “the apparent mass balance residual” ?
103-105: I am guessing that maybe this is discussed further below (?), but it seems like you are already potentially limiting the usefulness of this approach by restricting the momentum balance to SSA. I.e., if one of the main interests here is in improving the inference in regions of slower moving ice flow, which is presumably due to less sliding and more internal deformation, then SSA doesn’t seem like the right assumption to make for the model dynamics. I know that ISSM has higher-order approximations available (e.g., L1L2 or “Blatter-Pattyn”). Has that also been explored (acknowledging the obvious additional computational burden) and compared against the approach using SSA?
116: “… from THE regional climate model RACMO…”
Section 2.3: It sounds like the basal mass balance term in equation 2 is assumed to be ~0? If so, it would be good to note that explicitly here in the discussion of the apparent mass balance term.
152: “to prevent from taking” (omit “from”?); “or diving by zero” (“dividing by zero”)
176-177: Would it be worth commenting on the choice of median vs. mean? Is the median chosen because of the small number (5) of samples, such that the mean could be easily biased?
180: By “challenging to implement”, do you mean where the traditional / previous mass conservation approach does not perform well? Implementation sounds more like the approach is challenging, but I imagine the approach is just as easy to implement in these regions, it’s more the prior / baseline result that you are not happy with.
242: It’s not clear here exactly what “Fig.2(1)” is referring to.
Figure 5: In the caption for this and figure 4 it would be helpful to remind the reader which dataset is subtracted from which and shown in panels d-f (e.g., PINN minus original BedMachine product or vice versa?).
Table 4: It might also be useful to provide some percent / fractional metrics here? E.g., for the apparent mass balance RMSE, how does that number compare to the average apparent mass balance over the same area? Such a table could be added to the SI if it’s not deemed important enough for the main text.
3.2.2. – It’s left hanging a bit as to the significance of the differences in u, apparent mass balance, and sfc. elevation when using the different approaches. For example, how do these differences compare to those that arise when using the original BedMachine approach? Would it make sense to include those metrics (differences in u, apparent mass balance, and sfc elevation) somewhere here for comparison? It’s a bit unclear to me what the broader implications are of these secondary metrics w.r.t. using the derived datasets for modeling. If the authors have additional thoughts on this they would be welcome in the supplementary information.
326: If the discussion starting in 4.1 is intended to be specific to Deception, then maybe that should be noted earlier in this paragraph? Alternatively, if the discussion in 4.1 up to line 326 where Deception is mentioned is supposed to be generic, then perhaps line 326 should be something more like, “… a far more realistic bed topography map, particularly for Deception.”
329-330: Would “…slightly higher RMSE SUGGESTS …” be more appropriate here than “indicates”? I think the speculation in this sentence makes sense, but it seems like it is perhaps speculation as opposed to a concrete fact.
336-348: W.r.t. the prediction of thinner ice – is it also possible that this could be the result of the chosen stress balance model? E.g., in order for the SSA model to match surface velocities, it would need to assume a depth-averaged velocity profile that is larger than would be assumed in a model that allowed for internal deformation (E.g., L1L2), because SSA can only accommodate velocity via a change in the sliding component (unless I’m misunderstanding the model used here). If that is indeed the case, then it seems like the optimization process might necessarily bias the ice thickness on the thin side; if the depth averaged velocity is too large, the same flux (constrained by continuity equation and the apparent mass balance terms) can only be accommodated by reducing the ice thickness.
345: “These reasons imply that …” is a bit awkward. “This implies that …”? “These arguments imply that …” ?
353: “… and exceeds their INDIVIDUAL limitations …” ?
398: “We observe that the PINN better captures…” --> “We observe that the PINN captures observable features better with …”
403: “… state variable predictions”. (remove plural on “variable”)
415: “… mass-conserving approach, AS CONFIRMED BY THE DISCOVERY OF new bed features beneath Narssap and Deception.”
417: “… we recommend USING this approach …”
A last thought / general comment: The implied “geomorphology” of the three focus areas studied here look very different from one another. E.g., The Upernavik and Narsaap beds look very smooth when compared to Deception. In the areas where there are no troughs, they almost look like high-resolution DEMs from past, heavily glaciated regions of Canada. Is there any published work on previous Greenland glaciations that might provide some more insight into this? I’m not suggesting it should be part of this paper, but it could be interesting to look into whether or not the “smoothness” that your methods are implying about the bed in different regions is in line with current glacial geological / geomorphological understanding. It would seemingly be a further testament to the power of the methods used here if you were resolving that level of information about the bed through hundreds / thousands of meters of ice.