Scaling of instability time-scales of Antarctic outlet glaciers based on one-dimensional similitude analysis

Recent observations and ice-dynamic modeling suggest that a marine ice sheet instability (MISI) might have been triggered in West Antarctica. The corresponding outlet glaciers, Pine Island Glacier (PIG) and Thwaites Glacier (TG), showed significant retreat during at least the last two decades. While other regions in Antarctica have the topographic predisposition for the same kind of instability, it is so far unclear how fast these instabilities would unfold if they were initiated. Here we employ the concept of similitude to estimate the characteristic time scales of several potentially MISI-prone outlet glaciers 5 around the Antarctic coast. Our results suggest that TG and PIG have the fastest response time of all investigated outlets, with TG responding about 1.25 to 2 times as fast as PIG, while other outlets around Antarctica would be up to ten times slower if destabilized. These results have to be viewed in light of the strong assumptions made in their derivation. These include the absence of ice-shelf buttressing, the one-dimensionality of the approach and the uncertainty of the available data. We argue however that the current topographic situation and the physical conditions of the MISI-prone outlet glaciers carry the 10 information of their respective time scale and that this information can be partially extracted through a similitude analysis.


Introduction
Sea-level rise poses a future challenge for coastal regions worldwide (IPCC, WG II, 2014).The contribution from mass loss of the West Antarctic Ice Sheet has been increasing over the last two decades (Medley et al., 2014;The IMBIE team, 2018).Large parts of the Antarctic Ice Sheet rest on a retrograde marine bed, i.e., on a bed topography that lies below the current sea level and is down-sloping inland (Bentley et al., 1960;Ross et al., 2012;Cook and Swift, 2012;Fretwell et al., 2013).This topographic situation makes these parts of the ice sheet prone to a so-called marine ice-sheet instability (MISI, Weertman, 1974;Mercer, 1978;Schoof, 2007;Pattyn, 2018).Currently, this kind of instability constitutes a large uncertainty in projections of future sea-level rise (IPCC, WG I, 2013;Joughin and Alley, 2011;Huybrechts et al., 2011;Golledge et al., 2015;Winkelmann et al., 2015;DeConto and Pollard, 2016;Pattyn et al., 2018): If the grounding line (that separates grounded from floating ice) enters a region of retrograde bed slope, a positive ice-loss feedback can be initiated.Resulting self-sustained retreat, acceleration and discharge of the ice sheet can be hindered by the buttressing effect of ice shelves and topographic features (Dupont and Alley, 2005;Goldberg et al., 2009;Gudmundsson et al., 2012;Favier et al., 2012;Asay-Davis et al., 2016) or strong basal friction (Joughin et al., 2009;Ritz et al., 2015).
The ongoing retreat of two major outlets of West Antarctica, Pine Island Glacier (PIG) and Thwaites Glacier (TG), initiated by warm-water entrainment into their ice-shelf cavities and resulting basal melting (Jenkins et al., 2010;Pritchard et al., 2012;Paolo et al., 2015;Shean et al., 2018), is pointing toward a developing MISI (Rignot et al., 2014;Mouginot et al., 2014;Konrad et al., 2018;Favier et al., 2014;Joughin et al., 2014;Seroussi et al., 2017) with the potential of raising global sea level by more than 3 m (Bamber et al., 2009;Feldmann and Levermann, 2015).Whether the observed retreat indeed marks the start of a MISI or is a temporally limited response to oceanic variability still requires further studies (Hillenbrand et al., 2017;Smith et al., 2017;Jenkins et al., 2018).Similar warm-water exposition might apply to other Antarctic ice shelves in the future (Hellmer et al., 2012;Timmermann and Hellmer, 2013;Greene et al., 2017), bearing the potential to trigger unstable grounding-line retreat.It has been demonstrated that East Antarctica's Wilkes Subglacial Basin (WSB), would, once destabilized, contribute by 3 − 4 m to global sea level (Mengel and Levermann, 2014).The marine part of the adjacent Aurora Subglacial Basin stores ice of around 3.5 m sea-level equivalent (Greenbaum et al., 2015).Other studies investigated the response of the Filchner-Ronne tributaries to enhanced melting (Wright et al., 2014;Thoma et al., 2015;Mengel et al., 2016) and suggest that instability might not unfold there, possibly due to the stabilizing buttressing effect of the large ice shelf and narrow bed troughs (Dupont and Alley, 2005;Goldberg et al., 2009;Gudmundsson et al., 2012).In total, the marine catchment basins of Antarctica that are connected to the ocean store ice masses of about 20 m of sea-level equivalent (Fig. A1).Besides the question of whether local instabilities are already triggered, or will be triggered in the future, the question about which time scale a potentially unstable retreat would evolve if it had been initiated still remains.Naturally, numerical ice sheet models are used to investigate these instabilities.
Here we try to contribute information using a different approach which is based the concept of similitude (Buckingham, 1914;Rayleigh, 1915;Macagno, 1971;Szücs, 1980).To this end the presented method combines observational and model data with a similarity analysis of the governing equations for shallow ice-stream flow.This similarity analysis has been carried out in previous work (Feldmann and Levermann, 2016), yielding scaling laws that determine how the geometry, the time scale and the involved physical parameters for ice softness, surface mass balance and basal friction have to relate in order to satisfy similitude between different ice sheets.The approach is an analogy to, e.g., the derivation of the Reynold's number (Reynolds, 1883) from the Navier-Stokes equation (Kundu et al., 2012), which can provide a scaling law that assures similar flow patterns of a fluid (laminar or turbulent) under variation of its characteristic geometric dimension, velocity and viscosity.In our previous paper we have shown that the ice-sheet scaling behavior predicted according to the analytically derived scaling laws agrees with results from idealized three-dimensional, numerical simulations.In particular, this included the prediction of MISI-evolution time scales over a range of three orders of magnitude for different ice-sheet configurations respecting ice-sheet similarity.In the present study we take the next step and apply these scaling laws to the real world, aiming to infer relative response times to potential destabilization of eleven MISI prone Antarctic outlets, based on their dynamic and geometric similarity (ice streams on marine, landward down-sloping bed topography).This approach has the disadvantage of reduced dynamic complexity compared to numerical modeling, but the advantage of being based on the observed ice-sheet configuration which is generally not the case for an ice-sheet simulation.
Our approach is based on the fact that the currently observed shape of the considered Antarctic outlet glaciers represents the balance of the different forces that act on the ice.That is, the balance of the forces caused by the bed and ice topography, the surface mass balance and the frictional forces from the ice interaction with itself and its surrounding.This static situation carries the information for the initial ice-sheet response after a potential destabilization.The scaling laws used to infer these responsetime scales come out of a similitude analysis of 1) the shallow-shelf approximation (SSA, e.g., as described by Morland, 1987;MacAyeal, 1989) of the Stokes stress balance and 2) mass conservation (Greve and Blatter, 2009).
The SSA is a simplified version of the Stokes stress balance, accounting for the case of membrane stresses dominating over vertical shear stresses.It is thus a suitable representation of the plug-like ice flow observed close to the grounding lines of the ice streams analyzed here.The equation of mass conservation provides the time evolution of the ice thickness, being determined by the horizontal ice-flux divergence and the surface mass balance.
The derivation of the scaling laws from these two equations is presented in full detail in Feldmann and Levermann (2016) and thus shall be outlined only very briefly here: in their dimensionless form the SSA and the equation of mass conservation together have three independent numbers, analogous to the Reynolds number in the Navier-Stokes equation (Reynolds, 1883).
If each of these three numbers is the same for topographically similar situations (i.e., ice streams on retrograde bed) the solution of the ice-dynamical equations will be similar, satisfying three independent scaling laws (one for each dimensionless number).
These numbers and thus the scaling laws are combinations of parameters/scales of characteristics of the ice dynamics, physical ice properties, boundary conditions and geometry.In the present study we use scales for the bed elevation, b, the slope of the retrograde bed, s, the surface mass balance a and the bed friction C, obtained from observational and model data (for details see Sec. 2.2).
The procedure to infer these scales from the data is schematized in Fig. 1: For each of the tributaries considered in this study a transect is defined that represents a center line of the ice stream and covers its potentially unstable bed section close to the grounding line.Since this study aims at providing a scaling analysis the transect were not chosen following a strict definition that involves the velocity or topography fields.They were chosen as to be straight lines that allow the quantification of a generic slope of the topography.The results presented here do not depend on the precise choice of the position of these transects.Fitting a linear slope to the retrograde bed section (blue line), the bed elevation b at the starting point of the slope and the slope magnitude s serve as characteristic geometric scales for the sloping bed on which the ice sheet is grounded.Thus, b is a representative scale also for the thickness of the ice below sea level in the vicinity of the grounding line.Basal friction C is averaged along this section to represent local conditions and the surface mass balance a is obtained from averaging over the catchment that feeds the ice stream (Figs. 2 and S1-S9).For each tributary these physical parameters are taken relative to the reference values of PIG (primed vs. unprimed parameters).The resulting scaling ratios are then used to calculate time-scale ratios τ = T /T via two independent scaling conditions coming out of the similitude principle.We assume that the inferred time scales correspond to the outlet-specific initial response time to potential destabilization.We choose PIG as the reference as it is one of the most prominent tributaries of West Antarctica and it is relatively well observed.Due to the nature of the conducted scaling analysis the scales calculated here could also be expressed relative to any other of the examined oulets without changing the results.

