Similitude of ice dynamics against scaling of geometry and physical parameters

The concept of similitude is commonly employed in the fields of fluid dynamics and engineering but rarely used in cryospheric research. Here we apply this method to the problem of ice flow to examine the dynamic similitude of isothermal ice sheets in shallow-shelf approximation against the scaling of their geometry and physical parameters. Carrying out a dimensional analysis of the stress balance we obtain dimensionless numbers that characterize the flow. Requiring that these numbers remain the same under scaling we obtain conditions that relate the geometric scaling factors, the parameters for the ice softness, surface mass balance and basal friction as well as the ice-sheet intrinsic response time to each other. We demonstrate that these scaling laws are the same for both the (two-dimensional) flow-line case and the three-dimensional case. The theoretically predicted ice-sheet scaling behavior agrees with results from numerical simulations that we conduct in flow-line and three-dimensional conceptual setups. We further investigate analytically the implications of geometric scaling of ice sheets for their response time. With this study we provide a framework which, under several assumptions, allows for a fundamental comparison of the ice-dynamic behavior across different scales. It proves to be useful in the design of conceptual numerical model setups and could also be helpful for designing laboratory glacier experiments. The concept might also be applied to real-world systems, e.g., to examine the response times of glaciers, ice streams or ice sheets to climatic perturbations.


Introduction
In the fields of fluid dynamics and engineering, scaling laws are used to perform experiments with spatially reduced models in water channels or wind tunnels to predict the behavior of the associated full-scale system (e.g., Scruton, 1961;Li et al., 2013).Dimensional analysis and the principle of similitude allow us to derive such scaling laws analytically (e.g., Rayleigh, 1915;Macagno, 1971;Szücs, 1980).For instance, a dimensional analysis of the Navier-Stokes equation (Kundu et al., 2012) yields the Reynold's number (Reynolds, 1883) as one of the dimensionless parameters of the governing equation which characterize the dynamics of fluid flow.Under the assumption of the similitude principle the Reynold's number can provide a scaling law for the fluid's characteristic linear dimension, velocity and viscosity that assures similar flow patterns.The principle of similitude is applied well beyond the field of engineering, e.g., in zoology (land mammals move in dynamically similar fashion at equal Froude number; Alexander and Yayes, 1983) or biology (Stahl, 1962).
The similitude concept is also applied to some extent in the field of glaciology: in laboratory glacier experiments dimensionless numbers like the Reynold's, Froude and Ramberg numbers are used to check for (dynamic) similarity between the geometrically scaled model (based on the properties of the analogue ice material) and the real-world system (Burton et al., 2012;Corti et al., 2014).In Halfar (1983) and Bueler et al. (2005) the similitude principle is used to derive similarity solutions of the shallow-ice approximation (SIA; Hutter, 1983) of the full Stokes stress balance for the case of an isothermal, radially symmetric ice sheet.
Here we apply the concept of similitude to the dynamics of idealized ice sheets based on the shallow-shelf approximation (SSA; Morland, 1987;MacAyeal, 1989;Greve and Blatter, 2009).In particular we assume isothermal ice and a spatially uniform basal friction coefficient, conditions that have been used to analyze ice-sheet dynamics in a number of previous studies (e.g., Dupont and Alley, 2005;Goldberg et al., 2009;Gudmundsson et al., 2012;Pattyn et al., 2013;Asay-Davis et al., 2015).Neglecting the terms of vertical shearing in the stress balance and accounting for the small thicknessto-length ratio of ice sheets, the SSA represents the relevant dynamics of floating ice shelves and grounded ice streams.The capability of numerical SSA models to simulate these regimes that are characterized by fast plug-like flow has been demonstrated in various studies (e.g., Goldberg et al., 2009;Gudmundsson et al., 2012).The SSA can be complemented by the SIA (Huybrechts, 1990;Sato and Greve, 2012) to also include vertical shearing, which is dominant in the more stagnant interior parts of an ice sheet (Bueler and Brown, 2009;Pollard and DeConto, 2012;Thoma et al., 2015), whereas higher-order approximations (Schoof and Hindmarsh, 2010;Larour et al., 2012;Cornford et al., 2015) neglect less stress components in the full Stokes stress balance (Favier et al., 2012).The MISMIP3d benchmark revealed that numerical models applying the SSA can capture grounding-line dynamics comparable to more elaborate models in conceptual experiments (Pattyn et al., 2013;Feldmann et al., 2014).
A dimensional analysis of the ice-dynamic equations is often carried out to compare the magnitudes of the different acting forces and thus to derive physically motivated approximations, as done when deriving the SSA from the full Stokes stress balance (MacAyeal, 1989;Greve and Blatter, 2009).The non-dimensionalized form of the SSA itself and the involved dimensionless coefficients that result from the introduction of typical scales for, e.g., ice-sheet thickness and velocity, have been used to consider asymptotic limits of SSA ice flow in previous work (Schoof, 2007a;Dupont and Alley, 2005;Tsai et al., 2015;Haseloff et al., 2015).In the present study we utilize these coefficients to derive ice-sheet scaling laws for the geometry, response time and other physical icesheet parameters, a step that, to our knowledge, has not been taken before.The scaling behavior of ice sheets, which here is analyzed in a conceptual way, might be of use to better understand the large-scale evolution of the polar ice sheets.Of particular interest is the scaling of the ice-sheet response time (Levermann et al., 2013(Levermann et al., , 2014) ) against the background of Antarctic instabilities (Weertman, 1974;Schoof, 2007b;Rignot et al., 2014;Fogwill et al., 2014;Mengel and Levermann, 2014).The timescales of possible rapid ice discharge due to instability in the past (Pollard and DeConto, 2009;Pollard et al., 2015) and future (Favier et al., 2014;Joughin et al., 2014;Feldmann and Levermann, 2015b) are highly uncertain.
The paper is structured as follows: in the next section the governing equations in SSA are non-dimensionalized to derive ice-sheet scaling laws for one and two horizontal dimensions, respectively.Afterwards the analytically predicted icesheet scaling behavior is compared with results from numerical modeling.To this extent conceptual experiments are designed in two and three spatial dimensions.Steady states as well as the transient response to perturbation of the simulated ice sheet are analyzed for a systematic variation of the scaling parameters which are prescribed according to the scaling laws.We then examine analytically the implications of the scaling conditions for the response times of ice sheets considering the geometric scaling factors and basal friction parameter as independent variables.Eventually we discuss the results and conclude.

