Reducing uncertainties in projections of Antarctic ice mass loss

. Climate model projections are often aggregated into multi-model averages of all models participating in an intercomparison project, such as the Coupled Model Intercomparison Project (CMIP). The “multi-model” approach provides a sensitivity test to the models’ structural choices and implicitly assumes that multiple models provide addi-tional and more reliable information than a single model, with higher conﬁdence being placed on results that are common to an ensemble. A ﬁrst initiative of the ice sheet modeling community, SeaRISE, provided such multi-model av-erage projections of polar ice sheets’ contribution to sea-level rise. The SeaRISE Antarctic numerical experiments aggregated results from all models devoid of a priori selec-tion, based on the capacity of such models to represent key ice-dynamical processes. Here, using the experimental setup proposed in SeaRISE, we demonstrate that correctly representing grounding line dynamics is essential to infer future Antarctic mass change. We further illustrate the signiﬁcant impact on the ensemble mean and deviation of adding one model with a known bias in its ability of modeling grounding line dynamics. We show that this biased model can hardly be identiﬁed from the ensemble only based on its estimation of volume change, as ad hoc and untrustworthy parametrizations can force any modeled grounding line to retreat. How-ever, tools are available to test parts of the response of marine ice sheet models to perturbations of climatic and/or oceanic origin (MISMIP, MISMIP3d). Based on recent projections of Pine Island Glacier mass loss, we further show that excluding ice sheet models that do not pass the MISMIP benchmarks decreases the mean contribution and standard deviation of the multi-model ensemble projection by an order of magni-tude for that particular drainage basin.


