A Generalized Interpolation Material Point Method for Shallow Ice Shelves. 2: Anisotropic Nonlocal Damage Mechanics and Rift Propagation

Abstract Ice shelf fracture is responsible for roughly half of Antarctic ice mass loss in the form of calving and can weaken buttressing of upstream ice flow. Large uncertainties associated with the ice sheet response to climate variations are due to a poor understanding of these fracture processes and how to model them. Here, we address these problems by implementing an anisotropic, nonlocal integral formulation of creep damage within a large‐scale shallow‐shelf ice flow model. This model can be used to study the full evolution of fracture from initiation of crevassing to rifting that eventually causes tabular calving. While previous ice shelf fracture models have largely relied on simple expressions to estimate crevasse depths, our model parameterizes fracture as a progressive damage evolution process in three‐dimensions (3‐D). We also implement an efficient numerical framework based on the material point method, which avoids advection errors. Using an idealized marine ice sheet, we test the creep damage model and a crevasse‐depth based damage model, including a modified version of the latter that accounts for damage evolution due to necking and mass balance. We demonstrate that the creep damage model is best suited for capturing weakening and rifting over shorter (monthly/yearly) timescales, and that anisotropic damage reproduces typically observed fracture patterns better than isotropic damage. Because necking and mass balance can significantly influence damage on longer (decadal) timescales, we discuss the potential for a combined approach between models to best represent mechanical weakening and tabular calving within long‐term simulations.


Introduction
In Section S.1 of this supporting information, the early MISMIP+ creep damage accumulation for isotropic ( = 0) and mixed isotropic/anisotropic creep damage ( = 0.5) are reported at similar levels of rift propagation as given for fully anisotropic creep damage ( = 1) in Figure 2 of the main text. In addition, the damage field at calving is given for the fully anisotropic creep damage case where the entire Cauchy stress tensor is scaled by damage, rather than only scaling the deviatoric components. Further description and implementation details of the SSA zero-stress damage model (Sun et al., 2017) and the necking and mass balance modification (Bassis & Ma, 2015) are given in Sections S.2 and S.3, respectively

S.2 Description of zero-stress damage model
In the zero-stress criterion, closely-spaced crevasses are assumed to propagate to the depth at which the net longitudinal maximum principal Cauchy stress is zero (Nye, 1957). The net Cauchy stress at depth is parameterized as where ()) ( ) takes the same form as within the creep damage model from Equations (15)- (17).
We disregard the water pressure term for surface crevasses and assume dry surface conditions. A zero-stress isotropic damage variable was previously defined for SSA models as the ratio of the We restrict our zero-stress damage tests to the fully-isotropic and fully-anisotropic cases.
For full-anisotropy, the initial damage accumulation for the zero stress model occurs on a single plane aligned normal to the maximum principal stress of the undamaged configuration, as in the creep damage model. This plane subsequently rotates over time according to spin, as in Equation (8). However, unlike creep damage, anisotropic zero-stress damage accumulation must be restricted to this plane at subsequent time steps, and evolves according to the stresses normal to the plane because the zero-stress criterion assumes crevasses open in accordance with tensile (Mode I) fracture. Rifting is incorporated with the same 2-D critical damage rupture scheme from the creep damage model. For stability, we attempt to restrict the change in zero-stress damage between time steps using the same adaptive time-stepping scheme used for the creep damage simulations, but defining dD 7777 !"# = max; 4 ,-$ − 4 , <. Note that for zero-stress damage, we must eliminate the condition to restart the damage solution if dD 7777 !"# > 0.075 because damage is solved implicitly; however, we still achieve stable simulations with sharp rifting because the next time step is always adjusted with the goal of achieving dD 7777 !"# = . , where we set . = 0.05.

S.3 Description of damage modification for necking and mass balance
Necking describes the process in which basal crevasses widen under tension and the resulting feedback on crevasse evolution, where depending on strain-rates and crevasse-geometry, the ratio of crevasse penetration to ice thickness (i.e. damage) will either increase or decrease over time (Bassis & Ma, 2015). The ratio can increase due to greater thinning rates associated with the presence of crevasses. However, as crevasses grow, the local ice geometry simultaneously adjusts to hydrostatic equilibrium, and depressions fill with surrounding ice due to "gravitational restoring" stresses. If the system is dominated by these gravitational stresses rather than thinning, the ratio of crevasse penetration to ice thickness will decrease (i.e. healing). The ratio is further modulated by mass balance processes, such as melting and accumulation of snow or marine ice in crevasses. A previous study investigated this complex coupling of various processes, and an expression for large-scale ice flow was proposed using perturbation analysis that defines the rate at which damage is modulated according to necking and mass balance processes (Bassis & Ma, 2015). This model can be employed in conjunction with a mechanical damage model that tracks crevasse depths, but has not yet been tested to our knowledge.
When used in conjunction with the zero-stress model, this large-scale damage modification takes the form: where the first term in the parentheses describes the influence of necking on damage and the second term describes the influence of the melt rate, ̇. Within the necking term, parameter * is an effective flow law exponent and 0 describes the ratio of gravitational restoring stress to tensile stress. Derivation of these terms is non-trivial, and we direct the reader to the original publication for a detailed explanation. The expression is only valid in the long wavelength limit, which corresponds to the following assumptions: crevasses are wide compared to the ice thickness, perturbations are assumed to relax immediately to hydrostatic equilibrium, and the melt rate in crevasses is equivalent to the large scale melt rate. We solve (S2) immediately after completion of the SSA solution, and add the damage increment to the zero-stress damage calculated during the SSA.