Similitude of shallow ice-sheet dynamics
Here we derive scaling laws that determine how the geometry, response time 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.This is visualized conceptually in Fig. 1 for two ice sheets which differ in vertical and horizontal scale.Based on the governing equations in dimensionless form, we obtain dimensionless scale factors which depend on the scales of the geometric and physical parameters of the ice sheet.The requirement that each of these factors has to remain the same under a scaling of the parameters makes sure that the dynamic equations remain exactly the same.The resulting scaling laws thus put constraints on the parameter scaling, ensuring similitude between the different ice-sheet configurations.

Basic equations for similitude analysis
The problem addressed here is the one of an isothermal icesheet in SSA (Greve and Blatter, 2009).The two horizontal components of the stress balance in SSA with spatially uniform ice softness A are given by .Schematic of the similitude-analysis method carried out in this study.A reference system (blue ice sheet and bed topography) with geometric scales h and L, timescale T and physical parameters ice softness A, basal friction coefficient C and surface mass balance a is scaled in horizontal and vertical direction (red contours, primed system).The goal is to derive the scaled parameters of the primed system under which dynamic similarity between both ice sheets holds.A dimensional analysis of the governing equations yields dimensionless scale factors which have to remain the same under scaling to attain similitude.The resulting scaling laws determine the scaled (primed) parameters.
implies that the horizontal velocities do not vary in vertical (z) direction.The effective strain rate ˙ e (Greve and Blatter, 2009) can be written as We choose the basal shear stress in Eq. ( 1), τ b = (τ b,x , τ b,y ), to be given by a Weertman-type sliding law (Greve and Blatter, 2009): with constant friction coefficient C. The exponent m determines the particular type of the sliding law, including plastic (m = 0, magnitude of basal shear stress independent of velocity; Tulaczyk et al., 2000) and linear-viscous (m = 1, basal shear stress proportional to ice velocity; MacAyeal, 1989) behaviors.A value of m = 1/n = 1/3 is commonly assumed to represent sliding over rough bed (Schoof, 2007a;Joughin et al., 2009;Cuffey and Paterson, 2010).The evolution equation for the ice thickness, i.e., the ice thickness equation (ITE), which results out of mass conservation (Greve and Blatter, 2009) reads with horizontal ice flux Q = H v and surface mass balance a.Throughout the study we focus on the grounded part of the ice sheet and assume negligible melting/refreezing at its base.Hence the basal mass balance is not taken into account in Eq. (4).