Scaling laws and uncertainty criteria
The first of the two scaling laws used in this study is derived from the isothermal, two-dimensional SSA (Feldmann and Levermann, 2016, , Eq. 1).More precisely, the scaling condition comes out of the non-dimensionalized frictional term of the SSA (see Eqs. 9 and 11 in Feldmann and Levermann, 2016).In its general form it reads Levermann, 2016).The friction data applied in this study is obtained from an inversion model that assumes a linear relation between basal friction and ice velocity (m = 1) with spatially varying proportionality constant C (Morlighem et al., 2013).For this linear case Eq. ( 1) simplifies to As a consequence, for example, more slippery bed conditions (lower γ) and a steeper bed slope (higher σ) yield a shorter time scale (smaller response time).
The second scaling law (τ ) used here results from mass conservation (Feldmann and Levermann, 2016, , Eqs. 4 and 13).It is independent of basal properties but a function of the vertical scaling β and the surface-mass-balance ratio δ = a /a, i.e., τ = βδ −1 . (3) According to Eq. ( 3), an initially stable situation with stronger accumulation but initially thinner ice at the grounding line, for example, results in a faster response in case of destabilization.As detailed in (Feldmann and Levermann, 2016), the above scaling laws (Eqs. 2 and 3) are consistent with analytic solutions of the ice-dynamic equations (Schoof, 2007).
Calculation of the time-scaling ratios τ and τ provides two independent estimates for the response time of each outlet examined here, values which are relative to the reference PIG.We use the deviation between the two values as a measure for the uncertainty of the estimation, accounting for observational uncertainty and approximations in our approach.In order to consider the calculated response-time scaling of a tributary as of sufficient certainty we require the two calculated responsetime ratios to fulfill the following conditions: The first criterion compares the obtained qualitative scaling behavior which would be contradictory for c 1 < 0 (one timescaling ratio would indicate a faster, the other a slower response compared to the reference).The second criterion yields the relative error between the two calculated response-time ratios with respect to their mean and ensures that this uncertainty is at most 20 %.Ice streams with a larger uncertainty in their calculated time scaling are reported in the following, but considered unsuitable for the approach used here.

