Simulating the Antarctic ice sheet in the late-Pliocene warm period: PLISMIP-ANT, an ice-sheet model intercomparison project

. In the context of future climate change, understanding the nature and behaviour of ice sheets during warm intervals in Earth history is of

recent studies show that the oceanic melting and/or calving of the floating parts of the AIS, the ice shelves, is substantial (Pritchard et al., 2012;Rignot et al., 2013), and that an increase in sub-shelf melting can have a significant impact on grounded ice (e.g. Goldberg et al., 2012). Nonetheless, the sensitivity of the AIS to changes in ocean temperatures remains largely uncertain (Church et al., 2013).
Studying past intervals with a warmer than present-day (PD) climate can be used to gain a better understanding of the sensitivity of the AIS to predicted future climate warming. One such warm interval is the late-Pliocene warm period (also known as the PRISM interval, 3.264 to 3.025 million years before present; Dowsett et al., 2010), which can be considered as a possible analogue for future climate change at the end of this century. This warm period is a well studied interval in Earth's history using both models (e.g. Dolan et al., 2011;Haywood et al., 2013) and data (e.g. Dowsett et al., 2013). During the late-Pliocene, atmospheric CO 2 concentrations are estimated to vary between 300 and 450 ppm (e.g. Bartoli et al., 2011;Martinez-Boti et al., 2015), although an other estimate (Badger et al., 2013) also indicates lower concentrations close to the interglacial values found in ice cores (Lüthi et al., 2008).
Regardless of the rather large uncertainty of the atmospheric CO 2 concentration during this time period, multiple proxy estimates for temperatures show a clear signal of warming over the globe (Dowsett et al., 2010). In recent years the Pliocene Modelling Intercomparison Project (PlioMIP) has provided a framework for studying the Pliocene with climate models (Haywood et al., 2010. PlioMIP includes both atmosphere-only and coupled Atmosphere-Ocean General Circulation Models (AO-GCMs). CO 2 levels for the PlioMIP experiments were set to 405 ppm (Haywood et al., 2010 and further boundary conditions are based on the PRISM3 global reconstruction (Dowsett et al., 2010). For PlioMIP eight AO-GCMs were used to study both changes in sea surface and surface-air temperatures during the late-Pliocene warm period. The models simulated an increase of the global mean surface-air temperature between 1.84 and 3.60 • C compared to the pre-industrial . Although the range of the ensemble is quite large, the models appear to be able to reproduce the broad scale features of the late-Pliocene climate as evidenced in the proxy records Haywood et al., 2013;.
For the late-Pliocene warm period, estimates of the change in sea level are in the order of 10-30 m above PD (Raymo et al., 2011;Rovere et al., 2014;Rohling et al., 2014), although precise values remain uncertain and are possibly influenced by dynamic topography and glacial isostatic adjustment (Rowley et al., 2013). This requires a significant contribution from the Greenland  and the AIS. Large parts of the West Antarctic ice sheet (WAIS) are grounded on bedrock well below sea level (down to ∼ 1500 m). Hence, this ice sheet is more likely to disin-tegrate when ocean and atmosphere temperatures increase. This is shown by both sediment data (Naish et al., 2009) and modelling studies of the AIS for interglacials throughout the Plio-Pleistocene (Pollard and DeConto, 2009;. For the high estimates of sea-level increase, a contribution from the East Antarctic ice sheet (EAIS) is also required. However, thus far transient simulations with ice-sheet models have not been capable of reproducing a substantial retreat of the EAIS (Pollard and DeConto, 2012b;  as has been suggested based on marine sediments (Williams et al., 2010;Cook et al., 2013). A few modelling studies did show that the Wilkes and Aurora basins are sensitive to a particular climate forcing (Hill, 2009;Dolan et al., 2011). However, this result was based on a SIA-only ISM with imposed PRISM3 boundary conditions.
In this study we investigate the nature and behaviour of the AIS during the late-Pliocene warm period in terms of the Pliocene Ice Sheet Modelling Intercomparison Project (PLISMIP; Dolan et al., 2012). For this particular study of Antarctica, from now on referred to as PLISMIP-ANT, we use state of the art ice-sheet models (ISMs) that model both grounded and floating ice. Here, we aim to quantify the sea-level contribution from the AIS calculated with different ISMs forced by output of the HadCM3 AO-GCM Pliocene experiment (the same as Bragg et al., 2012). We have performed control experiments with pre-industrial and PD climate forcing to evaluate the equilibrium response of each model to the PD climate. Additionally we conducted two experiments forced with late-Pliocene HadCM3 climate forcing, one initialised with the PD ice sheet, and one that starts with the much smaller Pliocene ice sheet, as used by HadCM3 and PRISM3. A third Pliocene experiment has been performed with forcing the ISMs with a Pliocene HadAM3 climatology, with PRISM3 boundary conditions, but with a modern Antarctica. The five experiments are summarised in Table 1.

Methods
The basic setup of PLISMIP-ANT follows the experimental design outlined in Dolan et al. (2012). All experiments described here are steady state simulations for 100 000 years (100 kyr). We thus focus on the equilibrated response of the ice sheets to a particular climate forcing. Following Dolan et al. (2012) all models use the same climate forcing and same surface-temperature lapse rate correction of −8 • C km −1 : with T surf the temperature at the surface of the ice sheet and T GCM the temperature field of the climate model in • C, H surf the surface elevation of the ice sheet and H GCM the surface topography of the climate model in metres. Here, we firstly describe the experimental design as implemented specifically The Cryosphere, 9, 881-903, 2015 www.the-cryosphere.net/9/881/2015/  (Lythe et al., 2001) or Bedmap2 (Fretwell et al., 2013) and PRISM3 (Dowsett et al., 2010). PI: pre-industrial. for PLISMIP-ANT. Secondly, the different climatology of each experiment is described and compared. Lastly, we specify the setup of the ISMs.

Experimental design
To force the ISMs over Antarctica we use the monthly climatology obtained from simulations using the Hadley Centre Coupled Atmosphere Ocean Model version 3 (HadCM3; Pope et al., 2000;Gordon et al., 2000), which are set-up following the PlioMIP experimental design  and are comparable to those presented in Bragg et al. (2012). The participants use their models to run the ISM over Antarctica including ice shelves. All ISMs are initialised using a suite of boundary conditions. A complete description of PLISMIP within PlioMIP is given by Dolan et al. (2012); however, for PLISMIP-ANT some modifications have been made to the experimental design that we describe here. As shown in Table 1, for PLISMIP-ANT five experiments are performed. In addition to the atmospheric forcing of precipitation and temperature, we also use yearly averaged ocean temperatures from the ocean model of HadCM3 as input for the ISMs. The main setup of the experiments requires ISMs that are fully capable of simulating both grounded and floating ice, using the approximate stress balance equations for ice sheets; the shallow ice approximation (SIA) and for ice shelves; the shallow shelf approximation (SSA). Additionally, as a sensitivity experiment we also include three SIAonly models that are only capable of simulating grounded ice.

Model boundary conditions
The model spread between the eight climate models used in PlioMIP is quite substantial (see Fig. 3 in Haywood et al., 2013). However, the different models seem to largely agree that a significant polar amplification has occurred over the Antarctic region. For the current setup of the project, we use one AO-GCM climate forcing from PlioMIP (HadCM3). Our choice of HadCM3 as the initial climate model for use in PLISMIP-ANT is based on the long history of usage of HadCM3 within Pliocene climate modelling (e.g. Haywood and Valdes, 2004) and because it is a model with an average climate sensitivity, which simulates average Pliocene temperature changes in relation to predictions from the rest of the PlioMIP ensemble. For the first control experiment, we use a pre-industrial simulation of HadCM3 that includes the PD topography and ice sheets and a pCO 2 of 280 ppm. The yearly mean climatology is shown in Fig. 1a- We perform a second control simulation using PD reanalysis of ERA-40 for the surface climatology (Uppala et al., 2005) (Fig. 1e,f) and the World Ocean Database 2009 data set for ocean temperatures (Boyer et al., 2009) as illustrated in Fig. 1h. The climatology of ERA-40 is a bit warmer, averaged over the continent 4 • C and a bit wetter, around 0.2-0.5 m yr −1 more precipitation in coastal areas compared to the pre-industrial simulation of HadCM3 (Fig. S1). However, the largest differences occur over the interior of East Antarctica, where precipitation is up to a factor 5 lower. This has quite a significant influence on the reconstructed ice volume as will be shown later on. Although ERA-40 is also subject to uncertainty over Antarctica (Bromwich and Fogt, 2004;Bracegirdle and Marshall, 2012), we use here ERA-40 and the WOD-09 data sets as a secondary control test to simulate the present-day ice sheet, and to show the response of the ISMs to a different climatology for the late Holocene.
Two of the Pliocene simulations are forced with the Pliocene run of HadCM3, which uses the PRISM3 boundary conditions and a pCO 2 of 405 ppm, illustrated in Fig. 1i-l. Here ocean temperatures are depicted at the bottom of the PD ice-shelves of Bedmap1, which are horizontally extrapolated from the nearest ocean grid points since HadCM3 uses a modern land-sea mask, i.e. the alternate experimental design as given by Haywood et al. (2011). Outside the ice shelves, sea surface temperatures are shown. During model www.the-cryosphere.net/9/881/2015/ The Cryosphere, 9, 881-903, 2015 simulations, ocean temperatures are vertically adjusted according with the depth of the bottom of the shelves. Mainly due to the smaller AIS in PRISM3 the surface-air temperatures over Antarctica are warmer by about 7 • C on average compared to the pre-industrial climate. Similarly, the absence of ice in the Wilkes and Aurora basin results in an increase in annual total precipitation of about 0.4-0.6 m yr −1 over this particular region. Large temperature differences are also found in the ocean where sub-surface temperatures show a widespread increase of ∼ 2.6 • C on average below the ice shelves (Fig. S1c).
To investigate the significance of the imposed PRISM3 boundary conditions and starting conditions, we perform a third Pliocene experiment, Pliocene PD-Ant (Table 1). The ISMs are initialised with the equilibrium results of the Control HadCM3 run and forced with a Pliocene HadAM3 climate model simulation that has all the PRISM3 boundary conditions except using a modern Antarctic ice sheet (Phase 2 experiment in Dolan et al., 2012). The HadAM3 climatology is shown in Fig. S2. We use here the HadAM3 simulation since no such simulation is available for HadCM3, henceforth for ocean temperatures the HadCM3 Pliocene values are used. In comparison with the pre-industrial HadCM3 climate ( Fig. S2c-d) there is a significant increase in both temperatures and precipitation around the margins of the ice sheet.
For all ISMs we have provided monthly climatology of surface-air temperature and precipitation and yearly mean ocean temperatures at 19 depth levels for HadCM3 and 30 levels from the WOD-09 data set, ranging from the surface to ∼ 4.5 km depth. As a lower boundary condition for the 3-D ice-sheet temperature field, the preferred boundary condition is taken to be the heat-flux field from Shapiro and Ritzwoller (2004). For the initial ice-sheet thickness and bedrock topography we have used the Bedmap (Bedmap1) data set (Lythe et al., 2001) for the PD configuration and the PRISM3 ice sheet (Dowsett et al., 2010) for the Pliocene. We have also performed the same experiments using the recently published updated bedrock data set of Bedmap2 (Fretwell et al., 2013).
All ISMs are run on a 40 km by 40 km grid, 167 × 167 grid points. Climate fields are projected on this grid using a stereographic projection with OBLIMAP v2.0 (basic theory described in Reerink et al., 2010). For the projection, a central longitude was used of 0 • E, the central latitude was set to the south pole (i.e. hence in this case a polar stereographic projection) and the angle that defines the standard parallel was set to 24.7 • (for details see Reerink et al., 2010). The projection we use here requires a correction for the area of the grid points of the ISMs, for which we follow the methods described in Snyder (1987) with a latitude of true scale of 65.3 • . All volumes of the ice sheet and the contribution to sea level are calculated using the corrected area of each grid point.

Ice-sheet models
For simulating the Antarctic ice sheet over its complete domain in PLISMIP-ANT, we use ISMs that solve both ice flow for grounded and for floating ice. The models in this study include approximate equations of the Stokes equations of flow. The approximations are primarily based on the shallowness of a large ice body, with spatial scales that are much larger than the thickness of the ice. For grounded ice, the SIA (Appendix A1; Hutter, 1983) is used. The SIA only considers horizontal shear stresses and assumes the force of gravity to be the main driver of ice flow. On the other hand, for the ice streams and ice shelves, horizontal stretching, or a membrane-type flow is dominant, which is described with the SSA (Appendix A2; Morland, 1987). Both approximations are based on an isotropic, i.e. Glen's flow law, whereas in reality ice is a highly anisotropic material, which behaves differently for vertical shear stresses (SIA) and horizontal stretching (SSA) (Ma et al., 2010). Therefore, different enhancement factors for the flow parameter are used for SIA and SSA flow (see Appendix A). Recent developments in ISMs also include higher-order physics, or the full-Stokes solution of 3-D ice flow (see for example Pattyn et al., 2013). However, for paleoclimate applications that largely investigate the long term, on the order of 10-100 kyr, response of ice sheet, shallow models are still predominantly used (e.g. Huybrechts, 2002;Pollard, 2010;Pollard and DeConto, 2012b;Golledge et al., 2012;. All participants were asked to set up the ISM in a standard mode. In other words, use the models at their regular setup with their own parameter settings for the thermodynamics, mass balance and ice flow as would be used for regular simulations. The reasoning behind this is that we get an estimate in the differences in ice volume between different modelling groups that use their normal setups of the models, as they are used for other applications as well and that they simulate a more realistic modern ice sheet. By including the fixed lapse rate correction, Eq. (1), all ISMs are initially forced with the same surface temperatures and precipitation fields from the climate models.
All six ISMs that are used calculate ice-velocities with the SIA and SSA, see Table 3. Since it would be too exhaustive to describe here all aspects of the different models, we will provide a short description of each model and its specific methodology of calculating ice velocities, the surface mass balance and how the sub-shelf melting is included using the ocean temperatures from the climate forcing. The latter is described below, since this is generally a new aspect in most models. For a more detailed description of each ISM, the reader is referred to their respective references as included at the bottom of Table 3. All models incorporate a bedrock model, which is adjusted to changes in ice loading. For all models, the basic elastic lithosphere-relaxing asthenosphere (ELRA) model has been used (Le Meur and Huybrechts, 1996).  Table 1). Top panels shows results from a pre-industrial run of HadCM3. Middle panels is ERA-40 (1971ERA-40 ( -2000 (Uppala et al., 2005) and ocean temperatures from the WOD-09 data set (Boyer et al., 2009). Bottom panels illustrate the Pliocene HadCM3 run with full PRISM3 boundary conditions . From left to right, surface-air temperature in • C, Precipitation in m yr −1 water equivalent, sea surface temperatures and temperatures at the bottom of the PD ice shelves in • C and surface topography in the climate model in m. The black line in all panels represents the Bedmap1 outline of the grounding line.
A new aspect for most of the ISMs used in PLISMIP-ANT is the sub-shelf melting, or basal mass balance, which includes a parameterisation using ocean temperatures as climate forcing. For recent and future mass loss of the AIS, oceanic sub-shelf melting has been found to be significant (Pritchard et al., 2012;Rignot et al., 2013) and as such it is an important component to be included in the total mass budget of the ice sheet, especially for the much warmer ocean temperatures of the late Pliocene (see Fig. 1c and k). Most models use a parameterisation as described by Holland and Jenkins (1999) and Beckmann and Goosse (2003): with the different parameters as described in Table 2. T oc is the temperature of the ocean underneath the ice shelf, as vertically interpolated from the 3-D ocean temperature fields from the climate forcing. T f is the freezing temperature as given by Beckmann and Goosse (2003): with S 0 a mean value for the salinity of the ocean of 35 psu and z b the bottom of the ice shelf below sea level. The sub- Table 2. Physical parameters for the sub-shelf melt parameterisation.

Constant and description Value
shelf melt parameter F melt varies between ice-sheet models and is given in Table 3. Since the HadCM3 climate model does not resolve all points underneath the ice shelves, the ocean temperatures are extrapolated using a distance weighting scheme (similar to Maris et al., 2014). The SMB is largely calculated using the same method in all models. Precipitation is taken from the climate forcing and from this snow accumulation is determined depending on the surface temperatures. All models except ANICE determine surface melting with a positive degree-day (PDD) scheme (Reeh, 1991), with a PDD factor for ice melt of 8 and Table 3. Description of the ice-sheet models used for PLISMIP-ANT. All models apply the climatological forcing of temperature and precipitation with absolute values. Models are run on a 40 km by 40 km grid. For the bottom boundary condition of the ice temperature the heat flux field of Shapiro and Ritzwoller (2004) was used. The surface temperature is corrected with a surface lapse-rate of −8 • C km −1 . SMB: surface mass balance, fd: finite difference, SIA: shallow ice approximation, SSA: shallow shelf approximation, PDD: positive degree day, ITM: insolation-temperature melt, BG03: Beckmann and Goosse (2003).  to ice-shelf bottom.

AISM-VUB
The Antarctic ice sheet model (AISM) from the Vrije Universiteit Brussel (VUB) has been initially developed by Huybrechts (1990Huybrechts ( , 2002 and was further improved by Fürst (2013). For the present coarse resolution experiments, SIA and SSA are calculated separately for grounded and floating ice and coupled across a one grid cell wide transition zone. It is assumed that depth-averaged longitudinal deviatoric stresses only contribute to the effective stress at the grounding line (Huybrechts, 1990;Huybrechts and de Wolde, 1999). The ice shelf front coincides with the 150 m thickness contour. Sliding is calculated using a Weertman sliding relation inversely proportional to the height above buoyancy wherever the ice is at the pressure melting point. Surface melting is calculated with the PDD scheme, including meltwater retention by refreezing and capillary forces in the snowpack, driven by the surface temperature field of the climate forcing. Parameter settings are given in Table 3. Sub-shelf melting is parameterised as a function of local ocean-water temperature above the freezing point using Eq. (1). A distinction is made between protected ice shelves (Ross and Ronne-Filchner) with a melt factor of F melt = 5.2 × 10 −3 m s −1 and all other ice shelves with a melt factor of F melt = 21.8 × 10 −3 m s −1 . The parameters are chosen to reproduce observed average melt rates (Depoorter et al., 2013) under the Ross, Ronne-Filchner and Amery ice shelves for WOD-09 temperature observations and Bedmap2 shelf geometry. For the Control HadCM3 run and initial Bedmap1 geometry, average melt rates are a factor 2.0-2.5 too high for these three ice shelves. The two control and the Pliocene Ice-PRISM3 simulations are using an initial spin-up with fixed geometry for 10 kyr and consecutively for 40 kyr with fixed grounding line before the unconstrained 100 kyr simulations. The Pliocene Ice-PD simulation is integrated forward for 100 kyr restarting from the PD steady state of the Control HadCM3 simulation.

ANICE
The ANICE model is part of the IMAU-ICE package (Institute for Marine and Atmospheric research Utrecht), the icesheet model of Utrecht University. The package contains a range of ISM of different complexities, from shallow 1-D models to a full-stokes application. ANICE calculates both the SIA and SSA velocities for sheet and shelf ice, and add these together, no additional grounding-line parameterisations are included. Basal sliding is included as a Mohr-Coulomb plastic law, with basal stresses included in the SSA equations. The basal stress is calculated as a function of a till stress that spatially varies as a function of bedrock elevation de Boer et al., 2013). Surface melt-ing is calculated using an insolation temperature melt (ITM) model, using monthly values of the PD insolation at the top of the atmosphere and surface-air temperature (de Boer et al., 2013). The monthly precipitation field is adjusted with the change in surface temperature, the latter is adjusted according to Eq. (1). Furthermore, refreezing of rain and melt water is calculated using a potential retention fraction. Sub-shelf melting is calculated as described above, and only applied on grid points that are completely floating with the melt factor F melt = 2 × 10 −3 m s −1 and is combined with melt rates for exposed ice shelves and the deep ocean (Pollard and De-Conto, 2009;de Boer et al., 2013). No additional calving law is applied.

PISM
The parallel ice sheet model (PISM) used for this project is the most recent version v0.6 Feldmann et al., 2014). Velocities from the SIA and SSA are combined to yield total velocity . PISM v0.6 includes a sub-grid scheme described in Feldmann et al. (2014) that improves grounding line migration.
There are two components of this scheme: firstly, a subgrid interpolation system allows for the "correct" position of the grounding line to be estimated according to the flotation criterion. The ratio of grounded to floating area represented by the 2-D interpolated position is then used to modify the driving stress in that cell. Interpolated basal melt is neglected here. Secondly, surface gradients at the grounding line are calculated in one direction only (up-glacier in the last grounded cell, and down-glacier in the first floating cell). Basal sliding is included as a Mohr-Coulomb plastic law, with basal stresses included in the SSA equations ). An elevation-dependent prescription of the till friction angle is used (see Martin et al., 2011), ranging from 6 • for all areas of bedrock below 100 m elevation and linearly increasing to 15 • for all areas where the bed is above 1500 m elevation. Additionally, the subglacial till layer is also weakened by saturation of meltwater generated at the ice-sheet bed by geothermal, frictional and strain heating (Golledge et al., 2014).Variability in modelled ice volume in the PISM simulations arises from a thermodynamic feedback in which increased basal sliding (leading to volume loss) is the threshold response to a gradual saturation from meltwater saturation of the basal substrate layer. Surface melting is calculated with the PDD scheme. The subshelf melting rates are calculated with a modified form of the quadratic parameterisation of Holland et al. (2008): Here, T oc is used at a fixed depth of 600 m, as this was considered most representative of the water depth affecting most of the PD ice shelves. Additionally, two calving criteria are used: firstly, the eigen calving approach of Levermann et al.
(2012) that predicts calving losses according to horizontal spreading rates, and secondly a thickness limitation is imposed, such that shelves thinner than 250 m are automatically calved. The latter is a tuned value found through experimentation to yield ice shelf extents of reasonable fit to observed geometries.

PSU-ISM
The Penn State University (PSU) ISM has been widely used for paleoclimate applications (e.g. Pollard andDe-Conto, 2009, 2012a, b). The most recent version includes a grounding-line flux boundary condition as introduced by Schoof (2007), whereas a heuristic scheme is used to determine the transition from sheet to shelf ice flow (Pollard and DeConto, 2012b). Sliding is included as the standard Weertman sliding, but the basal sliding coefficients were tuned to minimise modern-day ice surface elevation errors (Pollard and DeConto, 2012b). The tuned coefficients are adopted from Pollard and DeConto (2012b), the tuning is not repeated in this study. Surface melting is included using a basic PDD scheme. The sub-shelf melt rates use the same Eq.

RIMBAY
RIMBAY is based on the 3-D ISM by Pattyn (2003) and a full description is given in Thoma et al. (2014). Here we use the shallow-approximation version of RIMBAY that combines SIA and SSA velocities in a similar way as PISM and ANICE. in RIMBAY the SSA and SIA velocities are added together with a smoothing gradient over the grounding line (Thoma et al., 2014), which mixes SIA and SSA velocities over 2 grid boxes, i.e. a distance of 80 km, to smooth the transition between SIA and SSA regions. Sliding is included with a basic Weertman sliding law, for which the sliding velocity is a function of the driving stress (see Thoma et al., 2014). Surface melting is calculated with a PDD scheme. Sub-shelf melting is calculated as described above with the melt factor F melt = 11 × 10 −3 m s −1 .

SICOPOLIS
Here we use SICOPOLIS (SImulation COde for POLythermal Ice Sheets) version 3.2-dev revision 498. The model calculates the SIA and SSA separately for sheet and shelf flow, respectively. The enhancement factor for ice flow on land are separate for glacial and interglacial ice. E SIA = 5 for glacial ice (older than 11 kyr for the Control simulations) and E SIA = 1 for interglacial ice, consistent with measurements from ice cores (NEEM community members, 2013). No additional grounding-line parameterisation or combinations are used. Sliding will initiate when ice at the base reaches its local melting point by applying a Weertman-type law in the form used in Sato and Greve (2012). This includes sub-melt sliding before reaching the melting point, when the ice is 1 • C below the melting point, it starts sliding with a fraction of the full sliding velocity, which is reached when the temperature is at the melting point. Surface melting is calculated with the PDD scheme, supplemented by the semi-analytical solution for the PDD integral by Calov and Greve (2005). Further, the model implements a retention model that takes into account the contribution from rainfall and surface melt to the formation of superimposed ice, for which a saturation factor of 0.6 is chosen (Reeh, 1991). The sub-shelf melting parametrisation is as described above, with different melt factors, F melt = 5 × 10 −3 m s −1 for protected, F melt = 5 × 10 −2 m s −1 for exposed and F melt = 5 × 10 −1 m s −1 for open ocean shelves. Melting at the grounding line points is included using the regression of Rignot and Jacobs (2002). This regression is only used for the control experiments. A simple ice thickness threshold method is used for calving, with a value of 50 m to enable calving, only applied at exposed calving fronts, the latter determined as in Pollard and DeConto (2009).

Results
All experiments are 100 kyr steady state runs; i.e. a constant climate forces the ISMs, for which only surface temperatures are adjusted with a constant lapse rate, Eq. (1), and ocean temperatures are adjusted according with the depth of the bottom of the shelves. Figure 2 shows the full 100 kyr simulated ice volume of all models for the first four experiments of PLISMIP-ANT, performed with all six ISMs. The model behaviour varies considerably due to differences in specifying initial conditions between the models, i.e. initial ice temperatures and differences in calculating velocities and the surface mass balance. In general, the models do come into an equilibrium state quite rapidly.

Modern control simulations of Antarctica
For PLISMIP-ANT two control simulations have been performed. The first simulation is the basic test for a comparison with the Pliocene HadCM3 forcing and uses a pre-industrial simulation of HadCM3 (Fig. 1a-d). Differences in the timeevolution of the models can be mainly ascribed to the initial ice-temperature distribution and consequently velocities of the ice (Fig. 2a). Most ISMs display the strongest increase in ice volume at the beginning of the 100 kyr simulation and the variability thereafter is quite limited. However, some models, such as the PSU-ISM (green), show a smaller initial increase in ice volume and a steady increase of volume. The small initial increase could be due to the initialisation procedure and the different basal-sliding parameterisation (Sect. 2.3.4). The different initial ice volume of AISM is caused by the initialisation procedure performed before the equilibrated 100 kyr run (see Sect. 2.3.1). Variability in ice volume in the PISM simulations arises from a thermodynamic feedback in which increased basal sliding (leading to volume loss) is the threshold response to a gradual saturation of the basal substrate layer by meltwater. Under a constant climate forcing, these glaciological feedbacks give rise to an ice sheet that is in a state of dynamic equilibrium (Van Pelt and Oerlemans, 2012; Golledge et al., 2014).

The
In general, the final grounded ice volume between models differs quite substantially, 25.9 to 29.6 × 10 6 km 3 for the Control HadCM3 simulation, whereas the Bedmap1 PD grounded ice volume is 26.6 × 10 6 km 3 (Fig. 4). Nonetheless, the topography and the extent of the ice shelves are similar compared to the PD initial ice sheet for almost all models (Fig. 3). Only the PSU-ISM does exhibit a retreated Ross ice shelf, largely induced by substantial sub-shelf melting. Although the setup of the individual ISMs is different, the difference of the final topography of the ISMs compared to the observed Bedmap1 surface elevation shows a rather similar pattern, i.e. a lower topography in the interior of the ice sheet and thicker ice around the edges, especially in the Lambert ice shelf, the Antarctic Peninsula and Coats land (see Fig. S3). The lower ice thickness in the centre and the thicker edges is a common feature in shallow ice-sheet models and can be ascribed to the course spatial resolution that does not capture the detailed topography and unknown basal conditions, especially for Bedmap1 (e.g. Martin et al., 2011;Maris et al., 2014;Briggs et al., 2014).
Differences between the ISMs are difficult to assess, but are a combination of changes in either the SMB, ice flux across the grounding line and basal melting, largely combined in a non-linear way. The larger volume of the PSU-ISM and AISM could be due to a slightly higher SMB, since the ice flux across the grounding line is quite large. The smaller volume and extent of SICOPOLIS (Fig. 4) is due to more basal melting at the grounding line relative to the other models. The Control Obs simulation with ERA-40 and WOD-09 climate forcing in general shows a smaller ice sheet compared to the pre-industrial control with HadCM3, with a range of grounded ice volume of 24.11 to 26.86 × 10 6 km 3 as illustrated in Fig. 4. The lower ice thickness in the interior of the EAIS is mainly due to lower precipitation over this area, which is known to be underestimated in ERA-40 and models of present-day climatology ( Van de Berg et al., 2005). Most ISMs do reconstruct an ice sheet that remains comparable to the PD ice volume and extent (Fig. S4). The extent of the ice shelves is simulated less well, in particular by PISM and SICOPOLIS, due to lower ice velocities across the grounding line and a lower SMB over the ice-shelves areas.
When comparing the two control experiments (Fig. 4), AISM simulates ice sheets that are both larger than Bedmap1, whereas SICOPOLIS simulates ice sheets smaller than Bedmap1, the latter with a smaller extent of grounded ice, mainly due to locally high rates of grounding line melting. RIMBAY and ANICE simulate ice volume closest to Bedmap1, but ANICE shows a smaller ice extent. The largest difference between the two control simulations is shown by PISM and the PSU-ISM (Fig. 4), which could be attributed to the difference in SMB between the two simulations over grounded ice that is larger for PISM and the PSU-ISM relative to the other ISMs. In general, for both control experiments and for all six ISMs a larger grounded ice area corresponds to a larger grounded ice volume (Fig. 4).

Antarctica in the late Pliocene
As shown in Table 1, we have performed two Pliocene experiments with the same HadCM3 climate forcing. Pliocene Ice-PD simulation uses the PD ice sheet as an initial state for the ISMs (Fig. 5a), whereas the Pliocene Ice-PRISM3 simulation is initialised with the much smaller PRISM3 ice sheet topography (Fig. 6a). For both simulations the response over the 100 kyr simulations is very different for the ISMs (Fig. 2c and d). For the Pliocene Ice-PD experiment the AISM, PSU-ISM and RIMBAY show an increase in ice volume, whereas ANICE, PISM and SICOPOLIS show an initial drop and then a recovery to a lower volume than initially (Fig. 2c). The three models with a smaller ice sheet behave in a similar way in the Pliocene Ice-PRISM3 simulation, as shown in  The final ice-sheet topographies are shown in Figs. 5 and 6. The warmer ocean temperatures in the late-Pliocene climate forcing (see Fig. 1k compared to Fig. 1c) result in complete disintegration of the ice shelves from the PD initial ice sheet for all models except RIMBAY (Fig. 5f), which might be partly due to the grounding-line physics included (Sect. 2.3.5). The differences in grounded ice volume between the other models are largely due to the differences in SMB.
For all six ISMs the ice sheet has a larger grounded volume in the Pliocene Ice-PD simulations compared to the Pliocene Ice-PRISM3 simulations (Fig. 7a). Moreover, when considering East and West Antarctica independently ( Fig. S5), there are some interesting features within the PLISMIP-ANT ensemble. Volume predictions for East Antarctica vary from 22.04 to 25.45 × 10 6 km 3 using the PD as an initial condition and from 21.01 to 24.08 × 10 6 km 3 using the PRISM3 ice sheet to initialise the ISMs. None of the models sustain the extent of retreat given as initial condition in PRISM3 (Figs. 5 and 6), due to a positive SMB over the Wilkes and Aurora basins. Furthermore, the largest difference between the two Pliocene simulations for the WAIS (see Fig. S5b) is simulated by RIMBAY, 1.79 to 6.42 × 10 6 km 3 and by the PSU-ISM 1.95 to 3.98 × 10 6 km 3 .
To test the impact of the imposed PRISM3 boundary conditions, the Pliocene PD-Ant experiment forced with Pliocene  HadAM3 climatology, with PRISM3 boundary conditions except a modern Antarctica, exhibits a good test if the models can also simulate a retreat from the present-day grounding line. Although we impose warm ocean and surface-air temperatures on the ice sheet, the ISMs do not exhibit a retreat of the ice sheet (Fig. 8). The PD extent of the ice sheet is maintained due to the imposed increase in precipitation, and ice fluxes across the grounding line that are larger due to thicker ice at the margins (Fig. S3). The grounded ice volume and extent of the Pliocene PD-Ant simulation is consistently larger for each ISM, except RIMBAY, relative to the Pliocene Ice-PD simulations ( Fig. 7a and c), the latter is forced with the full Pliocene HadCM3 climatology. Figure 9 provides an overview of the results from the four main experiments in terms of ice-sheet thickness. On a gridbox by grid-box basis, the ensemble of results has been sorted into order of thickness (thinnest to thickest members) and we have plotted the median (i.e. the mean of the 3rd and 4th member; Fig. 9, left panels) and the range (i.e. the difference between the largest and smallest ice thickness, divided

Intercomparison
The by two; Fig. 9, middle panels). Finally, we have also plotted ice-sheet presence that shows how many of the six ISMs predict ice of any thickness, both floating and grounded ice, in that particular grid box. The ice presence maps (Fig. 9, right panels) show that all models reconstruct an EAIS of nearmodern extent for the late-Pliocene, and that no ISM simulates the retreat in the Wilkes and Aurora basin as prescribed by the PRISM3 boundary conditions. The range of ice thickness in model predictions illustrates the degree of model agreement among the PLISMIP-ANT ensemble. The differences among the models are large, in particular for the Pliocene Ice-PD simulation over West Antarctica (Fig. 9h). In addition (Fig. 9i) shows where some models suggest a large-scale ice cover across West Antarctica and others only small ice caps. For the Pliocene Ice-PRISM3 simulation, initialised with the PRISM3 ice sheet, the median shows a smaller WAIS, whereas the EAIS is similar to that of the Pliocene Ice-PD simulation (Fig. 9j) warm ocean temperatures (Fig. 1k) and a negative SMB. The largest range in thickness for the Pliocene Ice-PRISM3 experiment is exhibited over the Wilkes and Aurora basin in East Antarctica (Fig. 9k). All ISMs simulate an advance over this area of the ice sheet relative to the initial PRISM3 ice-sheet topography (Fig. 9l), ice-sheet advance is caused by a positive SMB, whereas variability in extent and ice thickness is due to differences in ice velocities. In Fig. 10 we show cross sections through the ice sheets resulting from the six ISMs for three locations. The cross sections are shown in Fig. 9a and represent the regions with the largest spread among the models. Figure 10a and b show the cross sections through the Wilkes basin (A in Fig. 9a), for the Control HadCM3 and Pliocene Ice-PD simulations respectively. For the Control HadCM3 simulation (In Fig. 10a) the modelled topographies are largely similar, whereas the spread between the ISMs for the Pliocene Ice-PD simulation (In Fig. 10b) is notably larger. For the Lambert ice stream (B in Fig. 9a) the spread among the ISM simulations is large for the four experiments. All models generally show thicker ice for the Control HadCM3 simulation relative to PD (see also Fig. S3), whereas for the Pliocene Ice-PD the ice thicknesses vary but cluster around the initial PRISM3 surface elevation (In Fig. 10d). Similar for the Lambert ice stream, all ISMs simulate a higher topography over the cross-sectional area of the WAIS (Fig. 10e). As can be seen in Fig. 9h, the spread between the models is particularly large over West Antarctica for the Pliocene Ice-PD simulation, with some models simulat-ing a large ice sheet and others predicting ice-free conditions (Fig. 10f), due to a balance between grounding-line flux, subshelf and surface melting.

Sensitivity to initial bedrock topography
Recently, a new data set of bedrock topography, surface elevation and ice thickness for Antarctica (Bedmap2; Fretwell et al., 2013) has been published. To analyse the sensitivity of the modelled ice-sheet geometry to the new bedrock topography inferred from observations we have repeated the same experiments with the Bedmap2 data set, using exactly the same setup as was used for the original experiments. The Bedmap2 data set provides a significant improvement relative to the Bedmap1 data set used here, i.e. higher resolution, improved data coverage and precision (Fretwell et al., 2013). Moreover, Bedmap2 contains fewer inconsistencies between surface elevation, ice thickness and bedrock topography, which was a limitation in the Bedmap1 data set (Fretwell et al., 2013).
To repeat the experiments, a new initial Pliocene ice-sheet topography had to be generated for the Pliocene Ice-PRISM3 simulation. Here we have placed the PRISM3 ice-sheet configuration on the Bedmap2 bedrock topography. To account for the uplift of the bed due to the retreat of the ice sheet, relative to the Bedmap2 ice thickness, the ELRA bedrock model has been used to run the bedrock topography to isostatic equilibrium. The final bedrock topography and ice-sheet surface are then used as initial fields for the Pliocene Ice-PRISM3 exper-The Cryosphere, 9, 881-903, 2015 www.the-cryosphere.net/9/881/2015/ iment as shown in Fig. 11a. In general, differences with the original PRISM3 ice sheet are not large. However, bedrock elevation is significantly lower in the Wilkes and Aurora basin (see Fig. 11h). Naturally, some uncertainties are involved in this procedure such as the chosen bedrock model and its parameters and the accompanying uncertainties in the Bedmap2 data set (see Figs. 11 and 12 in Fretwell et al., 2013). However, we believe this is a reasonable first sensitivity test to identify how the ISMs will respond to a different initial bedrock topography. Figure 12 represents a comparison between the two simulations of the grounded ice volume and area. The Bedmap2 simulations for Control HadCM3 are also comparable to the PD ice-sheet extent and ice volume. Final ice volume for the Pliocene Ice-PD and Pliocene PD-Ant experiments are comparable for both bedrock topographies, with 24.7 to 31.5 × 10 6 km 3 and 26.1 to 31.3 × 10 6 km 3 for Bedmap1, and 23.1 to 30.3 × 10 6 km 3 and 26.7 to 29.4 × 10 6 km 3 for Bedmap2, respectively.
As is shown in Fig. 11, the final simulated surface topography for the Pliocene Ice-PRISM3 experiment shows a different result especially for the Wilkes and Aurora basin, where observations have improved considerably compared to Bedmap1. As shown in Fig. S6a, most models calculate an even lower volume than the initial PRISM3 ice sheet, also due to a reduced size of the central part of the ice sheet, whereas the area covered by ice is still larger (see Fig. S6c) and the Pliocene Ice-PD simulations results in a modern extent of the EAIS (Fig. S7). Similar results are obtained with using the SIA-only models ( Fig. S6b and d). The smaller volume in the Pliocene Ice-PRISM3 experiment with Bedmap2 hence display the difficulty of the ISMs used here to simulate a readvance of the EAIS to the present-day position in the Wilkes and Aurora basins, compared to the same experiments with Bedmap1.

Discussion
For the control simulations, all ISMs reconstruct an ice sheet close to its PD configuration and result in a smaller equilibrated ice sheet driven by the ERA-40/WOD-09 climate (Fig. 4), mainly due to the drier conditions across East Antarctica in ERA-40 relative to the pre-industrial simulation of HadCM3. The differences between the models for all experiments are rather small considering that all models are used with their own setup for determining ice temperatures and velocities. The variability within the model ensemble is a combination of the differences in SMB and ice fluxes. The average SMB for the six ISMs is 2113.3 ± 129.7 Gt yr −1 (Gt = 10 12 kg) and the ice flux across the grounding line is 346.5 ± 147.8 Gt yr −1 at the final step of each 100 kyr simulation. The ice fluxes are calculated afterwards, using final ice thickness and vertically averaged velocities, and are merely to give an indication of the range between models. Most ISMs for the two Pliocene simulations have a similar final steady-state topography (Figs. 5, 6 and 7) and show a retreat of the WAIS from its PD configuration in the Pliocene Ice-PD simulation due to higher ocean temperatures that enhance sub-shelf melting and ablation due to higher surface temperatures (Fig. S1a and c). Only RIMBAY simulates a WAIS that is larger than PD (Fig. S5) and the PSU-ISM simulates a rather large WAIS as well, although the ice shelves have completely disintegrated. For the EAIS all models produce a similar final surface topography and final volume for the two simulations. Here most ISMs do show an increase in ice volume for the Pliocene Ice-PD simulation relative to PD, mainly caused by a higher accumulation in the Wilkes and Aurora basins relative to the pre-industrial (see Fig. 1j). The largest difference in grounded ice volume is simulated by RIMBAY (Fig. S5). Although surface-air and ocean temperatures are largely the same for each ISM at the initial step of the 100 kyr simulations, the ISMs show quite a strong differences between the two Pliocene simulations ( Fig. 2c and d). The divergent behaviour within our intercomparison is largely due to differences in ice fluxes and sub-shelf melting (not shown), i.e. two features in the models that are not constrained in our experimental setup.
When comparing the bedrock topography of Bedmap1 and Bedmap2, the difference for the subglacial basins in the Wilkes and Aurora regions is in particular large (Fig. 11h) and of importance for ice-sheet growth and collapse. The re-advance and retreat of the EAIS in these regions involves marine ice-sheet instability (Weertman, 1974), related to grounding-line migration, hysteresis due to unstable retreat on reverse bed slopes, pinning on bedrock highs and variable buttressing due to ice shelves (Schoof, 2007;Jamieson et al., 2012;Gladstone et al., 2012;Parizek et al., 2013;Docquier et al., 2014). Although some models used here do include physical parameterisations to improve grounding-line migration (for example PISM, PSU-ISM and RIMBAY are used in Pattyn et al., 2013), the SIA-SSA models used here are not capable of simulating a retreat from modern extent for the given Pliocene climate forcing (Figs. 8 and S7). This could also be ascribed to a too course resolution at the grounding line, which is known to be of crucial importance for accurate simulations of the migration of the grounding line accurately (e.g. Vieli and Payne, 2005;Goldberg et al., 2009). A fixed grid size of 40 km was used here because of computational efficiency, whereas long-term simulations with higher-order physics are currently infeasible.

Comparison with SIA-only ISMs
The initial setup of PLISMIP was comprised of models that include the SIA only (Dolan et al., 2012), similar to the experiments performed for Greenland . Although a combination of the SIA and SSA is necessary to simulate the complete domain of the AIS, the main driver of ice flow for the EAIS is the SIA-based ice flow velocity. Here, we compare simulations with three SIA ISMs to the modelled EAIS with the SIA-SSA models. The three models are IcIES (Saito and Abe-Ouchi, 2004), BASISM (Hindmarsh, 2001) and IMAU-ICE, a SIA version of AN-ICE (de Boer et al., 2013). All three models use the SIA as described in Appendix A1 and use a Weertman type sliding law and have been used for the Greenland experiments Table 4. Contribution to the mean sea level (m) for all simulations relative to the initial PD ice sheet using Bedmap1 (B1) and Bedmap2 (B2). Contributions are calculated with a constant ocean area of 3.62 × 10 14 m 2 using only grounded ice above flotation and corrected for changes in bedrock elevation (see text for details). Differences for the three Pliocene experiments with Control HadCM3 are shown in brackets. as well, as described in Koenig et al. (2015). IMAU-ICE is largely similar to ANICE, only uses Weertman sliding. As is shown in Fig. S5c, final ice volume for the EAIS falls within the range of the SIA-SSA models, with IcIES on the low end and BASISM on the high end of the spectrum of SIA-SSA models. Similar to the six SIA-SSA models with Bedmap2, the three SIA-only models all show a smaller ice extent over the Wilkes and Aurora basins ( Fig. 11i-k). Also, all three models simulate a smaller ice volume using Bedmap2 ( Fig. S6b and d) relative to Bedmap1 ( Fig. 7b and d). The smaller EAIS as simulated in the Pliocene Ice-PRISM3 experiment with Bedmap2 does not simply reflect a decrease in SMB, but can largely be ascribed to the difference in bedrock topography and pinning on bedrock highs in the Wilkes and Aurora basin. The latter is also demonstrated by the fact that the SIA-only models show similar results as the SIA-SSA ISMs. Therefore, the initial state of the ice sheet is here of deciding influence and a retreat/re-advance could potentially be approximated in a more realistic fashion by models combining parameterisations of the above-mentioned grounding-line physics with sufficiently high spatial resolution (e.g. Cornford et al., 2013).

Contribution to late-Pliocene sea level
The contributions of Antarctica to sea level are shown in Table 4. All values are derived from the total ice volume at the last time step of each 100 kyr simulation relative to the PD mapped ice sheet on the 40 by 40 km grid, using ice thickness above flotation and a correction for bedrock change: where Hi 0af and Hi af are the ice thickness above flotation for the initial (either Bedmap1 or Bedmap2) and final modelled ice sheet in m water equivalent, respectively: Density of ice and seawater are taken as provided in Table 2 and an ocean area of O area = 3.62 × 10 14 m 2 is used. Hi is ice thickness (in m) and Hb is the bedrock topography (in m; negative below sea level). On average, the six ISMs calculate a sea-level contribution of −3.23 ± 2.93 m s.e. for the Control HadCM3 simulation relative to Bedmap1. On the contrary, for the Bedmap2 simulations the ISMs produce a difference with Bedmap2 ice of 1.11 ± 3.02 m s.e., i.e. on average too small for the Control HadCM3 simulation (  Table 3) was done different for each ISM. Surface-air and ocean temperatures have initially similar values for each ISM, but differences in final ice volumes calculated by the ISMs can be largely ascribed to differences in ice fluxes across the grounding line and subshelf melting. The PDD melt scheme seems less sensitive to initial conditions relative to the ITM scheme used by ANICE and leads to ice coverage over the WAIS land areas, also in the Pliocene Ice-PRISM3 experiments for the other five ISMs.
The small spread in final ice volume of the AISM can be ascribed to the spin-up procedure employed prior to the 100 kyr steady state simulations. Changes of ANICE are largely mass balance driven, whereas a smaller ice volume corresponds to more melt underneath the ice shelves and a larger volume because of increased precipitation. The results of PISM and ANICE are comparable, whereas a similar procedure is employed for calculating the sliding using the SSA. Sub-shelf melting of PISM is relatively more sensitive compared to the other models (see Fig. S4). The largest re-advance of all ISMs is shown by the PSU-ISM, due to the grounding line parameterisation including in the model. For the Pliocene Ice-PD experiment, RIMBAY does not show a large retreat of the WAIS compared to the other ISMs, possibly related to the similar enhancement factors used for SIA and SSA flow. SICOPOLIS includes an additional grounding-line melt parameterisation. Although the Control experiments show a volume close to that of the initial ice sheet, the grounding line has retreated significantly. Moreover the melt scheme also induces a much further retreat for the Pliocene PD-ANT experiment with Bedmap1.
Our simulations of the late-Pliocene warm period with Bedmap1 do not support a significant reduction of the EAIS across the Wilkes and Aurora basins as has been suggested by studies of marine sediments (e.g. Williams et al., 2010;Cook et al., 2013). Although the experiments using Bedmap2 with the initial PRISM3 ice-sheet geometry do simulate less ice cover over the Wilkes and Aurora basin, with imposing the Pliocene climate on a modern-day Antarctic ice sheet, a significant retreat from the current grounding line positions could not be realised (Figs. 8 and S7). Our experiments do show the importance of using an accurate data set of bedrock topography for ice-sheet models. Additionally, improved grounding-line and/or calving physics may be crucial to simulate the retreat of the grounding-line for subglacial basins of Antarctica (e.g. Fogwill et al., 2014;Mengel and Levermann, 2014;Pollard et al., 2015).
The spread in the PLISMIP-ANT model ensemble is quite considerable. Additionally, incorporating climate fields from different AO-GCMs  will allow evaluation of the uncertainties in climate forcing on the steady state response of the modelled AIS. For the Greenland ice sheet a previous intercomparison showed that this is important to take into account . The spread in surface-air and ocean temperatures between AO-GCMs over Antarctica is considerable  and it is likely that both calculated SMB and sub-shelf melting will also contribute to a large spread in the modelled AIS sealevel contribution. Moreover, it would be desirably that in a future study the design of the experiments is more restricted and that the ISMs uses the same SMB and sub-shelf melting scheme, such that model intercomparison can focus on ice-sheet dynamics alone.
In conclusion, results are depending on the ISMs, the imposed climate forcing and initial conditions of the ice sheet. The sea-level contributions from the Pliocene Ice-PD experiment relative to the Control HadCM3 simulations are on average 0.0 ± 6.1 m s.e. for Bedmap1 and −1.9 ± 4.9 m s.e. for Bedmap2. Although some models, e.g. ANICE and PISM, do exhibit a positive contribution for both the Bedmap1 and Bedmap2 experiments (Table 4). For the Pliocene PD-Ant experiment, all ISMs simulate an increase in ice volume, with an average sea-level drop of −3.7 ± 2.2 m s.e. for Bedmap1 and −1.5 ± 1.5 m s.e. for Bedmap2. On the other hand, positive sea-level contributions are simulated for the Pliocene Ice-PRISM3 experiments; however, these simulations are largely dependent on the initial ice-sheet topography. Hence, we cannot provide a conclusive sea-level contribution of the Antarctic ice sheet to Pliocene sea-level rise. All ISMs used within PLISMIP-ANT apply the shallow ice and shallow shelf (or shelfy stream) approximations to reduce computational time relative to solving the full Stokes equations of flow. Here we shortly describe the two approximations.

A1 The shallow ice approximation (SIA)
For modelling 3-D ice sheets over long time scales, the SIA is commonly used to calculate ice flow over land areas (e.g. Hutter, 1983;Huybrechts, 1990). For the SIA the normal, longitudinal, stresses are neglected relative to the horizontal shear stress. In this way, shearing stresses induced by vertical changes of the horizontal velocities are only balanced by the driving stress: ρ i gH ∇H s . The SIA velocities follow from an integral equation: Here, ∇H s is the horizontal surface slope, ζ the scaled vertical coordinate, ρ i = 910 kg m −3 the density of ice, g = 9.81 m s −2 the gravity acceleration and n = 3 the flow exponent in Glen's flow law. A(T * ) is the flow-rate factor (Pa −3 yr −1 ) depending on the ice temperature corrected for pressure melting dependent (T * ). E SIA is the flow enhancement factor (Ma et al., 2010), which is different for each ISM (see Table 3).

A2 The shallow shelf approximation (SSA)
To determine ice velocities for ice shelves, the approximate stress balance for the SSA includes longitudinal stress which are more dominant compared to the shear stress. The balance equations determine stretching velocities, i.e. the change of the horizontal velocities independent of depth in the horizontal plane. The SSA is largely used to calculate the velocities of ice shelves and ice streams (e.g. Morland, 1987;Bueler and Brown, 2009 Here, u and v are the SSA velocities in the x and y direction, respectively (in m yr −1 ) and H i is the ice thickness. For the SSA the stresses due to stretching are balanced by the gravitational driving stress and possibly the basal stresses τ b,x and τ b,y (in Pa) when applied on land. The parameter µ is the vertical averaged viscosity, a function of the strain rates and the vertical mean flow rate factor A(T * ) (e.g. Bueler and Brown, 2009 with A the vertical mean flow rate factor A(T * ) and E SSA the enhancement factor for the SSA velocities (Ma et al., 2010), which is different for each ISM (see Table 3).