Flow-line case
In the flow-line case the geometry of an ice sheet can be scaled in horizontal (x) and vertical (z) direction, using two scaling factors α and β, respectively (α, β > 1 for elongation and 0 < α, β < 1 for shortening).We define these as where the prime denotes the scaled system.In particular, Eq. ( 5) states that the ice-sheet length L scales according to L = αL.
Since we neglect the y direction here, we only have to consider the x component of the SSA (Eq.1a), in which all terms that include y drop out.The effective strain rate (Eq.2) thus simplifies to ˙ e = ∂v x ∂x and the SSA reads The ITE (Eq.4) in flow line is given by Now we bring these two equations into nondimensionalized form by introducing the dimensionless variables L , using the scales H, L and T for ice-sheet thickness, length and response time, respectively.We obtain www.the-cryosphere.net/10/1753/2016/The Cryosphere, 10, 1753-1769, 2016 and for the SSA and ITE, respectively.In Eq. ( 9) the two dimensionless constants θ and φ relate the different involved stresses to the driving stress.Extending θ with H/L and φ with L −1 we see that these scale factors relate the membrane stresses (Hindmarsh, 2006) and the basal stresses to the driving stress, respectively.In the floating ice shelf the driving stress is always fully balanced by the membrane stresses (no basal resistance, thus C = 0 in Eq. 9).Focussing on the grounded part of the ice sheet, we assume that its driving stress is balanced by a combination of membrane stresses and basal stresses.
The two governing Eqs. ( 9) and (10) of our problem remain exactly the same as long as each of the dimensionless factors θ , φ and ω remain the same.In other words, the icesheet dynamics are expected to be similar under a transformation that leaves these factors unchanged.Thus the scaling of the ice sheet's typical length and thickness scales according to Eqs. ( 5) and ( 6), i.e., L = αL and H = βH, in general requires (some of) the physical parameters a, C, A and its response time T to change in order to maintain similarity with respect to the unscaled ice sheet.We hence can infer three scaling conditions for the timescale ratio τ = T /T : with friction-coefficient ratio γ = C /C, ice-softness ratio ζ = A /A and surface-mass-balance ratio δ = a /a.This system of three equations has six unknowns from which four remain when we take α and β as given by the applied geometric transformation.Prescribing one of the three parameter ratios γ , δ or ζ hence determines the scaling of the other two parameters and the time scaling of the system.We can link the ratios of surface mass balance and ice softness by combining Eqs. ( 12) and ( 13), yielding a relation which is independent of the horizontal scaling factor α.For the case of a scaled ice-sheet geometry that is left unchanged in vertical direction (β = 1) ice softness and accumulation hence scale identically.Using Eqs. ( 11)-( 13) we can further express δ and ζ as functions of both geometric scaling ratios and the basal friction ratio: Inserting Eq. ( 14) into Eq.( 16) we also obtain a condition for the basal-friction ratio as a function of both geometric scaling parameters and the surface-mass-balance ratio: Results of an application of the derived scaling laws in numerical flow-line simulations are given in Sect.3.An alternative approach to derive the above scaling conditions based on boundary-layer theory is given in Appendix A.

Two-dimensional case with one time and one length scale
The two-dimensional SSA (Eq. 1) is derived from the full Stokes equation using a single horizontal length scale L and timescale T , respectively (Greve and Blatter, 2009).Continuing this line of thought, we introduce the dimensionless velocity in y direction, v * y = v y T L , in addition to the dimensionless variables from Sect.2.2 to non-dimensionalize the SSA equations.The dimensionless effective strain rate (Eq.2) then reads For the x component of the SSA (Eq.1a) we hence obtain The same coefficients and result from the y component of the SSA, which is not specified here.The nondimensionalized ITE (Eq.4) reads (20) Comparison between the flow-line and the two-dimensional SSA and ITE shows that we obtained the same number of dimensionless factors that appear at the same place and are identical to each other, i.e., θ = , φ = and ω = .Hence under the assumption of a single horizontal length scale the scaling relations for the two-dimensional SSA are the same as in the flow-line case.

Two-dimensional case with time and length scales for both horizontal directions
Starting again from the two-dimensional SSA (Eq. 1) we now make the less-constraining assumption of two horizontal The Cryosphere, 10, 1753-1769, 2016 Table 1.Parameter values as prescribed in the unscaled reference simulations for the flow-line setup (2-D) and the three-dimensional channel setup (3-D), respectively.For the scaling experiments the bed geometry (b x and b y ) and the parameters a, A and C are multiplied with the scaling ratios from to Table 2.The terms "BC-300" and "BC0" refer to the bed geometries described in Feldmann and Levermann (2015a) and b y,G refers to y component of the bed topography used in Gudmundsson et al. (2012).Table 2. Scaling ratios as used for our numerical simulations.Prescribed scaling ratios are highlighted in bold and the others result from Eqs. ( 11)-( 17), ( 22) and ( 23).Each row corresponds to a scaling experiment that is carried out in flow line (2-D) and in a threedimensional channel setup (3-D).The parameter values prescribed in the simulations are obtained by multiplying b x , b y , C, a and A (see Table 1) with the given ratios α, β, γ , δ and ζ .The analytic values for χ are used to fit the sections of linear grounding-line retreat in Figs. 5 and 6.
length scales L x and L y and accordingly two timescales T x and T y , yielding the dimensionless velocities v * In this case the effective strain rate (Eq.2) does not simplify to a single term as in the previous sections but consists of several mixed terms.The SSA thus expands to a much longer expression which we detail in the Appendix B. Although we obtain a multiple of dimensionless coefficients that need to remain the same for the ice sheet to fulfill similarity under scaling, the resulting scaling laws are identical to the ones derived above (see Appendix B).This implies that our requirement of similarity results in the constraint that the ice sheet can have only one timescale T = T x = T y and one length scale L = L x = L y as opposed to our initial assumption of distinct scales for each horizontal direction.
We investigate ice-sheet scaling also in a threedimensional setup in the next section.