Data
The locations of the chosen transects are motivated by the bed topography and the flow field of the analyzed outlets to capture their retrograde bed section directly upstream of the grounding line (e.g., see Fig. 2).The ice-stream characteristic parameter values are all obtained from datasets that represent present-day conditions of the Antarctic Ice Sheet.The bed topography stems from the BEDMAP2 dataset (Fretwell et al., 2013) which is the most recent continent-wide compilation of Antarctic ice thickness and basal topography and involves data from various sources, including satellite, airborne radar, over-snow radar and seismic-sounding measurements that total 25 million data points.The basal-friction dataset is taken from Morlighem et al.
(2013) and is a product from inversion of observed present-day Antarctic ice-surface velocity (Rignot et al., 2011) with the Ice Sheet System Model (Larour et al., 2012).The thermomechanical, higher-order model uses anisotropic mesh refinement that allows for horizontal resolutions down to 3 km along the Antarctic coast, i.e., our region of interest.The field of basal friction is consistent with results from other inversion models (Joughin et al., 2009;Morlighem et al., 2010;Pollard and DeConto, 2012).
Surface mass balance is obtained from the Regional Atmospheric Climate Model (RACMO version 2.3p2, van Wessem et al., 2014) and is averaged over the period from 1979 to 2016.RACMO is forced by ERA-Interim reanalysis data at the lateral boundaries, simulating the interaction of the ice sheet with its atmospheric environment, involving relevant processes such as solid precipitation, snow sublimation, and surface meltwater run-off.For most of the tributaries the surface mass balance a is obtained from averaging over the entire feeding catchment basin (Zwally et al., 2012).For some basins that have a long coast line and are drained by several major tributaries the averaging area is constrained to a region upstream of the ice stream of interest (Figs.S3-S6 and S9).