Introduction
The contribution of the Antarctic ice sheet to sea-level rise (SLR) has steadily increased during the last two decades. In the early 1990s, the amount of snow falling over the ice sheet was more or less balanced by the total coastal discharge. Today, the ice sheet loses mass at a rate of ≈ 80 Gt yr −1 , equivalent to ≈ 0.2 mm yr −1 of the global eustatic SLR (Shepherd et al., 2012). Proximal geological evidence shows that the western part of the Antarctic ice sheet may have collapsed during warm periods of the late Pleistocene (Scherer et al., 1998). Such collapses have probably been driven by an unstable retreat of the marine-based regions (i.e., underlying bedrock below sea level), characterized by a retrograde bed slope. The underlying process, marine ice sheet instability (MISI), is supported by theoretical (Weertman, 1974;Schoof, 2007a) and numerical results (Durand et al., 2009;Pattyn et al., 2012). After MISI has been initiated, it could lead to a collapse of the contemporary west Antarctic ice sheet and have the potential to cause sea level to rise by ≈ 3.3 m (Bamber et al., 2009), leading to a drastic impact on human societies (Nicholls and Cazenave, 2010). However, conditions for the initiation of such a collapse, and rate of retreat remains poorly known (Church et al., 2013). The potential for MISI underscores the urgent need for reliable projections of Antarctic mass balance in order to conceive efficient regional and global adaptation strategies.
Current projections for mean sea-level rise in 2100 range from 0.28 to 0.98 m depending on the Representative Concentration Pathways (RCP) scenarios, and the contribution of ice sheets represent about a third of the total projected SLR (Church et al., 2013). However, this likely range excludes the possibility of a collapse of west Antarctica. Since the latest Published by Copernicus Publications on behalf of the European Geosciences Union.
IPCC Assessment Report (AR5), new modeling initiatives tend to show that both Pine Island and Thwaites glaciers may have initiated MISI (Favier et al., 2014;Joughin et al., 2014). Significant progress in the ability of marine ice sheet models to reproduce observed dynamical changes in coastal regions has led to these novel results . Still, ice sheet models have not reached the level of development that other models, simulating other components of the climate system, have reached. Antarctic and Greenland ice sheet model ensembles, in particular, remain in their infancy. Only one attempt has been produced so far, namely the SeaRISE initiative, which has been extensively reported in three pivotal papers by Bindschadler et al. (2013) and Nowicki et al. (2013a, b). Results of all participating models were aggregated into unweighted model averages to produce SLR projections. However, the confidence in related projections remains low because of the unproven ability of many participating models to cope with coastal dynamics (Church et al., 2013).
Parallel to the SeaRISE initiative, specific model intercomparison exercises (MISMIP and MISMIP3d 1 ) have been designed to improve our understanding of grounding line dynamics (i.e., dynamics of the limit between the grounded ice sheet and the downstream floating ice shelf). These initiatives have led to the formulation of requirements regarding physics and numerical approaches to adequately simulate the flow of coastal outlet glaciers in contact with the ocean (Pattyn et al., 2012. In this respect, Favier et al. (2014) proposed a multi-model intercomparison to evaluate the response of Pine Island Glacier (PIG) to changes at the grounding line, based on models meeting these MISMIP and MISMIP3d requirements.
In this paper, we assess the origin of uncertainty in recent ice sheet model projections of Antarctic sea-level contribution for PIG, based on SeaRISE and results due to Favier et al. (2014), Seroussi et al. (2014) and Joughin et al. (2010), guided by MISMIP and MISMIP3d. We further evaluate the potential bias introduced by models that are limited by marine ice sheet physics and reassess SLR projections for Pine Island Glacier based on MISMIP-tested models. We clearly demonstrate the effect of abandoning the "one model, one vote" approach (Knutti, 2010). A brief inventory of the physics implemented in common ice sheet models is presented in Sect. 2. In Sect. 3, we compare the SeaRISE sample with a simplified ice sheet model (SISM) and demonstrate that a proper representation of grounding line dynamics is essential in reducing uncertainties. Finally, a global ensemble analysis for PIG basin is presented. 1 MISMIP: Marine Ice Sheet Model Intercomparison Project.

Stokes equations and approximations
The basic problem in ice sheet modeling is to solve the gravity-driven flow of an incompressible and nonlinear viscous ice mass, further extended with a constitutive equation relating stresses to strain rates, i.e., where τ is the deviatoric stress tensor and D ij are the components of the strain rate tensor. The effective viscosity η is then expressed as where D e is the strain-rate invariant. Models use a temperature-dependent coefficient A, and set n = 3, according to Glen's flow law. A is possibly adjusted through an enhancement factor E that classically accounts for the anisotropic behavior of ice (Ma et al., 2010). The velocity (and pressure field) of an ice body is computed by solving the Stokes problem: where p is the isotropic pressure and g the gravitational acceleration. Apart from the boundary conditions, which are discussed below, this model represents the most complete mathematical description of ice sheet dynamics and is commonly called a full-Stokes model. Owing to the considerable computational effort, approximations to these equations are often used, such as higher-order, shallow-shelf and shallow-ice approximations. These approximations involve dropping terms from the momentum balance equations as well as simplifying the strain rate definitions and boundary conditions. Higher-order Blatter-Pattyn-type models consider the hydrostatic approximation in the vertical direction by neglecting vertical resistive stresses (Blatter, 1995;Pattyn, 2003). A particular case of this type of model is a depth-integrated hybrid model, combining both membrane and vertical shear stress, and of comparable accuracy to the Blatter-Pattyn model (Schoof and Hindmarsh, 2010;Cornford et al., 2013). Vertical shearing terms are included in the calculation of the effective viscosity, but the force balance is simplified. A further approximation, known as the shallow-shelf approximation (SSA), is obtained by neglecting vertical shear (Morland, 1987;MacAyeal, 1989).
However, the earliest and most common approximation in large-scale ice dynamics simulations is the shallow-ice approximation (SIA). This approximation incorporates only vertical shear stress gradients opposing the gravitation drive, The Cryosphere, 9, 2043-2055, 2015 www.the-cryosphere.net/9/2043/2015/ which is valid for an ice mass with a small aspect ratio (i.e., thickness scale much smaller than length scale) in combination with a significant traction at the bedrock. Its main advantage is that all stress and velocity components are locally determined. The approximation is not valid for key areas such as ice divides and grounding lines (Hutter, 1983;Baral et al., 2001), since it excludes membrane stress transfer across the grounding line (Pattyn et al., 2012). The non-validity of SIA is overcome by some models through the use of grounding line flux or grounding line migration parametrizations based on solutions derived from matched asymptotics (Schoof, 2007b(Schoof, , 2011.

Boundary conditions
We will not list all boundary conditions of thermomechanically coupled ice sheet models, but focus on those that are of importance for grounding line migration. These pertain to the initialization of the ice sheet and conditions at the contact of the ice sheet with the ocean boundary.

Initialization
Initialization of ice sheet models to reproduce the current ice sheet state is commonly done through long-term paleo simulations (paleo spin-up). This has the advantage of establishing a reasonable temperature regime within the ice column (Rogozhina et al., 2011). However, reproducing current ice sheet geometry and velocities remains of limited accuracy. To circumvent this pitfall, inverse methods have recently been introduced. Basal drag (or ice viscosity) is inferred by minimizing the misfit between observed and modeled surface velocity (e.g., Morlighem et al., 2010;Gillet-Chaulet et al., 2012) or observed surface elevation (Pollard and DeConto, 2012b). For thermomechanically coupled simulations, a steady-state temperature field is used (Pattyn, 2010), which nevertheless has a limiting impact on shortterm projections .

The marine boundary
It has long been hypothesized that grounding line migration may provoke unstable behavior when the ice sheet rests on a retrograde bed slope below sea level (Weertman, 1974), leading to a potential collapse of marine-based areas (Mercer, 1978). However, despite major developments in numerical ice sheet modeling, the majority of state-of-the-art models in the 1990s and 2000s did not exhibit significant grounding line retreat when ocean forcing (through ice shelf melting) was applied. The response of Antarctic ice sheet models to a warmer climate led to a higher volume due to atmospheric precipitation increase, but not to a mass loss due to inland migration (Huybrechts and de Wolde, 1999;Ritz et al., 2001;Houghton et al., 2001;Solomon et al., 2007). The incapacity of numerical ice sheet models to cope with grounding line migration has been seriously challenged by remote sensing observations at the end of the 1990s showing grounding line retreat and substantial mass loss in the Amundsen Sea sector (e.g., Rignot, 1998). However, ice sheet models were still unable to cope with such rapid changes when IPCC AR4 was released (Solomon et al., 2007). A verification of ice sheet models became feasible due to a boundary layer theory developed by Schoof (2007a), who showed that in the absence of lateral buttressing, grounding line positions are unique and stable on a downward-sloping bedrock. The theory also confirmed the MISI hypothesis in the case of a retrograde bed (unstable grounding line positions), although marine ice sheets are not unconditionally unstable in two horizontal dimensions as buttressing may stabilize grounding line over a retrograde bed . Two model intercomparison exercises have subsequently been organized to verify whether numerical ice sheet models produce results in agreement with the boundary layer theory (MISMIP, MISMIP3d; Pattyn et al., 2012. It further allowed requirements for numerical models to be identified, in order to cope with grounding line migration. Basic conclusions of both intercomparisons are that in order to resolve grounding lines, the inclusion of membrane stresses across the grounding line is required at a sufficiently small grid size (< 500 m), even when a subgrid interpolation of the grounding line (< 5 km) is preferred. An exception is the use of a grounding line parametrization based on the boundary layer theory due to Schoof, which has been successfully implemented in a series of models (Pollard and DeConto, 2012a; Thoma et al., 2014). This works well for coarse spatial resolutions, but the short-term transient response remains questionable when compared to other approaches (Drouet et al., 2013), especially since the theory has been developed for the steady-state case.
In a more recent paper,  further scrutinized the results from the ice2sea MISMIP3d intercomparison to demonstrate that a clear distinction in the response to marine forcing could be related to the complexity of the model physics. The study shows that at least higherorder approximations are necessary to accurately simulate the flow across the grounding line, as the presence of vertical shearing in the force budget softens the effective viscosity at the grounding line, leading to a faster response on short timescales. The result is clearly different for SSA models that are stiffer at the grounding line and seem to overestimate the contribution to SLR. This has also been confirmed by a model intercomparison of Pine Island Glacier (Favier et al., 2014).

Description of SISM (simplified ice sheet model)
To demonstrate the importance of the proper inclusion of a marine boundary in large-scale ice sheet models, we developed a simple (but in terms of marine conditions -wrong) ice sheet model. The simplified ice sheet model (SISM) is a numerical ice sheet model based on the physics inherent www.the-cryosphere.net/9/2043/2015/ The Cryosphere, 9, 2043-2055, 2015 to well-known ice sheet models (Huybrechts, 1990;Fastook and Holmlund, 1994;Saito and Abe-Ouchi, 2004, e.g.,). It is a two-dimensional vertically integrated model, solving the Stokes equations according to SIA. The time-dependent evolution of the ice sheet is based on mass conservation and the model is initialized through 1000 (and 100) years of surface relaxation, starting from the BEDMAP2 data set (Fretwell et al., 2013). Details of the model and the model runs are given in Appendix A. Contrary to the SeaRISE experiments, climate forcing is not applied and the present climate conditions are retained during the whole run (see Appendix A for details). Note that models with a similar degree of physical complexity in the description of ice flow have been included in the SeaRISE multi-model ensemble (Bindschadler et al., 2013).
Grounding line dynamics are not explicitly included in SISM. However, melting at the grounding line is introduced by subtracting the amount of basal melt from the surface mass balance at the last grounded grid point. Ice thickness becomes zero when the ice thickness in that grid point -determined from ice advected from upstream and the local mass balance -becomes zero (or negative). Therefore, grounding line retreat is purely due to melting and not due to any physical process operating at the grounding line. Hence, a marine ice sheet instability (retreat of the grounding line on a retrograde slope in absence of melt perturbation and significant buttressing) is not simulated with this simplified model (Pattyn et al., 2012).
Using SISM, we perform a number of the Antarctic SeaRISE experiments and investigate the impact of including a model with a known bias on the ensemble projection. The SeaRISE initiative led to the first attempt to evaluate multi-ice sheet models ensembles. It is important to note that at the time the experiments were designed, circa 2008, SeaRISE's primary goal was to investigate the sensitivity of ice sheet models to external forcing. Its baseline hypothesis presumed that there was no "best" ice sheet model around and that ensemble modeling would potentially lead to a better understanding of ice sheet models (Bindschadler et al., 2013). Six models participated in the SeaRISE modeling of the Antarctic ice sheet, with a large variety of approximations to the Stokes equations and different treatments with respect to implementing grounding line dynamics (see Table 1). More details on the physics and numerics of SeaRISE models can be found in Bindschadler et al. (2013).
The SeaRISE experiments all start from an initial presentday ice sheet, which is built up using either a paleo spinup or assimilation methods. Perturbations in boundary conditions are then imposed for 500 years and compared to a control run to remove the long-term drift. Climate forcing experiments refer to the ensemble mean of AR4 A1B changes in temperature and precipitation being imposed for 94 years and held constant at the values of year 94 for the remainder of the 500-year runs. An amplification factor of 1, 1.5 and 2, respectively, is applied in order to simulate warmer climate scenarios (experiments C1, C2 and C3, respectively). Subsequently, basal sliding perturbations are implemented through a uniform increase of basal sliding (amplified by a factor 2, 2.5 and 3, respectively, for experiments S1, S2 and S3, respectively) and the sensitivity of Antarctic ice shelves to sub-ice shelf melt was performed by applying a uniform melt rate at the base of floating ice (2, 20 and 200 m yr −1 , experiments M1, M2 and M3, respectively). These sensitivity experiments (or combinations of them) were further used to evaluate the dynamic contribution of Antarctica to sealevel rise for the 21st century under the various RCP scenarios (Bindschadler et al., 2013;Levermann et al., 2014), with estimated median contributions ranging from 0.07 m for the low-emission RCP2.6 scenario and 0.09 m for the strongest RCP8.5 (Levermann et al., 2014). While the reliability of such projections has been questioned (Church et al., 2013), this has not been further evaluated and discussed so far.

PIG ensemble
While SeaRISE focussed on modeling the whole Antarctic ice sheet, a number of studies have simulated the effect of ice shelf melting at the basin scale. Pine Island Glacier (Joughin et al., 2010;Favier et al., 2014;Seroussi et al., 2014) and Thwaites Glacier (Joughin et al., 2014) in particular have shown considerable contemporary grounding line retreat and thinning . Joughin et al. (2010) present a first comprehensive modeling of Pine Island Glacier (PIG) based on a SSA model and using assimilation methods for initialization. Although this particular model did not participate in any MISMIP intercomparison, its physics (SSA) and spatial resolution around the grounding line (down to 140 m) make this model compliant with MISMIP recommendations. Favier et al. (2014) propose a model intercomparison of PIG based on three models of varying complexity, i.e., an SSA model, a higher-order model and a full-Stokes model. These three models also took part in the ice2sea MISMIP3d intercomparison  and produced verified results at the spatial resolution used in the PIG intercomparison. Their approach consists of computing an initial state as close as possible to the current geometry and surface velocities using assimilation methods. Melting and calving perturbations are further applied. They show that the response of PIG is mainly driven by the bedrock topography rather than the type and the amplitude of the perturbation and further conclude that PIG is probably already engaged in a MISI. The study finally estimates the contribution of PIG to SLR over the next 20 years ranging from 3.5 to 10 mm. Finally, Seroussi et al. (2014) use a higher-order model over PIG to simulate its dynamical response to marine forcing over the next 50 years using a higher-order model, with spatial resolutions down to 500 m at the grounding line. This particular model did not perform the MISMIP experiments, but as with the Joughin et al. (2014) model, physics and numerical implementation are conform to MISMIP and MISMIP3d recommendations.
All models presented above will be compared below for the PIG basin. However, we start the analysis with an evaluation of the importance of marine processes on ice sheet response on a pan-Antarctic scale. Figure 1 displays the contribution to SLR after 200 years as a function of the change in grounded area. It clearly highlights the fact that large contributions to SLR always go along with significant changes in the extension of the grounded ice sheet. In other words, having models able to cope with grounding line dynamics is a prerequisite before establishing projections of upper bound dynamic contribution of the Antarctic ice sheet to SLR. The evolution of the grounded area as a response of SeaRISE experiment M3 is presented in Fig. 2. Despite the drastic and unrealistic perturbation (200 m yr −1 melt rate, designed to approximate a sudden collapse of all ice shelves), the response of the participating models varies widely, from a limited grounding line retreat to almost a complete collapse of all the marine sectors within a period of 200 years. Moreover, amongst the models presenting a significant retreat, the impacted regions are different with significant differences in grounding line retreat rates.

Grounding line migration in the SeaRISE ensemble
Large differences in model response are essentially due to two factors: models that correctly implement melting under the ice shelves will fail to produce a significant retreat if the grounding line area is not properly sampled (spatial resolution below 500 m), when using a physical approximation based on SSA, or lacking a parametrization of grounding line dynamics based on the boundary layer theory due to Schoof (2007a). This failure has been clearly illustrated by Vieli and Payne (2005) and Docquier et al. (2011). On the other hand, models that implement melting at the grounding line, i.e., the last grounded grid point, melt grounded ice away, thereby mimicking grounding line retreat. The result is unphysical in both implementation (since melting occurs under the ice shelves) and reaction (spatial resolution and/or physical model are unappropriate). The SISM model illustrates this perfectly, as grounding line retreat in this model is not due to ice-dynamical processes at the grounding line, but due to ice being melted away at the grounded line. The retreat rates produced by this model are within the range produced by SeaRISE, due to the fact that several models within SeaRISE implement grounding line melt in a similar fashion.

Impact of SISM on the SeaRISE ensemble
Since the number of models participating in the Antarctic SeaRISE experiments is rather limited, we may expect that adding a model (e.g., SISM) to the sample will significantly impact on the ensemble mean projections, thereby questioning its relevance. Its effect is illustrated in Fig. 3 and Table 2. Compared to other models, the contribution to SLR with SISM is close to the SeaRISE ensemble mean for sliding experiments and is amongst the largest for melting perturbations, but it is not a striking outlier. As a reminder, SISM is based on simple model physics, isothermal and a surface mass balance set to present-day conditions and not evolving following any RCP scenario. Taking SISM into account in the ensemble unweighted mean leads to two distinct impacts. When considering melt perturbation, adding SISM to the ensemble usually increases both the mean and standard deviation of the ensemble projections. The increase in mean is substantial, up to 20 % for experiments M2 after 100 years. We can anticipate that adding a biased model which would present a limited capacity of grounding line retreat would lead to a decrease of the ensemble mean contribution to SLR The Cryosphere, 9, 2043-2055, 2015 www.the-cryosphere.net/9/2043/2015/  together with an increase in the related standard deviation, as the sample size increases. The particular case of sliding experiments and experiment M3 is instructive: the projected contribution of SISM is fortuitously close to the SeaRISE ensemble mean. Including the SISM in the ensemble mean projections slightly affects the mean but also decreases the standard deviation. Ironically, in the particular situation where a biased model projection is coincidentally close to the ensemble mean, introducing such a model may be wrongly interpreted as improving the confidence in the ensemble projection.
The SeaRISE experiments were rerun with SISM, starting from a different spinup (100 years instead of 1000 years; Fig. 3). Despite significant differences (several hundreds of meters in ice thickness in coastal areas) between both geome-tries, the prognostic runs are hardly affected in terms of SLR contribution over the next 200 years, indicating that an initial state close to the observed ice sheet geometry is a weak constraint when large perturbations are applied.

Ensemble analysis on PIG
In view of the small SeaRISE sample, we extended the sample with recent regional studies (basin-scale), focused on Pine Island Glacier (Joughin et al., 2010;Favier et al., 2014;Seroussi et al., 2014). Since the most significant changes in grounding line position and mass balance are currently observed over PIG , this drainage basin appears to be the most appropriate region to evaluate the impact of model physics/numerics on SLR projecwww.the-cryosphere.net/9/2043/2015/ The Cryosphere, 9, 2043-2055, 2015 tions. Amongst these studies, Favier et al. (2014) argue that PIG is already experiencing MISI and that forthcoming mass change projected by models is relatively similar irrespective of the perturbation amplitude, even for an almost complete collapse of the current ice shelf. Therefore, this implies that comparison with SeaRISE experiments is feasible, despite the difference in melt perturbations between the various studies. However, this comparison remains qualitative as most of the models we compare use different setups that are certainly responsible for a part of the spread between models. Figure 4 presents the evolution of the cumulated contribution of PIG drainage basin to SLR for the period 2000-2050 according to the SeaRISE M3 experiment and according to Joughin et al. (2010), Favier et al. (2014 and Seroussi et al. (2014). As mentioned above, estimations from SeaRISE range from a very limited retreat of the grounded line (e.g., Fig. 2a) and relatively low contribution to SLR (below 5 mm cumulated in 2050) to an extremely high discharge of 3 mm yr −1 and a collapse of the entire drainage basin within a century (e.g., Fig. 2d). As expected, SISM is amongst the models predicting the highest contribution for PIG, but this model result stays within the whole sample range.
A striking feature of Fig. 4 is that all projections due to Joughin et al. (2010), Favier et al. (2014 and Seroussi et al. (2014) occupy a limited range compared to the full range of the SeaRISE sample, with SLR contribution between 2.3 and 18.8 mm by 2040, compared to 2.8 and 146.4 mm, respectively.

Discussion
Most models in the SeaRISE sample have a coarse spatial grid size (> 10 km, see Table 1), which -despite the physical approximations -do not sample grounding line dynamics as stipulated in the MISMIP intercomparison (Pattyn et al., 2012. Only one model uses the parametrization based on the boundary layer theory due to Schoof (2007a), and has been tested against MISMIP, exhibiting a behavior similar to models that do capture grounding line dynamics at high spatial resolution. While such parametrized models are probably less reliable for transient effects, they capture the essence of grounding line migration and stability . The basin-scale simulations for PIG are performed with models that capture grounding line dynamics at the spatial resolution required and with appropriate physics. Most models did participate in the MISMIP intercomparison at the spatial resolution used in the PIG analysis.
Models that capture grounding line dynamics are within the dark gray envelope in Fig. 4 (associated mean of 8.6 mm and standard deviation of 4.9 mm in 2040), compared to the light-gray envelope that represents the SeaRISE sample (mean of 50.1 mm and standard deviation of 67.5 mm in 2040). Not only do the MISMIP-verified models occupy  (Bindschadler et al., 2013). Starting times of experiments computed by Joughin et al. (2010) and Seroussi et al. (2014) correspond to the acquisition year of the surface velocities used for inversion, 1996 and 2008, respectively. As detailed in Favier et al. (2014), the starting time of their experiments corresponds to the last grounding line measurements available, i.e., completed in 2011 (Park et al., 2013). Cumulated contribution was offset to zero in 2009. The light gray (dark gray) envelope encompasses the SeaRISE ensemble (MISMIP-verified ensemble, respectively). a smaller range than the full sample but a distinction between the physical representation of each of the models tends to appear. It is not expected that a model according to full Stokes would show the same results as a SSA model, since the physical model is different. Such a distinction clearly appears for models that capture grounding line migration as was already shown in  based on an ideal ice sheet geometry. For instance, SSA grounding lines move faster in reaction of a perturbation, because they are stiffer at the grounding line (i.e., the viscosity of the grounding line is higher). Higher-order models (and in the limit full Stokes models) produce a slower response, as the viscosity at the grounding line is lower due to the inclusion of vertical shearing in the stress tensor. Such differentiation is not captured whenever the spatial resolution at the grounding line is too coarse and obliterates the characteristics of the physical model. Results computed by Seroussi et al. (2014) with a higher-order model are in agreement with the estimation computed with the higher-order model used in Favier et al. (2014) Favier et al. (2014). Finally, the model based on the parametrized approach has the highest contribution to SLR of the sample of models that capture grounding line migration. Although this comparison remains qualitative because boundary conditions and perturbations differ from one modeling experiment to the other, this result is in line with the results of MISMIPs, which tend to demonstrate that application to a "real" case seems to endorse the conclusions of "idealized" case simulations. However, it would probably deserve specific controlled experiments (i.e., similar forcing, domain, initial geometry . . . ) to better apprehend the intermodel spread.
It is peculiar to note that the models due to Favier et al. (2014) exhibit marine ice instability (and presumably also the case for Joughin et al., 2010 andSeroussi et al., 2014) and their response is to a large extent indifferent to the amplitude of the perturbation applied. Yet, their contribution to SLR on the timescales considered is smaller than the majority of the models that were used in SeaRISE that did not capture any MISI. This poses serious questions as to whether the inherent complexity of an ice sheet model (thermomechanics, sliding, surface mass balance) is decisive in the process of representing ice sheet response to marine forcing. This issue will definitely become important when ice sheet models will be fully coupled to ocean models at a spatial resolution that complies with both systems. Furthermore, coupling could exhibit a series of other feedbacks between the rate of sub-shelf melting and changes in the sub-shelf cavity shape, which are currently unknown.

Conclusions
The SeaRISE initiative has been the first multi-ice-sheet model ensemble projection to evaluate the future contribution of Antarctica to SLR. Results of all participating model results were taken into account, irrespective of the inherent difference in complexity between the models. A similar approach is used in AOGCMs (Atmosphere-Ocean General Circulation Models) community (Knutti, 2010). However, it is probably simpler to evaluate the ability of ice sheet models to perform adequately when compared to AOGCMs, as fewer key processes are at play. In any case, and whatever the component of the modeled climate component, first order processes must be taken into account (Knutti, 2010). As an example, it sounds particularly inappropriate to use an atmospheric model unable to compute a radiation balance to make any projections of mean surface temperature. Similarly, to compute projections of Antarctic contribution to SLR, ice flow models have to be evaluated on their ability of incorporating grounding line dynamics. If this process is not implemented within a model, any of its projections pertaining to coastal regions are unreliable, even on decadal timescales. Furthermore, solely based on the evolution of the modeled ice sheet, it may be hard to discriminate whether the projection is reasonable or not. Indeed, ad-hoc parametrizations can force any model to retreat, but the lack of physics makes any projection of the retreat and contribution to SLR untrustworthy. Owing to the small number of ice sheet models, including biased models has a strong effect on the mean, as well as on the dispersion of the results.
Benchmark experiments to evaluate the ability of models to cope with grounding line dynamics have been recently developed and others will emerge (MISMIP+). Ice sheet models should be evaluated using these benchmarks before being applied on actual cases. Such an approach has been followed by Favier et al. (2014) for PIG. Taking into account only a selection of models with appropriate physics and numerics to compute grounding line dynamics significantly reduces the spread of the projected contribution to SLR, reinforcing our confidence in projecting future glacier evolution. Initiatives to produce new multi-ensemble models will undoubtedly be launched in the near future (e.g., ISMIP6 2 ). Their ability to decrease uncertainties will strongly depend on whether inappropriate models (i.e., unvalidated grounding line dynamics) will be included or not.

Appendix A: SISM description and model setup
SISM (simplified ice sheet model) is an isothermal twodimensional ice sheet model based on the shallow-ice approximation (SIA), using Glen's flow law as a constitutive equation. The vertical mean horizontal ice velocities in the grounded ice sheet are calculated from the local ice geometry (e.g., Huybrechts et al., 1996) u = u b − 2 n + 2 A(ρg) n H n+1 |∇h| n−1 ∇h, where u = (u, v) are the horizontal flow velocities, u b is the velocity due to basal sliding, n = 3 is the flow-law exponent, A = 5 × 10 −17 Pa −n yr −1 is the flow-rate factor, ρ = 910 kg m 3 is the ice density, H (m) is the ice thickness and h (m) is the surface elevation. Basal sliding follows a typical Weertman sliding law, and is defined as where τ d = −ρg H ∇h is the driving stress (Pa), A s = 3 × 10 9 (m yr −1 Pa −m ) a sliding factor and m = 2 (Pollard and DeConto, 2012a). Ice-sheet evolution is based on mass conservation, written as a diffusive equation, i.e., whereȧ is the surface mass balance (m yr −1 ), and D spatially varying diffusion coefficients, defined by D = −u H /∇h. The model is numerically solved on a finite-difference grid with a spatial resolution of 15 km. BEDMAP2 data (Fretwell et al., 2013) are used as basic input for the model. The ice sheet equation (Eq. A3) is solved as a diffusion equation and diffusion coefficients are calculated on a staggered Arakawa B-grid. Since it is an isothermal model, sliding was allowed on the whole domain, but given the very low driving stresses in the interior, its impact on overall ice dynamics remains limited. Initialization of the model is based on a relaxation of the surface elevation, starting from the BEDMAP2 present-day ice sheet geometry, for a period of 1000 and 100 years, respectively. Surface mass balance is obtained from Van de Berg et al. (2006) and Van den Broeke et al. (2006), based on the output of a regional atmospheric climate model for the period 1980 to 2004, and calibrated using observed mass balance rates. During these runs, Eq. (A3) was solved for the grounded ice mask. Starting from this initialization, the model was run forward another 200 years to determine the model drift for the forcings.
Similar to the SeaRISE experiments, the model is forced according to three basal sliding and three basal melt scenarios and for a period of 200 years, during which the basal sliding factor is multiplied by 2, 2.5 and 3, respectively, and basal melting at the grounding line is set to 2, 20, and 200 m yr −1 , respectively. Since the model does not have ice shelves, melting was applied at the last grounded grid point that is in contact with the ocean. Whenever ice thickness in that particular grid point becomes afloat or zero for a bedrock below sea level, the grid point becomes ocean and is not taken up in the grounded mask anymore. This simple, and profoundly wrong, mechanism allows the grounding line to retreat. Grounding line advance is, however, not possible and not anticipated given the significant melting factors applied. Similar mechanisms of grounding line retreat were also applied for some of the models within the SeaRISE intercomparison.