Comparison with simulations
We compare our analytical findings with results from numerical simulations applying the Parallel Ice Sheet Model in conceptual geometric setups.The model is the same as used in (Feldmann and Levermann, 2015a) but here run in SSA-only mode.We define a reference topographic geometry which is prescribed in an unscaled reference experiment (indexed as "ref") along with the parameter values shown in Table 1.The scaling experiments use geometrically scaled versions of the reference bed topography and the physical parameters are modified according to the scaling laws derived in Sect.2.2.
Halving the horizontal and/or vertical length scales of the reference topography we obtain three geometric configurations which are shortened in vertical (α, β) = (1.0,0.5), horizontal (α, β) = (0.5, 1.0) or both directions (α, β) = (0.5, 0.5), respectively.To be able to calculate the other physical parameters a, A, C that apply to the scaling experiments according to the three scaling relations (Eqs.15-17) we need to prescribe one more scaling ratio in addition to α and β.Setting γ = 1 (identical basal friction) and δ = 1 (identical surface mass balance), two sub-sets of simulations are thus generated.The resulting scaling ratios which determine the parameter values are shown in Table 2 for each of the seven experiments.We apply the described procedure using (1) a flow-line setup (one horizontal and one vertical direction, bed topography in black in Fig. 2) and (2) a three-dimensional channel-flow setup (flow-line setup extended into second horizontal direction with valley-shaped bedrock in this direction to form a bed trough, see Fig. 4).Details are given in Appendix C.
The experiments are designed to perturb an ice sheet in equilibrium, triggering a marine ice-sheet instability that unwww.the-cryosphere.net/10/1753/2016/The Cryosphere, 10, 1753-1769, 2016 2 ,γ =1 (continuous) and 2-D ref (dashed).Initial steady-state before perturbation in grey; profile at the end of perturbation in red and final steady-state in blue.Output of the reference simulation is scaled by factor 0.5 in both horizontal and vertical direction to allow for comparison of shapes between the two simulations.2 ,γ =1 (continuous) and 3-D ref (dashed).Initial steady-state before perturbation in grey and final steady-state in blue.Output of the reference simulation is scaled by 0.5 in both horizontal and vertical direction to allow for comparison of shapes between the two simulations.
folds unaffected by the ceased perturbation.The speed of unstable grounding-line retreat and the equilibrium ice-sheet profiles before and after the instability serve as a measure to compare the scaling of the dynamic response and the steadystate geometry, respectively.

Comparing timescales of instability
All of our simulations show a similar pattern of groundingline evolution after perturbation (Figs. 5 and 6): after a phase of little to negligible grounding-line retreat the retreat rate increases (grounding line passes the coastal sill and enters the retrograde slope), reaching its maximum value around the minimum of the bed depression before declining to zero (grounding line stabilizes on inland up-sloping bed).The initial and final grounding-line positions of comparable setups (continuous lines) match or are close to each other.The similarity of ice-sheet shapes between different geometric configurations becomes apparent when laying the modeled steadystate ice-sheet profiles on top of each other and scaling the spatial axes according to α and β (shown exemplarily in Figs. 2 and 3).
The simulations clearly differ in the timescale of the evolution of the marine ice-sheet instability (MISI) which can  3 for a cross section along y = 0).Steady-state grounding-line positions for simulations 3-D α= 1 2 ,β= 1 2 ,γ =1 (continuous) and 3-D ref (dashed).Grey lines mark the position of the coastal sill and the bed depression, respectively.Output of the reference simulation is scaled by 0.5 in horizontal direction to allow for comparison of shapes between the two simulations.
be measured by the grounding-line retreat rate ẋgl = ∂x gl ∂t .
To compare different simulations we introduce a retreat-rate scaling ratio: Dependent on which additional parameter is prescribed to be identical under geometric scaling, we replace the timescale ratio using Eqs.( 11) or ( 13) to obtain scaling laws for the retreat rate as functions of the geometric scaling ratios only: We can thus calculate the retreat-rate ratios for all considered geometric configurations (Table 2).The grounding-line curves of our simulations are approximately linear over the time period during which the grounding line passes the bed depression and its retreat rate is largest.We fit a slope to the linear section of the unscaled simulation (purple slope fitted to black curve in Figs. 5 and 6) to obtain our reference retreat rate.Using the calculated retreat-rate ratios from Table 2 we can predict the grounding-line retreat rates for the scaled setups.Superimposing the linear sections of the scaled experiments with the respective analytically calculated slope (Figs. 5 and 6) gives a good match between numerical results and theory (see Figs. D1 and D2 for scaled versions of the time series).Our simulation ensemble of scaled ice sheets thus exhibits similarity as predicted from theory, regarding transient ice-sheet dynamics and steady-state geometry.
4 Implications for the response times of ice sheets Based on the scaling laws derived in Sect. 2 we explore analytically the implications of a scaling of ice-sheet parameters and geometry for the response-time scaling.Making the assumption that the basal friction parameter stays the same (γ = 1) while allowing different values for the surface mass The Cryosphere, 10, 1753-1769, 2016 balance and the ice softness, we are able to calculate the response-time ratio τ (Eq.11) as a function that only depends on the geometric scaling (α and β) and the friction exponent m: Using this equation in combination with Eqs. ( 12) and ( 13) we obtain contour maps for the ratios τ , ζ and δ in the α-β phase space (Fig. 7a-c for the common choice of an exponent value of m = 1/n with n = 3; e.g., Schoof, 2007a;Greve and Blatter, 2009;Cuffey and Paterson, 2010).Therein the blue and red areas correspond to the regimes of an increasing and decreasing parameter value under geometric scaling, respectively, which are separated by a white curve along which the considered parameter remains the same.