Results
The examined individual retrograde bed slopes vary in magnitude by a factor of six and bed elevation differs by up to 1000 m (Figs.3a and A2a).The lengths of the slope sections indicate how far the considered retrograde slopes reach inland before an (intermediate) section with vanishing slope or up-sloping bed follows (Figs. 2 and S1-S9).We scale these bed sections with respect to the reference according to the obtained vertical and horizontal length-scale ratios (Table 1).The resulting bed profiles collapse towards the reference showing similar down-sloping, while still exhibiting their characteristic pattern (Figs.3b and   A2b).This assures that the chosen geometric measures adequately reflect the characteristics of the retrograde bed slopes.
Within the ensemble of analysed outlets TG has the smallest response time (Fig. 4), being between 1.25 to almost 2 times as fast as PIG, which is found to be the second fastest.For the chosen outlet of WSB the two calculated ratios indicate a response that would be twice as slow as PIG.East Lambert Rift (ELR) and Institute Ice Stream (IIS) that feed Amery Ice Shelf and Filchner-Ronne Ice Shelf, respectively, are found to be 4 to 6 times slower than PIG.By far the slowest response shows Support Force Glacier (SFG) which is (more than) 10 times slower than PIG.
For each of these six tributaries (Table 1, blue) both imposed quality criteria (Eqs.4 and 5) are fulfilled with a maximum error of c 2 = 20 %.The other five regions are discarded since one or both of the criteria are not met, with c 2 ranging from 46 to 98 % (Table 1, grey).Each of the discarded ice streams exhibits extremely slippery bed conditions, with basal-friction values of up to three orders of magnitude smaller than the reference (Table 1).On the other hand these outlets also have the smallest bed slopes of the ensemble.Both slope and friction are hence subject to large relative uncertainties that mutually amplify in the calculation of τ and hence amplify the mismatch with τ .
Throughout the successfully examined outlets, diagnosed basal friction and surface accumulation each differ by one order of magnitude.Most slippery conditions and highest snowfall are found in the Amundsen Sea Sector, both favoring the short time scales of PIG and TG.We also infer the scaling of the ice softness A using a third scaling law coming out of the similarity analysis from Feldmann and Levermann (2016) (see Appendix A).The resulting ratios differ by two orders of magnitude.For most of the ice streams discussed here a lower softness value corresponds to a larger time scale.In other words, stiffer ice tends to slow down the response to destabilzation.

