A recent paper by Reference CalovCalov and others (2010) presents an intercomparison of several numerical ice-sheet models where the development of large-scale thermo-viscous instabilities is compared. In all but one of the models, the mechanics are approximated using the shallow-ice approximation (SIA). My purpose in writing this letter is to point out the dangers of using this approximation in this kind of analysis.

Modelling any system that contains physical instabilities immediately invokes the issue of whether the results are real or numerical artefacts. This issue hit glaciology in the 1990s, when ice-sheet modellers were seeking to understand the formation of ice streams and the possibility of ice-stream surges. One postulated mechanism was thermo-viscous instability (Reference PaynePayne, 1995). This arose from the fact that ice softens as it warms, meaning that any frictional heating can allow the ice to flow faster, leading to greater dissipative heating, engendering a positive feedback. Researchers using the simplest mechanical approximation, the SIA, produced some spectacular and broadly consistent results when modelling flow in two dimensions (Reference PaynePayne, 1995; Reference Greve and MacAyealGreve and MacAyeal, 1996; Reference PattynPattyn, 1996). However, once workers moved to three dimensions (3-D), while the models appeared to be able to spontaneously develop ice streams – the fast-flowing corridors of ice observed in nature – it was seen that results were inconsistent between workers even in terms of the number of streams generated (Reference PaynePayne and others, 2000). This was a serious obstacle to developing ice-sheet models that could explain, for example, ice-sheet surges, which are believed to have played a significant role in modulating palaeoclimate.

A further investigation into this is found in Reference CalovCalov and others (2010), where several models carry out an inter-comparison of thermomechanically coupled models, using parameter settings that encourage thermo-viscous instabilities. They find instabilities developing, stating

1. 1. usage of non-symmetrically discretized, but mathematically correct, numerical schemes, whose small asymmetries can grow to macroscopic size;

2. 2. usage of the usual symmetrically discretized numerical schemes allows asymmetrical order of floating-point evaluation, and these slight asymmetries can grow to macroscopic size;

3. 3. real bugs in the coding.

I would like to point out a fourth possibility. In Reference HindmarshHindmarsh (2004) I carried out a small-amplitude stability analysis of thermo-viscous flow, looking at the specific issue of the spatial scale at which these instabilities were generated, reaching the conclusion that using the SIA did not give a unique or correct estimate of the spatial scale of these instabilities in 3-D. Plane flow does not suffer from this problem, nor is there any evidence that working in 3-D without thermo-viscous coupling introduces instability. This analysis immediately explained why different well-constructed SIA-based models were obtaining entirely different answers. A similar small-amplitude perturbation analysis (Reference HindmarshHindmarsh, 2006a) indicated that including more advanced mechanical models would solve the problem by virtue of the way they smooth the dissipative heating horizontally in a physically appropriate fashion. In Reference HindmarshHindmarsh (2009) I demonstrated how consistent patterns of ice streams could be generated, with different grid sizes, using a numerical model containing horizontally acting stresses. A rule of thumb seems to be that once the grid size is comparable with the ‘membrane coupling length’ (Reference HindmarshHindmarsh, 2006b), around 10–30 km, spurious effects can be expected. Any instability phenomenon that requires resolution better than this is likely to run into serious problems when modelled using the SIA.

My point of view is that modelling thermo-viscous coupling using the SIA is doomed to failure owing to illposedness, and this is the real explanation of the discrepancies between the models. This is not numerical error; the ill-posedness (of the third kind according to Reference HadamardHadamard’s (1902) definition) is a spurious sensitivity to initial conditions at short transverse wavelengths of less than a few tens of kilometres. No numerical scheme, no matter how well designed, can prevent this. Ironically, good numerical schemes will represent this ill-posedness better than poor ones as their increased accuracy at short wavelengths makes them more susceptible to correspondingly sized instabilities. Any noise, such as the sources of numerical error mentioned by Reference CalovCalov and others (2010), can excite instabilities at shorter wavelengths in the SIA-based equations, but not in more complex models. In consequence, I do not believe that different numerical schemes for solving the SIA will remedy the situation.

There is a broader issue of whether, in the presence of instabilities, we should expect model results to be consistent or symmetric in the sense discussed by Reference CalovCalov and others (2010). This comment refers in particular to larger-scale instabilities that might cause a substantial part of an ice-sheet basin to collapse. An instability is, by definition, sensitive to initial conditions. Whether this is ‘ill-posed’ or not is quite a subtle question. In weather prediction this reduces to the difference between ‘weather’ and ‘climate’; the former becomes ill-posed after a few days owing to the loss of predictability, while the latter by construction refers to the well-posed properties of the system. We can only investigate this issue of comparability by experiments with and intercomparisons between models that are not ill-posed in the way SIA models are. I personally do not see a problem in principle with asymmetry developing about the centre line in a symmetrically forced developing surge; this is quite common in other cases of instabilities (e.g. turbulent flow around a cylinder). Here the idea that large-scale basin collapses can occur is the well-posed part of the problem, while I expect the details (e.g. the exact geometry) to be sensitive and in a formal sense ill-posed.

17 September 2011