Linking horizontal and vertical scales
To be able to follow physically motivated curves through the phase space we link the horizontal and the vertical scale.In idealized flow-line experiments (Feldmann and Levermann, 2015a) it has been shown that the Vialov ice-sheet profile (Vialov, 1958;Greve and Blatter, 2009), though derived under the assumption of the SIA, can be used to also approximate SSA ice-sheet profiles.Motivated by the Vialov profile, for which the central (maximum) ice-sheet thickness is proportional to the square root of the ice-sheet length, we assume a relation between the ice-thickness scale H and the length scale L of the form With α = L /L and β = H /H it follows that for the postulated ice-sheet proportion the two geometric scaling factors are linked such that and Eq. ( 24) then reads We are interested in finding a critical value of the exponent in Eq. ( 24) which determines a threshold in the α-β phase space between the two regimes of increasing (τ < 1) and decreasing (τ > 1) ice-sheet response time under applied geometric scaling.Assuming horizontal elongation (α > 1), which according to Eq. ( 26) implies also vertical elongation (β > 1, see Fig. 7d), it follows that τ < 1 only if the exponent in Eq. ( 27) is negative.Hence there exists a critical threshold with m ∈ (0, 1] and thus q c ∈ ( 1 2 , 1], above which the scaled, i.e., elongated, system responds faster compared to the unscaled system.This is visualized in Fig. 7a for m = 1/3.The area between the dashed (q = q c = 2/3) and the continuous (q = 1) curves is in the regime of τ < 1 for α > 1. Vice versa, for a shortened ice sheet (α < 1 and hence β < 1) in the area between these two curves holds τ > 1.The same qualitative scaling applies to the ice softness whereas the surface mass balance scales oppositely (Fig. 7b and c).
An exponent of q = 1/2 which represents Vialov proportions constitutes the lower asymptotic limit of the domain of all possible q c (limit m → 0; Eq. 27 requires m > 0 for α to remain finite).Thus a Vialov-shaped ice sheet exhibits a response-time scaling oppositely to the scaling explained above (the dotted Vialov curve in Fig. 7a lies always outside the region between continuous and dashed curve, independently of m).
Assuming Vialov conditions under identical friction, the scaling of the response time (Eq.24), surface mass balance (Eq.15) and ice softness (Eq.16) becomes independent of m, which is visualized in Fig. 8. Evaluating the curves in the left vicinity of α = 1, meaning a small reduction in both vertical and horizontal ice-sheet extent, yields a plausible scaling of the ice-sheet parameters in a warming atmosphere: rising atmospheric temperatures cause an increase in surface mass balance (δ > 1 in Fig. 7c; Frieler et al., 2015) and also lead to a softening of the ice (ζ > 1 in Fig. 7b; Cuffey and Paterson, 2010).The response time then decreases (τ < 1 in Fig. 7a).In this picture a warming-induced ice-sheet retreat would hence shift the ice sheet into the regime of faster response to perturbation, tending to accelerate potential further retreat.