Discussion
According to observations, PIG and TG show the largest ice discharge rates in the Amundsen Sea Sector of West Antarctica since the last four decades, including phases of rapid ice speedup, dynamic thinning and grounding-line retreat (Rignot et al., 2014;Mouginot et al., 2014;Konrad et al., 2018).PIG's most recent acceleration was found to coincide with the ungrounding of a region of only lightly grounded ice and the grounding line that is now at the upstream end of this area stabilized during the last few years (Mouginot et al., 2014).Meanwhile PIG's deceleration in ice discharge is more than compensated by TG's speedup which can be attributed to reduced ice-shelf buttressing that previously stabilized the upstream grounded ice (Rignot et al., 2014).It thus remains an open question whether the current destabilization of TG will result in retreat rates larger than PIG's recent retreat rates and hence confirm our finding of a faster time scale of TG compared to PIG.Regional model simulations of these regions (Favier et al., 2014;Joughin et al., 2014) indicate that TG's discharge rate is up to a factor of two larger than PIG's after the basal-melt perturbation has ceased and the glaciers have relaxed into a configuration of constant ice discharge.
The relatively large time scales of SFG, IIS and ELR, calculated here, are consistent with results from dynamical modeling that find these outlets of Filchner-Ronne Ice Shelf and Amery Ice Shelf to be far less sensitive to instability (Gong et al., 2014;Wright et al., 2014;Thoma et al., 2015;Ritz et al., 2015;Mengel et al., 2016) than WSB (Mengel and Levermann, 2014) and PIG and TG (Favier et al., 2014;Joughin et al., 2014;Feldmann and Levermann, 2015;Seroussi et al., 2017).In these cases the buttressing effect of the large ice shelves and the narrow, channel-type bed topography confining these outlets as well as relatively high basal friction (Table 1) might prevent MISI initiation and limit grounding-line retreat.
Our results are obtained from the application of scaling laws, a procedure commonly pursued in hydrodynamics and engineering (e.g., Scruton, 1961;Li et al., 2013) and also to some extent in the field of glaciology (e.g., Burton et al., 2012;Corti et al., 2014).The two scaling laws used here (Eqs. 2 and 3) are derived from an approximation of the Stokes stress balance (Feldmann and Levermann, 2016).The approximation accounts for the dynamics of shallow ice streams that are characterized by dominating longitudinal stresses and rapid sliding and typically have an onset of several 10 to 100 km upstream of the grounding line.Our analysis focuses precisely on such regions.The examined sections are all situated in the proximity of the grounding line and indeed exhibit slippery bed conditions and fast ice flow as indicated by inferred basal friction and velocity data, respectively (Figs. 2 S1-S9 f and g).
The proposed method is based on the assumption that the ice-stream specific properties can each be represented by a single scale.This is realized by averaging spatially varying, complex fields along a representative transect (ice softness, basal friction) or over the entire catchment basin (surface mass balance).While such an approach means a substantial idealization of the complex nature of ice dynamics it allows for a timescale analysis that is based on representative ice-property scales and considers the relevant physics through the (simplified) stress balance.Our analysis will only be valid as long as this representation is sufficient.
Furthermore, our approach assumes idealized conditions of unstable retreat.Though taking into account bed friction as a possible stabilizing factor (Ritz et al., 2015) it is strongly limited by not explicitely accounting for the stabilizing effect of ice-shelf buttressing (Dupont and Alley, 2005;Gudmundsson, 2013).That is, the approach estimates the time scale of retreat only after no significant ice shelf is present anymore or for the situation of an unconfined ice-shelf that does not provide any backstress (e.g., as might be the case for the present-day TG ice shelf, Rignot et al., 2014).Furthermore, our investigation is constrained to the short initial time period (decades to centuries) after destabilization during which the grounding line passes the investigated retrograde bed-slope section near the coast.
The used datasets of observed bed and ice geometry still show substantial gaps in coverage, though the regions analyzed here lie at the lower end of the uncertainty range (Fretwell et al., 2013, their Figs. 11 and 12).Measurements of the Antarctic surface mass balance are in general sparsely distributed over the continent.In contrast, results from numerical modeling are known to have biases with respect to observations.To reduce uncertainty we use a surface-mass-balance dataset obtained from a regional climate model that was calibrated by the available observational data (van de Berg et al., 2006;van Wessem et al., 2014).