Role of basal friction exponent m
The response-time scaling considered here is a function of the basal friction exponent m (Eq.24) and the visualization of the response-time ratio in the α-β phase space (Fig. 7) accounts for only one value of m.To examine the influence of m on the scaling we cut several hypersurfaces through the phase space, sampling the domain of the exponent.
Fixing the horizontal scale, i.e., going along α = 1, causes vertical elongation (shortening) to always result in a shorter (longer) ice-sheet response time (Fig. 9a).In this case the parameter choice of m only determines the curvature of τ (β).Fixing the vertical scale (β = 1) results in opposite behavior of τ , i.e., horizontal elongation (shortening) always yields a longer (shorter) ice-sheet response time (Fig. 9b).Equal geometric scaling of the two directions (α = β) gives a similar picture as obtained for α = 1 (the magnitude of the negative β-exponent is always larger than the α-exponent), with the difference that here the time scaling becomes independent of the geometric scaling for m = 1 (Fig. 9c).
Requiring the response-time scaling law (Eq.24) to be independent of m yields the relation (with k a real number) and thus τ = α k .In general, a negative (positive) value of k then results in a faster (slower) response when elongating (shortening) the ice sheet horizontally.The  In the range of ±50 km around the minimum depression, the slope of the curve of the unscaled simulation is fitted to obtain a reference retreat rate of 0.54 km/yr (purple slope fitted to black curve).This slope is used to predict the slopes, i.e., retreat rates, for the scaled experiments (other curves overlaid by purple lines with predicted slopes) according to Table 1.case of k = 0 yields the same timescale (τ = 1), independent of the α value (Fig. 9d).The case of k = 1 corresponds to the Vialov case for which the timescale ratio increases linearly when elongating the ice sheet horizontally.

Discussion and conclusions
Carrying out a dimensional analysis of the stress balance in SSA and the equation of mass conservation we derive icesheet scaling conditions for the vertical and horizontal length scales, the response time and the relevant physical parameters which determine ice-sheet behavior .Specifically, we find that the scaling relations derived for the SSA in flow line (Eqs.11-13) also hold for the SSA in two horizontal dimensions under the assumption of a single horizontal timescale and a single length scale.Only the twodimensional SSA accounts for stress components that allow for horizontal shearing and hence the effect of buttressing.
Our analysis also shows that although the full SSA accounts for both horizontal dimensions there can only exist one timescale T and one length scale L, as opposed to one for each dimension (T x , T y and L x , L y ) under the principle of similitude.
To non-dimensionalize the SSA stress balance we introduce scales for ice-sheet length, thickness and response time without assuming typical numerical values for these scales.We thus do not neglect further terms in the SSA stress balance by comparing orders of magnitudes of acting stresses (as done in, e.g., Schoof, 2007a) but consider the general case in which both membrane and basal stresses balance the . Scaling of response time τ , surface mass balance δ and ice softness ζ under the assumption of Vialov-type geometric scaling (β = α 1/2 ) and identical basal friction (γ = 1).The resulting scaling conditions are independent of m and given in the legend (n = 3).
driving stress (Eqs.9 and 19).In other studies the scales for ice-sheet thickness, length and time are used to express the friction parameter C, resulting in a dimensionless SSA stress balance that is characterized by a single scaling parameter (often denoted as , Schoof, 2007a;Tsai et al., 2015).In the present study we consider C to be an independent parameter/scale and thus obtain two scaling parameters θ and φ in the SSA stress balance.The resulting scaling laws hence involve the scaling of the basal roughness explicitly.The same holds for the scaling of the surface mass balance a.
The presented scaling conditions can provide rules in the design of model setups for numerical simulations as well as laboratory experiments to obtain parameter sets that leave the ice-sheet geometry (shape and extent) unchanged.For instance, a doubling of the basal-friction parameter under identical surface mass balance requires the ice softness to be reduced to 1 / 8 (Eqs.14 and 15), or a doubling in surface mass balance under identical basal friction requires a doubling of the ice-softness value (Eq.15).Note that these equations by no means make a statement about the physical dependency between ice softness, basal friction or surface mass balance.Our results, derived under the principle of similitude, provide conditions that have to be fulfilled in order to respect self-similarity of idealized ice sheets.In other words: if two (idealized) glaciers, which are in equilibrium with their environment and the SSA equation, have the same qualitative shape but differ, e.g., in surface mass balance, then they also www.the-cryosphere.net/10/1753/2016/The Cryosphere, 10, 1753-1769, 2016 Figure 9. Response-time scaling for hypersurfaces through the α-β-m phase space according to Eq. ( 24) for (a) α = 1, (b) β = 1, (c) α = β and (d) the constraint that the response-time scaling is independent of m.In each panel the legend gives the scaling law for τ that results from the applied constraint.
need to differ in basal friction and their specific relation is given by the scaling conditions.For the numerical simulations conducted in this study we apply parameter configurations that halve the geometric scale in horizontal and/or vertical direction with respect to the reference.The resulting ice-sheet response times range over 3 orders of magnitude (see Table 2).Irrespective of whether they are in a two-or three-dimensional setup, the modeled ice-sheet dynamics, represented by the rate of unstable grounding-line retreat (Figs. 5 and 6), as well as the geometry, represented by ice-sheet shape and grounding line position in equilibrium (Figs.2-4), exhibit the scaling behavior predicted from the analytical calculations to a good approximation.For the flow-line setup three scaled parameter sets show different qualitative ice-sheet evolution compared to the reference, while still complying with the expected response-time scaling.This difference is attributed to the design of the reference setup, i.e., the closeness of the initial steady-state grounding line to the point of instability (local bed maximum).Very small deviations from this position trigger unperturbed instability or prevent landward-induced instability in the scaled setups (see Appendix C).
In contrast to the flow-line configuration the threedimensional setup inherently accounts for the buttressing effect in the initial steady-state simulation due to the presence of a confined ice shelf (Dupont and Alley, 2005;Gudmundsson et al., 2012).However, the ice shelf is removed in the course of perturbation to prevent scale-dependent influences that would originate from a forcing through sub-shelf melting, surface accumulation or ice softness.Thus the speed of grounding-line retreat (and hence ice-sheet response time) is only indirectly affected by the former buttressing effect.An investigation of the response-time scaling under direct influence of ice-shelf buttressing requires a carefully designed experimental setup that maintains the ice shelf during perturbation (as in Asay-Davis et al., 2015) and accounts for the scaling also in the applied forcing.
To analytically investigate the implications of geometric scaling for the ice-sheet response time we make the simplifying assumption of identical basal friction (γ = 1).Though the response-time scaling law then still depends on the sliding exponent m (Eq.24) the qualitative response-time scaling (shorter or longer response time) turns out to be independent of the choice of m (Fig. 9): vertical ice-sheet elonga-tion (shortening) leads to a faster (slower) ice-sheet response and the opposite holds for the horizontal direction.In other words, thicker or shorter ice sheets tend to respond faster than thinner or longer ones.Equal scaling in horizontal and vertical direction (α = β) causes larger ice sheets to respond faster than smaller ones.
Assuming a relation between the horizontal and vertical scale of the form β = α q with 0 < q ≤ 1, we find a critical m-dependent threshold q c for the exponent (Eq.28) above which larger (smaller) ice sheets always exhibit a shorter (longer) response time.The case of q = 1/2 represents the lower asymptotic limit for all possible q c and corresponds to an ice sheet with Vialov-type proportions for which the central ice thickness is the square root of the horizontal extent.Conceptual flow-line experiments similar to the ones conducted here (Feldmann and Levermann, 2015a) revealed that the Vialov profile, which results under simplified conditions from the SIA of the full Stokes stress balance in flow line, can also reasonably approximate the ice-sheet shape in SSA.It is important to note that the SIA and SSA with their distinct underlying assumptions account for two very different types of flow (dominant vertical shearing vs. plug flow).The above-mentioned comparability of ice-sheet profiles resulting from these two approximations is thus not expected to apply in general but might be a product of the simplified, idealized nature of the setup and the choice of parameter values (e.g., relatively high basal friction).In the same study a comparison between steady-state ice-sheet profiles before and after collapse suggested a scaling of β = α 1/2 .For such an ice sheet the time scaling is identical to the scaling of its length, i.e., elongation (shortening) results in slower (faster) response which is opposite behavior than for the above discussed case of q > q c > 1/2.A thought experiment that is consistent with the scaling behavior derived for this kind of profile reveals that in the course of an ice-sheet retreat that is triggered by atmospheric warming the ice-sheet response would become faster, with self-accelerating effect on further retreat (Fig. 8).Note that all the considerations made above presume self-similarity of ice sheets and are only valid for a fixed basal-friction parameter.
In place of prescribing basal friction, the assumption of identical surface mass balance (δ = 1) or ice softness (ζ = 1) results in a more trivial response-time scaling which either equals the vertical scaling (Eq.13) or depends on the vertical scaling via a power-law relation with exponent −n (Eq.12), respectively.Since n is always positive (Cuffey and Paterson, 2010) in the latter case the response time decreases with increasing vertical extent.There are several other ways to analyze the implications of the scaling conditions derived here on ice-sheet dynamics that are not covered in this study.
Our approach includes several assumptions (shallow stress balance, isothermal ice, choice of sliding law, parameter constraints) and thus simplifies the problem of ice flow.At the same time it allows for the fundamental scaling analysis conducted here which incorporates the relevant physics of fast ice flow and results in scaling conditions that relate important physical parameters of an ice sheet to each other.A similitude analysis based on a less simplified stress balance than the one used here would certainly better account for the complexity of real-world systems, but that is beyond the scope of the current study.All statements on the ice-sheet scaling behavior made here therefore need to be considered in light of the idealized character of the underlying simplified SSA stress balance.
The SSA is of vertically integrated form and thus in particular does not account for variations of ice-sheet velocity within the ice column.The assumption of uniform ice softness further reduces complexity, neglecting the dependency of the ice softness on ice temperature, which typically varies in horizontal and vertical direction.The applied Weertmantype sliding law (Eq.3) is a common choice (Fowler, 1981;Schoof, 2007a;Pattyn et al., 2013) amongst others used to describe the sliding of ice sheets over bedrock (Greve and Blatter, 2009;Cuffey and Paterson, 2010;Tsai et al., 2015).Though we prescribe a uniform basal friction coefficient the resulting basal stress field that enters the SSA can vary spatially and temporally.The sliding law covers diverse types of sliding behavior depending on the sliding exponent m in Eq. ( 24).Except for the plastic limit (m = 0), it relates the scale of basal stress to the scale of velocity, resulting in a scaling law which links the scaling of ice-sheet geometry, friction and response time (Eq.11).
Our analytic exploration of the derived ice-sheet scaling behavior applies several constraints to the parameter space and is thus far from being holistic but is aimed to allow for (simplified) statements on the influence of geometric scaling on response time.The set of scaling conditions presented here shall provide a model which allows for a fundamental comparison of the large-scale scaling of the geometry and relevant parameters that determine ice-sheet dynamics.In particular the response-time scaling conditions might be suitable to analyze speed of the transient response to climatic perturbations of the polar ice sheets that took place in the past or might become relevant for the future.