Conclusions
Combining observations and modeling data with a similarity analysis of the governing, simplified ice-dynamic equations we apply the similitude principle (Fig. 1) to eleven MISI-prone Antarctic outlets to infer their relative time scales (Fig. 4).
Assuming potential destabilization of these outlets (MISI triggered) in the absence of ice-shelf buttressing our results suggest that PIG and TG have the smallest response time to destabilization (Fig. 4, Table 1), with TG responding between 1.25 to almost 2 times as fast as PIG.Further considered ice streams draining East or West Antarctica are found to respond twice (WSB) to 10 times (SFG) as slow as PIG.This also suggests that the dynamic regime of PIG is more similar to TG than to any other of the outlets analyzed here.This way, the proposed scaling approach may help field workers in their decision on which glaciers to observe to be able to study a broad and thus insightful spectrum of different (dissimilar) glacier dynamics.
The presented analysis can not make any statements about absolute time scales of the instabilities.It is constraint to relative statements between the different outlet glaciers, based on their specific geometry and physical conditions which carry information of the specific dynamic balance within the respective outlet glacier.Regarding absolute response times of individual tributaries, continued observation of the retreat of PIG and TG that extend the already available data might allow to deduce long-term trends, thus giving estimates for the absolute time scales of instability.These could be used to calibrate the method presented here.Our method comes with a large number of assumptions and limitations.Uncertainties will reduce as observational data will further improve and the relevance of the results will increase with sophistication of the method (e.g., explicit Blue variables are obtained from observations and model data.The equations for τ (Eqs. 2 and 3) result from a similarity analysis of the stress balance and conservation of mass (Feldmann and Levermann, 2016).shows the bed geometry of the retrograde slope section, scaled with vertical and horizontal scaling factors (Table 1) and normalized to the dimensions of Pine Island Glacier.Maps of (d) surface mass balance (van Wessem et al., 2014), (e) bed topography (Fretwell et al., 2013), (f) surface velocity (Rignot et al., 2011) and (g) basal friction (Morlighem et al., 2013).Drainage basins are obtained from Zwally et al. (2012).
The transect as shown in panels a and b are provided as black and yellow lines in panels d-g.(b) Bed topography from panel (a) scaled according to similitude theory using the horizontal and vertical length-scale ratios obtained from observations (see Table 1) and normalized to the dimensions of the PIG slope (that has the same aspect ratio in both panels).For better visibility the x and y axis are cut (see Fig. S4 for full extent of the ELR bed section).Marine bed topography in blue (colorbar), ice shelves in yellow and continental shelf in grey are taken from BEDMAP2 (Fretwell et al., 2013).
Table 1.Scaling ratios for all examined Antarctic outlets, clockwise, starting with the reference Pine Island Glacier.The ratios for the vertical scale (β), the retrograde-slope magnitude (σ), the basal friction (γ) and the surface mass balance friction (δ) are inferred from observational and model data (see Methods).The horizontal scaling ratio α, used for the scaling of the retrograde bed topography (Figs. 3 and A2), is calculated via α = βσ −1 .The scaling of the ice softness ζ is obtained from Eq. (A2).The response-time ratios τ and τ are calculated according to the scaling laws given by Eqs.
(2) and (3).The last two columns show the resulting values from application of the two quality criteria (Eqs.4 and 5) according to which tributaries are either accepted (blue) or discarded (grey).Abbreviations for the outlets as in Fig. 4. involving the ice softness ratio ζ = A /A. Due the absence of Antarctic ice-softness data we cannot utilize this scaling law to calculate a third independend time-scale ratio.However, combining this third scaling law Eqs. 2 with the first one and A1 we eliminate the time-scale ratio to calculate the ice-softness ratio via Being dependent on the scaling of the geometry and the basal friction, ζ can be considered as the ratio between ice-softness 10 scales that are representative for the analyzed outlets cross-sections.
being a function of the slope scaling ratio σ = s /s, the vertical scaling factor β = b /b and the scaling ratio of basal friction γ = C /C.The friction exponent m stems from the Weertman-type sliding law that enters the SSA (see Eq. 3 in Feldmann and

Figure 1 .
Figure1.Schematic visualizing the presented approach to obtain the response-time ratio τ between different Antarctic outlets as a function of their specific physical properties, here exemplified for Pine Island Glacier (unprimed reference) and Thwaites Glacier (primed scales).

Figure 2 .
Figure 2. Scaling of bed geometry and maps of data involved in the response-time calculation (here shown exemplarily for Thwaites Glacier (TG) with respect to the reference Pine Island Glacier (PIG)).(a) and (b) Elevation of ice surface (grey) and bedrock (black) with retrograde slope section (colored) that is used to infer the vertical scale and the slope magnitude for the transects shown in panels (d) -(g).Panel (c)

Figure 3 .
Figure 3. (a) Unscaled bed topography of retrograde slope sections of the six examined Antarctic outlets that fulfill the time-scaling criteria.

Figure A1 .Figure A2 .
Figure A1.Ice thickness above flotation of marine regions (bedrock below sea level) that are connected to the ocean (colorbar).Ice shelves in blue and continental shelf in grey.Ice and bed topography are taken from the BEDMAP2 compilation(Fretwell et al., 2013).