Introducing the dimensionless velocities
L y the non-dimensionalized form of the effective strain rate (Eq.2) reads Insertion into Eq.(1a) yields the following expression for the x component of the two-dimensional SSA: The Cryosphere, 10, 1753-1769, 2016 with the dimensionless coefficient x which has the same form as (Eq.19) but is specific for the x direction.The terms I , II = IV and III are evaluated in the following to obtain dimensionless factors for the SSA equation.The first expression I reads In order to obtain the same equations independent of an applied ice-sheet scaling, the dimensionless coefficients need to remain the same.We start with the first set of coefficients: We immediately see that Eq. (B14) gives the same time scaling (in x direction) as derived for the more constraint cases, i.e., in flow-line (Eq.12), as well as for the two-dimensional case that assumes only one time and length scale (Sect.2.3).Comparison of Eqs.(B14) and (B17) directly yields the condition α x = α y .Furthermore, replacing β and ζ in Eq. (B15) using Eq.(B14) we obtain τ x = τ y .These two conditions can also be deduced by comparing scaling relations that are derived from different coefficients, e.g., I,1 , I I,2 and I I I,1 .
The same procedure can be carried out for the y component of the SSA, leading to the same outcome due to the symmetry of both horizontal components of the SSA.Applying our findings it follows x = in Eq. (B2) and the dimensionless ITE is identical to Eq. ( 20) with the same coefficient .We thus found that in order to fulfill the required scaling similarity in the considered two-dimensional SSA-case there can only exist one horizontal length scale and one timescale (as opposed to one in each horizontal direction, as assumed initially).All the scaling relations derived for the flow-line SSA case (Eqs.11-13) hold here.
Figure1.Schematic of the similitude-analysis method carried out in this study.A reference system (blue ice sheet and bed topography) with geometric scales h and L, timescale T and physical parameters ice softness A, basal friction coefficient C and surface mass balance a is scaled in horizontal and vertical direction (red contours, primed system).The goal is to derive the scaled parameters of the primed system under which dynamic similarity between both ice sheets holds.A dimensional analysis of the governing equations yields dimensionless scale factors which have to remain the same under scaling to attain similitude.The resulting scaling laws determine the scaled (primed) parameters.

Figure 4 .
Figure 4. Bed topography of the three-dimensional channel setup, here shown in the scaled version with α = β = 0.5 (see Fig. 3 for a cross section along y = 0).Steady-state grounding-line positions for simulations 3-D α= 1

Figure 5 .
Figure5.Time series of grounding-line position for the reference and three geometrically scaled flow-line experiments with (a) identical basal friction and (b) identical surface mass balance, respectively.Grey horizontal lines indicate location of the minimum of the bed depression for both the scaled and unscaled case.In the range of ±50 km around the minimum depression, the slope of the curve of the unscaled simulation is fitted to obtain a reference retreat rate of 0.54 km/yr (purple slope fitted to black curve).This slope is used to predict the slopes, i.e., retreat rates, for the scaled experiments (other curves overlaid by purple lines with predicted slopes) according to Table1.

Figure 6 .
Figure6.Time series of centerline grounding-line position (along y = 0) for the reference and three geometrically scaled 3-D channel experiments with (a) identical basal friction and (b) identical surface mass balance, respectively.The fitting method is the same as described in Fig.5with a reference retreat rate of 1.18 km yr −1 .

Figure D1 .
Figure D1.Scaled time series of grounding-line position for the reference and three geometrically scaled flow-line experiments with (a) identical basal friction and (b) identical surface mass balance, respectively.For better comparison the grounding-line curves are shifted along the time axis to overlap where the grounding line passes the minimum of the bed depression (grey horizontal line).Around this point the scaled curves of unstable grounding-line retreat approximately collapse into a single curve, indicating that the scaled speed of the grounding-line instability is approximately the same throughout the ensemble of scaling experiments.See Fig. 5 for unscaled version of this figure.

Figure D2 .
Figure D2.Scaled time series of centerline grounding-line position (along y = 0) for the reference and three geometrically scaled 3-D channel experiments with (a) identical basal friction and (b) identical surface mass balance, respectively (analogous to Fig. D1).See Fig. 6 for unscaled version of this figure.