As a conveyor belt transferring inland ice to ocean, ice shelves shed mass through large, systematic tabular calving, which also plays a major role in the fluctuation of the buttressing forces. Tabular iceberg calving involves two stages: first is systematic cracking, which develops after the forward-slanting front reaches a limiting extension length determined by gravity–buoyancy imbalance; second is fatigue separation. The latter has greater variability, producing calving irregularity. Whereas ice flow vertical shear determines the timing of the systematic cracking, wave actions are decisive for ensuing viscoplastic fatigue. Because the frontal section has its own resonance frequency, it reverberates only to waves of similar frequency. With a flow-dependent, nonlocal attrition scheme, the present ice model [Scalable Extensible Geoflow Model for Environmental Research-Ice flow submodel (SEGMENT-Ice)] describes an entire ice-shelf life cycle. It is found that most East Antarctic ice shelves have higher resonance frequencies, and the fatigue of viscoplastic ice is significantly enhanced by shoaling waves from both storm surges and infragravity waves (~5 × 10−3 Hz). The two largest embayed ice shelves have resonance frequencies within the range of tsunami waves. When approaching critical extension lengths, perturbations from about four consecutive tsunami events can cause complete separation of tabular icebergs from shelves. For shelves with resonance frequencies matching storm surge waves, future reduction of sea ice may impose much larger deflections from shoaling, storm-generated ocean waves. Although the Ross Ice Shelf (RIS) total mass varies little in the twenty-first century, the mass turnover quickens and the ice conveyor belt is ~40% more efficient by the late twenty-first century, reaching 70 km3 yr−1. The mass distribution shifts oceanward, favoring future tabular calving.
Ice shelves are thermal and dynamic protectors of inland sheet ice. The predominance of ice shelves around the Antarctic Ice Sheet (AIS; Figure 1) is a manifestation of the delicate balance between inland ice stream discharges and atmospheric/oceanic interactions. The retreat of the AIS’s calving fronts, under warming climatic conditions, may cause adjustments of the flow speeds of inland ice streams, and force the ice geometry to reach a new balance (e.g. Mercer 1978; MacAyeal et al. 2003; Rott et al. 2002; Rignot et al. 2008; Scambos et al. 2004; Cook and Vaughan 2010; Joughin and Alley 2011). In certain ice–bedrock–ocean configurations (e.g., the inland-slanting bathymetry of the West Antarctic Ice Sheet), a positive feedback forms and a new balance is achievable only when the bulk of the ice is consumed (Hughes 2011; Joughin and Alley 2011; Schoof 2007; Hindmarsh 1998).
Ice-shelf dynamics is a major challenge in ice-sheet modeling (Bindschadler 2006; Martin et al. 2011) for a physical process-based assessment of the sea level rise (SLR) over the twenty-first century. Despite recent increased understanding of ice-sheet physics (Joughin and Alley 2011), the dynamic processes of calving, especially the large, systematic tabular ice attrition, still is not fully understood (Benn et al. 2007; T. Hughes 2011, personal communication). There remains a gap before a full projection capability is possible of future ice mass loss through systematic calving. There are several prevailing hypotheses concerning calving. Glaciological stresses are the major control on rift opening (Joughin and MacAyeal 2005; Alley et al. 2008; Larour et al. 2005). Scambos et al. (Scambos et al. 2005) proposed that incision from ponded water is responsible for calving of shelves residing in relatively warm places. Similarly, many researchers (Shepherd et al. 2003; Glasser and Scambos 2008; Vieli et al. 2007; Rignot and Jacobs 2002) also identified warm ocean water incision as a contributor. In situations of significantly fatigued cracks, periodic forcing at the calving front by storm-caused infragravity waves (Bromirski et al. 2010) serves as a trigger for calving of peripheral ice shelves. Modeling of rift generation and development has evolved from the elastic (Reeh 1968) to viscoplastic mechanical modeling (MacAyeal et al. 2003; Pralong and Funk 2005; Albrecht et al. 2011), focusing on different aspects of the ice fatigue process. Unfortunately, none includes sufficient mechanisms to simulate the entire life cycle of calving. Moreover, no existing scheme has explicit dependence on ice flow, although strain rate arising from spatial variations in velocity, which determines the location and depth of surface crevasses, is believed to be the “first order” control on calving (Benn et al. 2007).
None of the proposed calving mechanisms answers the question of “when calving occurs” (D. MacAyeal 2012, personal communication). Both the elastic calving mechanism (Reeh 1968) and the growth of ductile cracks (Kenneally and Hughes 2006) examine calving as local processes near the calving front. Tabular calving is a nonlocal process involving ice flow vertical shear within the shelf. This study proposes a calving scheme that reconciles existing schemes and explains the multifrequency nature of tabular calving, verified by observed calving-frequency data. Critical for this research is a recently available 1-km-resolution ice thickness and elevation dataset for the peripheral shelves (Griggs and Bamber 2011). The proposed scheme is primarily for large, systematic tabular calving. It does not indicate that smaller-scale calving does not occur in between. Essentially, the occurrence of smaller-scale tip-end calving is inevitable because of random collisions with floating ice driven by ocean currents and large-amplitude waves. These harbingers postpone and may even completely inhibit the occurrence of the systematic calving.
In section 2, there is a brief description of the ice model and the calving formula. Verification of the calving scheme is presented in section 3. In section 4, a projection is made of the future systematic calving of Ross Ice Shelf (RIS), the largest of the Antarctic ice shelves. The role of ocean waves in tabular calving is estimated objectively.
2. Ice model and input data
The Scalable Extensible Geoflow Model for Environmental Research-Ice flow submodel (SEGMENT-Ice) is an extensively evaluated ice dynamics model (Ren et al. 2011a; Ren et al. 2011b; Ren and Leslie 2011; Ren et al. 2013a). It handles various ice–bedrock–ocean configurations existing over Antarctica. By solving the 3D conservation equations for mass, momentum, and energy under multiple rheological relationships, SEGMENT-Ice provides prognostic fields of the driving and resistive forces and describes the flow fields and dynamic evolution of ice thickness profiles.
Details of the model numerics and physics are provided in Ren et al. (Ren et al. 2013a). In particular, SEGMENT-Ice is projected to a spherical coordinate system. Because the simulation domain contains the South Pole, to avoid the singularity, a rotation is made to shift the pole out of the calculation space (Ren et al. 2013b). The vertical sigma coordinate is stretched to better represent the near surface and the bottom ice layers. Vertical stretching also allows an explicit representation of the air–ice–water interface for ice shelves. The freeboard for most ice shelves is well less than 100 m. If the ice layer is uniformly divided, as in most peer models, it is difficult to set the lateral boundary conditions according to different configurations of ice, water, and air. SEGMENT-Ice also has a four-dimensional variational data assimilation (4DVAR) ability to better use the remote sensing data (e.g., surface ice velocity) and to better parameterize the basal sliding conditions (Ren et al. 2013b). The model supports multiple rheologies (of different materials). If a grid point is next to air, free-stress lateral boundary conditions apply. If it is next to water (ocean), no shearing stress is from its water neighbor, because the viscosity of water is many orders of magnitudes smaller than solid ice. The imbalance in normal (hydrostatic) pressure (Figure 2b), however, is the driving mechanism for shelf advancing. The following sections focus on AIS applications and provide a context for the numerical experiments. Only the calving scheme is detailed here and interested readers are referred to Ren et al. (Ren et al. 2013b) for details of the general model configuration.
2.1. Calving is a fatigue process
The AIS ice flow field consists a slow-moving sheet flow (a viscoplastic flow < 5 m yr−1), channelized streamflow (a combination of sliding and the deformation of underlying basal sediments; 100–500 m yr−1), and the free creeping of the fringe ice shelves (reaching 1000 m yr−1 at the frontal edge; Figure 2a). Many studies have found that transitions among flow regimes are highly correlated with the strength of ice-bed coupling. For example, the “grade glacier” theory (Alley et al. 2003) generalizes silt production and transportation as an integrated component of the ice erosion on glacier beds. Hindmarsh (Hindmarsh 1998) analyzed the stability of a thin till layer overlaid by a deep ice layer, using a perturbation method based on realistic assumptions on viscosities for ice (as viscous fluid) and till (power-law dependence on stress and effective pressure). His model provides an elegant explanation of the ice–bed coupling as an important mechanism for deforming-bed features of stream ice flow (e.g., formation of drumlin under ice and causing intermittent acceleration and deceleration of ice streams over the Siple Coast). Hughes (Hughes 2011) found that the transitions from sheet flow, through streamflow, to shelf flow are based on a progressive reduction in the strength of ice-bed coupling. This indicates the importance of pressure meltwater and basal granular material. SEGMENT-Ice has an ice–granular material–bedrock configuration for regions with granular material present. As ice–ocean contact provides an important path for climate warming to influence ice loss, Ren and Leslie (Ren and Leslie 2011) proposed a grounding-line dynamics scheme.
Ice-shelf calving is essentially a fatigue process of viscoplastic ice. Antarctic snow falls under a wide range of temperature, humidity, and winds. Ensuing conversion to firn further complicates the polycrystalline structure so that natural glacier ice has strong heterogeneous rheological properties (Budd and Jacka 1989). SEGMENT-Ice enhances ice viscosity according to the age of ice and topographical slope. Brittle ice has a very short elastic range (up to 0.3% deformation). Immediately after, there is primary creep and nonrecoverable deformation. Material breaks when stress in one principal direction is not properly balanced in the other two principal directions (the von Mises yielding criteria),
where is principal deviatoric stress (full stress adjusted by the hydrostatic stress) and is fracture strength. For ice, it is a function of ice strain history and temperature, which is ~2 MPa for the RIS.
SEGMENT-Ice systematically simulates the evolution of ice geometry, temperature, and ice flow. The ice rheological relationship is applied so that the deviatoric stress fields are estimated as strain rate multiplying the ice viscosity. Full stress is the sum of deviatoric stress and the hydrostatic stress tensor. Because the model equation is not constructed on principal stress coordinate, a linear transformation is needed when checking the yielding criteria [Equation (1)]. The underlying physics is contained in the appendix, where a cantilever beam assumption illustrates the mechanics of the limiting length of ice-shelf extension. The nonlocal calving scheme in SEGMENT-Ice (Figure 2b) emphasizes the role of ice geometry (e.g., the forward-slanting calving front profile; Figure A2b) in tabular calving of the embayed Antarctic ice shelves (MacAyeal et al. 1996) and reconciles existing ice-shelf calving hypotheses.
Because parts of the calving front have different ice thicknesses, temperatures, and flow speeds, it is unlikely that a through cut (e.g., a “big” crack through the entire depth of shelf and all across the ice front arch) is completed at once. Instead, one sector breaks first, and then the crack propagates laterally to neighboring areas (Figures 2 and 3). From Griffith’s energy criterion (to be detailed soon), further cracks develop if the elastic energy released from further cracking exceeds the required surface energy, set here at 30 mJ m−2, to propagate the crack.
2.2. General cracking procedure of brittle matter (Griffith’s law)
It is well known that attempts to calculate the tensile strengths of simple crystals give results very much higher than those observed experimentally. Griffith hypothesized this is because of the presence of a large number of minute cracks, known as Griffith cracks, inside the material. Here, it is not intended to repeat the sophisticated derivations, which include algebraic skills. The essence of the Griffith theory is that particles are confined in an energy trap within a certain range of the equilibrium position. Depending on the bonding, this can be determined by magnetic/electronic forces and/or van der Vaals–London forces. Very similar to the free-length theory in turbulent fluid dynamics, once a particle can move away from the equilibrium position farther than the interparticle average spacing, fluidity appears and macroscopically manifests itself as a fracture or yield. Thus, there are two rival players in the further development of an existing crack: one is the release of potential energy [W; as in Equation (A1)] but the other is the increase in surficial energy (creating new interfaces). Assume, on a flat plate, there is an elliptical hole of major axis 2l and subject to an average tensile stress σ (Figure 3). The concentration of stress at the acute tip is
where r is the inverse curvature at the ends of the major axis. Taking r as the intermolecular spacing a, the crack will spread if reaches , the maximum tensile stress that can be sustained by the material without cracking. To estimate the value, the process of cracking must be considered: this produces two new surfaces within the material with a distance apart of the order of the intermolecular spacing a and with each possessing surface energy J per unit area that is an intrinsic and measurable property of the material. The released strain energy must exceed this surface energy for cracking to proceed. Using the same concept as in Equation (A1) and applying a linear strain–stress relationship, the amount of strain (mechanical) energy is . The critical requirement for cracking to proceed is
where D being the uniform depth of the plate (or ice-shelf thickness). That is,
Equating in Equation (4) with the expression in Equation (2), the uniaxial tension in the direction perpendicular to the crack is obtained as , which is on the same order as results from experiments with ceramics and glasses. This expression also explains the accelerating characteristic of crack propagation before being arrested by another crack.
From Griffith’s energy criteria, further cracks develop if the elastic energy released from further cracking exceeds the required surface energy (e.g., 30 mJ m−2 for solid ice) to propagate the crack. Application of this criterion to the problem at hand with fractured ice shelves is complicated by it being a multiple-body problem. It must be assumed that an ice shelf contains a large amount of randomly oriented incipient cracks and that the cracks proceed in the direction of the highest local stress (tip ends of longest cracks). The most dangerous configuration is that of a crack heading to and merging with another “parallel” crack. At the other extreme, resting of a crack happens when two perpendicular cracks meet and release the concentrated stress, as depicted in Figure 3.
For the case of ice-shelf calving, the energy input for propagating the existing crack can be from many sources. However, if the force exerted to the shelf always goes out of phase with its motion direction, negative work is done and this does not contribute its distortional energy (W) buildup. For this reason, periodic forcing that can build up coherence with the motion of shelf itself transfers energy more efficiently than the random form of forcing. For periodic wave forcing, wavelength, amplitude, and the frequency difference of wave and resonance frequency of ice shelf all are important factors affecting the energy transfer/coupling between shelf and waves. For ocean waves, wave speed is determined (to a large degree) by water depth (ocean bathymetry); thus, wavelength and frequency are related.
Next, examine the lateral propagation of an existing crack on an ice shelf. Suppose the crevasse depth is H and the crack surface is a plane surface without any curvature. Crevasse depth is a function of longitudinal strain rate. In SEGMENT-Ice, it is a function of vertical flow shear and ice thickness (or it is a function of ice geometry and flow structure). Actually, the vertical penetrating of the crack ceases shortly after the frontal section is in hydrostatic with the ocean water. According to Griffith’s theory, the increment in the newly exposed surface area is related to the input of energy by
where is the length increment of the crack at ice surface and is energy input for each periodic forcing. To estimate the total number of periodic forcing that is needed before the length of crack is accumulated to finish a complete calving, the above expression can also be expressed as
In practice, as the length of the crack accumulates, the number of periodic forcing is recorded simultaneously.
2.3. One shelf’s noise is another shelf’s music: The resonance frequency of ice shelves
Because cracks that do not cyclically open and close do not grow rapidly, tides and other randomly oscillating factors play roles in ensuing crack development. An example is the wintertime refreezing of percolated surface melting water in the crevasses of ice shelves. Expansion associated with refreezing enlarges the crevasses and is effective in calving of lower-latitude, warmer ice shelves (J. Jacka 2012, personal communication) and explains the domino-like calving scenario reported in MacAyeal et al. (MacAyeal et al. 2003).
Tidal amplitudes around Antarctica generally are small compared with shelf ice thickness. Only tides with resonance frequencies of the ice shelves can effectively transfer energy into shelves and influence crack-opening rates. To evaluate the role of waves in the calving process, the resonance frequencies of peripheral ice shelves are estimated by solving the deflecting equation using a Newton–Raphson iteration scheme, with realistic ice-shelf configurations and distributed ice density and temperature. The ice shelf is partially immersed, so hydrostatic-pressure-caused frequency decreases are taken into consideration,
where is the resonance frequency without tip-end pressing, p is hydrostatic pressure, and is buckling load. From Equation (6), taking into account that the discord of forcing with the resonance frequency of the ice shelf, the crack tearing rate is expressed as
where a is the crack length (in lateral direction parallel to the calving front) starting from 0, N is the fatigue cycle, and is the Paris coefficient (Timoshenko and Gere 1963). Apparently, J and H in Equation (6) are ascribed into coefficient . The is a negative number indicating the exponential damping of the tidal tearing when it is not synchronized with the natural frequency of the ice shelf, is the tidal frequency, is resonance frequency of the ice shelf [i.e., in Equation (7)], and is the range of stress intensity change (proportional to the tide amplitude squared). Ice fatigue is fundamentally different from metal fatigue as the hydrogen bonding is weaker and the energy trap is more easily overcome, so phase changes occur and repair the break (damage). Thus, for regions with seasonal melting, a [in Equation (8)] is modified annually, to represent the refreezing healing contribution.
As a possible cause of ice-shelf disintegration and fracture, mechanical interaction of ocean waves with ice shelves was explored first by Holdsworth and Glynn (Holdsworth and Glynn 1978) and has more recently been reevaluated by a number of researchers (Sergienko 2010; Sergienko 2013; Bassis et al. 2008; Bromirski et al. 2010). Sergienko (Sergienko 2010) established a theory of normal modes based on an elastic ice-shelf assumption and explained how the flexural-gravity waves depend on the coupling between the elastic ice shelf and the underlying ocean flow. Practical applications remain uncommon. Quantifying the wave effects on ice-shelf fatigue is possible from Equation (8). From Equation (8), tidal fluctuations on time scales much shorter than needed for “ice to adjust rheologically” (Vaughan 1995) have limited effects on calving. Equation (8) also explains why sea ice dampens the storm-surge-caused waves most significantly. If there are no waves, there is no apparent limitation of the horizontal dimension of sea ice. Wave energy breaks the ice into sizes that have the resonance frequency of the driving wave. Thus, the average size of sea ice is controlled by prevailing wave activity. Apparently, the ice pieces, as a group, are effective in damping storm surges and dissipate the wave energy into heat. Thus, sea ice is a critical buffer for sea waves to fatigue ice shelves.
The fatigue process details affect not only the calving period but also the shape of the newly formed calving front. Sometimes it can form a backward-slanting interface with the ocean (Kenneally and Hughes 2006), an unstable configuration that will eventually restore forward-slanting configuration. Mechanical energy from waves is dissipated through turbulent eddies. Waves of small amplitude are unimportant in progressing the cracking of an ice shelf because 1) the energy it transfers to shelf is proportional to the square of the amplitude, so smaller wave transfer energy at a much lower rate—if the transferring rate is close to the dissipation rate, there is little damage to the ice shelf—and 2) the residual elastic range may not be overcome by small waves and contribute little to the fatigue.
The calving scheme is schematically presented in Figure 2b. Ice is brittle and wherever there are strong concentrations of strain rate, it breaks (C1). The top panels of Figure 2b are schematic diagrams of ice profiles, showing the different flow regimes. Ice shelves thin by creep thinning. The negative vertical strain rate (compression) causes a horizontal divergent (positive) strain rate. In the diagram, the white bulk arrows are stress (hydrostatic pressure) on the right side of the calving front exerted by the ocean, decreasing to zero at sea level. The red bold arrows are static stresses exerted on the left side of the front, decreasing linearly to zero at ice upper/subaerial surface. The red curve is the net horizontal stress, which reaches maximum at sea level. The vertical profile of the horizontal ice velocity field determines that there will be a “mushroom”-shaped outspread section that is out of hydrostatic balance with the ocean water. There is a limit to the length of this section before it breaks off from the main body of the shelf (bk in the diagram). The limiting length is regularly approached to finish a systematic calving [Equations (A7) and (A8)]. There are random components in the calving processes, such as hydrofracturing (Doake and Vaughan 1991). In the bottom row of Figure 2b, the two stages of tabular calving are illustrated. In the first stage (), lateral crevasses can arrest horizontal extension of flow-transverse crevasses (those parallel to the front). Only the latter part of the first stage sees rapid crack vertical extension. The second stage dominates the calving irregularity and, in this stage, the front-parallel crevasses merge and extend horizontally up to hundreds of kilometers. The end of this stage signifies the birth of giant tabular icebergs. Resistive stress reduces significantly and shelf ice flows faster, beginning another calving cycle. As a concluding remark, in the proposed two-stage calving scheme tsunamis can cause calving but only if the shelf has already been preconditioned by ice flow. Put in the Benn et al. (Benn et al. 2007) convention, the von Mise criteria–governed, ice flow–determined first stage of calving (illustrated in the appendix using a cantilever beam model) is the first-order calving mechanism. The wave activity, like other cyclic forcings, is of secondary importance.
2.4. High-quality input data from SeaRISE project
The particular cryospheric data fields of AIS used for the tabular calving experiments are quality-controlled data from the Sea Level Response to Ice Sheet Evolution (SeaRISE) project (http://websrv.cs.umt.edu/isis/index.php/data). The adaptive design of SEGMENT-Ice makes it possible simultaneously to use 5-km-resolution ice geometrical data for the slow creeping inland portion and 1-km-resolution data for the shelf regions. The unique frontal geometry (see appendix) and oceanic erosion both can vary on a 10-km scale. Thus, the 1-km ice-shelf data (Griggs and Bamber 2011) are critical for the success of this study. The ocean bathymetry is originally from Nitsche et al. (Nitsche et al. 2007). The basal (ground) heat flux is from Fox Maule et al. (Fox Maule et al. 2005). For modeled ice velocity verifications, we used Making Earth System Data Records for Use in Research Environments (MEaSUREs; Joughin and Alley 2011) and Modified Antarctic Mapping Mission (MAMM; Jezek 2003) measurements. The National Ice Center (NIC; http://www.natice.noaa.gov/products/antarctic_icebergs.html) records are used to verify the modeled statistics of calving events. For East Antarctic coasts, acoustic records of about 8 years (2004–12) also are used.
In addition, SEGMENT-Ice needs atmospheric and oceanic forcing as input. For the verification period, we used the National Centers for Environmental Prediction–National Center for Atmospheric Research (NCEP–NCAR) reanalyses (Kalnay et al. 1996), employed widely by climate researchers. The reanalyses were merged with atmospheric forcing provided by the Geophysical Fluid Dynamics Laboratory (GFDL) High Resolution Atmospheric Model, version 2.1 (HIRAM2.1; Zhao and Held 2010; Putman and Lin 2007) for 1948–2009. For the projection of tabular calving timing, location, and iceberg size, atmospheric forcings are from HIRAM2.1 model run. The HIRAM simulations, specifically the C180 HIRAM2.1 simulations, are chosen over other Program for Climate Model Diagnosis and Intercomparison (PCMDI) climate model simulations for the Fourth Assessment Report (AR4) primarily because of its finer horizontal (~50 km) and vertical grid spacing (32 vertical levels) in the version described by Zhao and Held (Zhao and Held 2010). For the twenty-first-century projection, only the Special Report on Emissions Scenarios (SRES) A1B emission scenario runs are used. Oceanic data are used in the chemical-potential-based grounding-line scheme. Estimates of the future ice-shelf calving are based on the International Panel on Climate Change (IPCC) SRES A1B scenario.
Unique to the nonlocal tabular calving scheme is that it depends on 3D ice flow (see appendix). Accurate, full 3D velocities are required for estimating systematic tabular calving. The model-simulated surface ice flow fields compare well with remotely sensed measurements (Ren et al. 2013b) for all of the three flow regimes over the entire AIS. There are few ice velocity observations at depth, but the model-estimated AIS total mass loss over 2003–08 agrees with satellite observations from the Gravity Recovery and Climate Experiment (GRACE; Ivins and James 2005) to within 2%, despite the very different horizontal resolution, as both are about −180 km3 yr−1. As the dynamic component of mass balance is a major term in the total AIS mass balance, this agreement between SEGMENT-Ice and GRACE is encouraging. Ice flow fields are directly related to the dynamic mass balance component (i.e., mass changes due to flow convergence and divergence) of the total mass balance. For the total mass balance to fit the observations for the correct physics, the thermal regime within the ice also needs to be accurately simulated.
3.1. Thermal regime inside the Ross Ice Shelf
Ice temperature affects every aspect of the total mass balance and thus the ice thickness H. Ice flows faster as the temperature increases, because of reduced ice viscosity; dynamic mass balance is thus sensitive to ice temperature. The bottom of an ice shelf is at pressure melt with the ambient salty ocean water. Whether the ocean water freezes to the ice shelf or the shelf ice melts into ocean depends on the net energy gain of the bottom ice layer. In this sense, ice temperature also is critical for the net surface mass balance. It is therefore important to determine which processes dominate ice temperature evolution.
Local rate of change of temperature depends on advection, conduction, and heat sources. In SEGMENT-Ice, by switching on and off certain terms in the thermal equation, the relative importance of the three factors can be compared. Flow shear and hence strain rate within the shelf are small and the internal friction as a heating source is minimal. Over the RIS, the surface ice temperature varies by over 30°C. The latitudinal distribution is affected strongly by the confining embayment. Temperature gradients are large along the west coast (the toes of the Transantarctic Mountains). Because of the geographical distribution of the ice temperature, ice thickness, and ice flow patterns, it is impossible to make general statements about the relative importance of the three mechanisms without detailed numerical simulations. The relative contributions to the general heat balance are examined for the five locations selected (the markers in Figure 4, covering different ice flow regimes over RIS). In the insets of Figure 4, the advection heating rates are plotted for the five selected locations. Because temperature decreases upward within ice shelves, flow convergence regions experience vertical advective heating. Because the horizontal poleward temperature gradients, southerly horizontal flow means cooling effects. The net advective effects contain both vertical and horizontal advective contributions and sometimes they cancel each other (among the five, three grids experience advective heating and two experience advective cooling). The sections close to the calving front of the RIS experience strong advective cooling because southerly ice flow is divergent. Consequently, the vertical advection works in synergy with the horizontal advection in causing local cooling. In the central stagnant region (e.g., 80.13°S, 179.5°E), the vertical advection (heating) cancels horizontal advection (cooling by southerly flow). The net effect at the selected location is warming but of lower magnitude (1.1 × 10−4 K yr−1) compared with a more southerly grid downstream of an ice rise (81.9°S, 175.77°W), where horizontal advection effects are minimal (velocity υ is only ~61.8 m yr−1). Over the RIS, magnitudes of advective heating/cooling rates are generally less than 2 × 10−4 K yr−1. The conduction-caused cooling rate, on the other hand, is about 0.02–0.15 K yr−1, which is more than two orders of magnitudes larger than the advective cooling rate. Ice-shelf temperatures primarily are controlled by conductive cooling in the vertical direction and secondarily by the advection of cold upstream temperature but only minimally by the strain heating.
The color shading in Figure 4 is the model-simulated twentieth-century climatological oceanic erosion rates (m yr−1). There are almost no latitudinal features in the geographical pattern. Rather, there are strong correlations with ice thickness profile. This is because the longitudinal gradients of oceanic temperature are smaller than the decrease in the pressure melting point with increasing ice thickness. Thus, beneath the thicker ice, the differences between ocean temperature and melting point temperature are larger. At the same latitude, thinner ice regions have a smaller melt rate also because the strong vertical temperature gradients help remove heat from the bottom ice layer. All else being equal, the oceanic melt is thus smaller. Strong oceanic melt centers coincide with regions of strong temperature advection, because these locations also have larger ice thickness or have a tendency for thickness increases.
Ice shelves are integral components of the entire ice sheet, and its creeping velocity is not locally determined and many other features are nonlocal (MacAyeal et al. 1996). Although we presented only the RIS, the actual numerical experiments are carried out over the entire AIS and neighboring oceans. For the grid points close to borehole measurement sites [e.g., Mill Island close to Shackleton Ice Shelf (Roberts et al. 2013); Bruce Plateau, Antarctic Peninsula (Zagorodnov et al. 2012)], the simulated ice temperature profiles agree with borehole measurements. The differences between modeled and observations, especially beneath the active layer, are much smaller than instrument error.
Accurate thermal regime simulation also warrants a satisfactory simulation of the ice flow field over the entire Antarctic Ice Sheet (Ren et al. 2013b). Zoomed in on the RIS region, there is a stagnant area on the Ross Ice Shelf (<100 m yr−1) of thicker ice at the central part off the steep Transantarctic Mountain coast. Outside this stagnant region, the ice speeds accelerate toward the ice/ocean front and can reach 2000 m yr−1 at the calving front. The discrepancies are generally less than 50 m yr−1 near the calving front, which is insignificant compared to the absolute magnitude of the flow (Figure 4).
3.2. Calving frequency, average tabular iceberg width, and resonance frequencies in the present climate
Calving is a low-frequency phenomenon. The historical record of calving is very short compared with the low frequency of such events (S. Jacobs 2011, personal communication). Direct observations of tabular icebergs are available only in the remote sensing era. The NIC provides the locations and times of all major Antarctic icebergs. The data were sorted to obtain birthplaces of the icebergs of size > 40 km2. The coast of Antarctica is divided into 20° longitudinal bins. Calving events in each bin are counted for 1996–2009 and occurrence frequencies are the number of calving events divided by the number of years (14). The model-simulated calving frequency is compared with observations (Figure 5). As the observational period is <20 years, bins with calving periods all longer than 20 years have no observed calving (e.g., 160°–140°W, 100°–80°W, and 40°–60°E). Thus, observations likely have a 20° neighborhood with at least one ice shelf with a calving period less than 14 years. This is close to the SEGMENT-Ice-estimated calving frequency. Precise calving timing and sizes of the icebergs are not well simulated at present.
The resonance frequencies (Figure 6) of the peripheral ice shelves cover a wide range (5 × 10−5 to 2.5 × 10−6 Hz). Antarctica is geologically quiescent and the primary periodic forcing is from oceanic waves. The influence of sea swells on ice-shelf flexure has long been known (e.g., Holdsworth and Glynn 1978; Williams and Robinson 1979) and recent proposed mechanisms show that ocean gravity waves can potentially affect ice-shelf stability (MacAyeal et al. 2009). Without estimating the ice-shelf resonances, previous conclusions are primarily conceptual and qualitative. From Figure 6, the RIS and the Filchner-Ronne Ice Shelf (F-RIS) both have frequencies of 60–80 × 10−6 Hz, much lower than the shoaling storm surge waves and even the longer infragravity waves. The smaller ice shelves, such as the Shackleton Ice Shelf over East Antarctica (90°E; about 34 000 km2), have resonance frequencies that overlap the storm surge frequency (>240 × 10−6 Hz). Sensitivity of the Shackleton Ice Shelf to wave activities is shown in the acoustically recorded calving events (A. Gavrilov 2012, personal communication). Calving of ice shelves with resonance frequencies overlapping less well with sea waves but with well-fractured surfaces also is observed, coinciding with the arrival of long-period oceanic infragravity waves (e.g., Wilkins Ice Shelf; Bromirski et al. 2010).
Following the above overview peripheral ice-shelf calving, the focus is on the RIS, the largest AIS ice shelf. Using twentieth-century atmospheric/oceanic forcing, model simulations of rift geometries are compared with interferometric synthetic aperture radar (InSAR)-detected time series of rifts on the RIS (Joughin and MacAyeal 2005). When ice flows through Framheim Island (southwest of Whales Bay), there are numerous transverse crevasses. Horizontal extensions of these crevasses are limited and usually are arrested by flow-parallel crevasses before the outspread section reaches length bk, the along- (ice) flow direction dimension of a tabular iceberg (see Figure 2 and the appendix). Then the unbalanced gravitational driving from the outspread section greatly accelerates the rifting process [Equation (6)]. Thus, for the eastern section of the RIS, the rift starts north of the Roosevelt Island and extends westward, eventually calving and forming a tabular iceberg. As the RIS coastline of RIS is convex, a single rift seldom runs through the entire ice shelf. Once the eastern sector is calved, the western part protrudes more and calves easier.
3.3. Future calving of the Ross Ice Shelf
Because the fracture strength of ice depends only weakly on ice temperature, climate warming influences calving events primarily from varying the ice geometry. At present, the RIS’s upper surface net snow precipitation is one order of magnitude smaller than dynamic mass balance and the lower surface oceanic melt (or accretion). The total mass balance is essentially the interplay of the latter two. Close, and beneath, the calving front, intrusions of warm current causes strong oceanic erosion bands (Figure 4). Most likely, ice necking occurs in these bands. Effects of climate warming on ice flow speeds and on surface (basal) mass balance work in synergy in certain locations and cancel each other in other places. Numerical modeling therefore is necessary to determine the net effects. Projections of future calving are made with HIRAM2.1 provided atmospheric forcing parameters. Despite the apparent increase in ice flow speeds, net surface precipitation, and oceanic erosion, the total ice within the RIS did not experience any apparent change (retaining a ~9.38 × 104 km3 level). After adjusting the sawtooth fluctuations associated with tabular calving, the shelf mass experiences a slight increase of ~0.45 km3 yr−1, which resulted from the hysteresis of the oceanic erosion. Although the systematic increase of total ice mass is negligible, the mass distribution becomes more toward the calving front, a geometry favoring tabular calving (Figure 7).
In the A1B scenario of atmospheric forcing, SEGMENT-Ice simulates calving rates and the possible tabular iceberg width (bk) after the year 2050 (Figure 8). In general, large ice shelves (e.g., F-RIS and RIS) and ice shelves at cold locations (e.g., Queen Maud Land, Amery Ice Shelf) calve slowly, at ~20–25 years on average over the RIS, about 5 years shorter than around the year 2000. Thinner ice shelves in relatively warmer environments calve faster but the icebergs are smaller. The model output is an average calving width. For the RIS, the average width is ~18.5 km. Each series of systematic calving discharges ~800 km3 of ice into the ocean. In contrast, for smaller ice shelves fringing the Amundsen Sea coast calving is more frequent but of smaller size.
Cracks that eventually make a tabular calving usually do not start locally. Instead, they are advected from upstream (Figure 2a). The calving always occurs within 70 km of the coastline. Stacking calving locations along the coastline make them hard to discern. A better way to present future calving locations is by tracing back the future calving fronts’ present locations. Numbers in Figure 2a is produced for such purpose. The contour lines are labeled by years. The information should be interpolated as “crevasses around line 50 contribute to tabular calving ~50 years later, when it is advected downstream to the coastline.” From Figure 2a, because of the cold trace left by B15, tabular calving along the RIS in the next 20 years will primarily be in the western sector (Ross Island and Franklin Island side). Calving is imminent along the Franklin Island sector (within 5 years). Then, systematic calving would resume a pattern similar to B15, with a systematic crack starting from Roosevelt Island and developing to the centerline. The ice overturning pattern generally reflects the surface flow pattern.
3.4. Wave effects on RIS calving
The RIS is assumed representative of the Antarctic ice shelves. When the extension of the frontal part is close to the limiting length or the outspread section is sufficiently necked, the stress perturbations associated with waves may lead to decisive big cracks and trigger a calving. However, for waves with amplitude two orders of magnitude smaller than shelf thickness, shortening of bk is insignificant. For example, the calving front ice thickness is 200–400 m for RIS. The M2 tides, which have magnitudes less than 20 cm around the Ross Sea, generally cannot force the shelf to respond rheologically. For extremely high waves from tsunamis, depending on the timing bk can shorten by up to 14%. As the RIS calving frequency is ~25 years and as the most recent calving was in 2000, tsunamis occurring after 2015 should not be disregarded in causing calving of RIS or at least enhancing the formation of a tabular iceberg. The big crack in an ice shelf serves only as a discontinuity for the tensile stress and the compression stress still is continuous. Thus, the buttressing effects still are not removed until the front section is completely separated from the main shelf. That occurs only after the birth of the tabular iceberg; then the buttressing effects are reduced and ice flow accelerates, preparing for next calving.
Tsunami occurrence is irregular and this study does not predict changes in its occurrence frequency. Instead, a set of sensitivity experiments is performed to investigate tsunami damage to ice shelves. Because of the convergent tectonic plate and the associated subduction zone, the Pacific Ocean has most (~80% of the global total) tsunami events. All tsunamis affect Antarctic ice shelves. Because of the bathymetry and geographical location, there are seafloor funnels surges into the mouth of Ross Sea, and the embayment’s configuration magnifies them; the RIS experiences wave-amplitude enhancement from tsunamis originating from the Indian Ocean coasts. If the wave amplitude is 15 m, the wave frequency is close to the resonance frequency of the Ross Ice Shelf, 5 × 10−5 Hz, and the limiting length of 30 km is approached. Only 70 such waves (about four tsunami events) are needed to fully fatigue and complete a tabular calving of the RIS.
Calving frequencies estimated in Figure 5 are upper limits, without considering random perturbations. In reality, calving will be finished prematurely because the fatigue from wave actions and random collisions with passing icebergs. The calving frequency also is for large, systematic calving of tabular icebergs. It by no means indicates that smaller-scale calving does not occur in between. Actually, the occurrence of smaller-scale tip-end calving is inevitable because of random collisions with floating ice driven by ocean currents and large-amplitude waves. These harbingers postpone and may even completely inhibit the occurrence of the systematic calving. For example, the 2011 Japanese tsunami and the 2012 Indonesian tsunami caused ~110 calving events of total volume of 21 km3 for the east Antarctic coast off the Cook and Shackleton Ice Shelves (A. Gavrilov 2012, personal communication). These events confirm that sea swells may be important for smaller-scale tip-end calving.
In addition to maintaining the mass balance of AIS, large tabular calving events leave strong signatures on ocean circulation and sea ice formation, affecting inland ice discharge through the buttressing effect (Thomas 1973). Understanding the calving process has long been a major unsolved problem in glaciology (Benn et al. 2007; Reeh 1994). Although several empirical calving relations have been proposed, all are tied to specific datasets and/or specific regions. A single, all-embracing calving law has proved elusive (Meier and Post 1987; Reeh 1994; van der Veen 2002). Here, a study is made of the calving of the Antarctic ice shelves in the twentieth and twenty-first centuries and the underlying viscoplastic ice fatigue mechanisms. The proposed calving scheme reconciles existing calving mechanisms and, because it uses large-scale ice-shelf properties, it is readily implementable in ice models coupled with climate models.
The primary cause of tabular calving is the vertical direction force imbalance between gravity and buoyancy, provided by ocean water for the protruding section. Crevasses that can grow and eventually generate a tabular iceberg coincide with the gravity-determined limiting extension length of the forward-slanting ice-shelf front. Tabular calving seldom is accomplished by one through-cut crack because, as the frontal part lowers, the buoyancy increases and a new balance is gained before the leading edge of the crack can penetrate the entire ice thickness. The ensuing periodic fatigue, primarily from oceanic waves, determines when the tabular iceberg occurs. As compression stress still can be transferred across large cracks and only the tensile strength diminishes, the buttressing effect lingers until complete separation of the tabular iceberg with the ice shelf. The big cracks and subsequent wave fatigue together determine the long-term calving frequency (Figure 2). Fatigue effects from wave activities are quantified by a formula involving the resonance frequencies of ice shelves. Forcings at frequencies close to the shelf resonance frequency transfer energy to ice shelf more efficiently. Irregularity of the large waves associated with tsunamis and storm surges contributes to the multiperiodicity nature of tabular calving. Over longer time scales on the order of centuries, climate control dominates the calving activities. Because shelf ice geometry is controlled by total mass balance (net snow precipitation, oceanic melt/freeze, and flow divergence/convergence), intrusion of warm currents, changed atmospheric parameters, and ice flow regimes as climate warms all affect tabular calving. Numerical modeling under realistic emission scenarios show that ice geometry adjusts so more mass is distributed close to the front, favoring tabular calving.
We thank Professor P. Bromirski for discussions about tides and swells near the Antarctic coast, T. Hughes for references he recommended, and Professor D. MacAyeal for suggesting important resources. Useful discussions/communications with Professors R. Thomas, J. Bamber, S. Jacobs, R. Bindschadler, J. Bassis, and H. Conway improved this manuscript.
Cantilever Beam Approximation for Ice-Shelf Attrition
Because of the vertical flow shear within an ice shelf, the upper layers will be stretched more than lower layers (Reeh 1968), so that the front of the ice is slanted forward, overhanging more and more with time. The section toward front thus displaces less water than its weight. This motion cannot proceed without a simultaneous downward bending of the frontal part. This bending will partially compensate the imbalance of gravity and buoyancy of the frontal section. However, if there is no pinning point, the bulk of this imbalance is supported by the inland neighboring section by submerging more than enough to support its own weight. The submerging geometry of ice shelf is thus flow dependent. We here propose a nonlocal calving scheme that considers calving as a function of ice-shelf thickness, temperature, and flow vertical shear.
When put under stress, the behavior of general material has the following stages: elastic region, permanent damage starts after reaching the yield point, plastic deformation, and necking and separation. Within the elastic range, the fatigue caused by external stress is minimal and reversible after the load is removed. The maximum strain when entering from the elastic region to the plastic region is a good measure for brittle and ductile material. The limit usually is set at 10−3. After each such strain cycle, the unit cells’ orientation gets “trained” by the external stress and becomes better organized. The following cycle usually means a raised yielding point (so-called strain hardening). Many measures such as Young’s modulus (E), Poisson ratio (ν), modulus of rigidity (G), and incompressibility (K; bulk modulus) are elastic moduli (and indicating nature of material in the elastic range).
Yielding/fracture is actually a measure of how much mechanical energy per unit of material can store, an energy form of measure of fatigue resistance capability. As shown in the simplest one-dimensional form (Figure A1), if an elastic splinter (of Young’s modulus E and length L) is deflected upward by an amount H, the work needed to be done (also the mechanical energy stored in the system) can be obtained by a simple integration,
where W is mechanical energy, F is restoring force, and the other variables are as in the illustration. Per unit length, the density of mechanical energy is . This energy intensity expression is applicable for a three-dimensional form. More interestingly, this expression has a formal similarity for fluids and can be formally extended to describe granular material as well (to be shown later). For example, in the case of fluids, E is in the form of viscosity and is the strain rate expressed in the form of flow shear. We return to the discussion of a simple yielding model, the Cantilever beam model, as indicated in the top panel of Figure A2.
As shown in Figure A2, it is assumed that the material is stiff and the deflection in the z direction is sufficiently small that linear deformation theory for elastic material is valid. It further is assumed that the cross section is rectangular to simplify the derivation by removing variations in the y direction and works only in the x–z plane (the convention for moment and torque are all in a right-handed coordinate system, as in Figure A2). The beam material has an elastic modulus E and density . Other geometrical parameters are as labeled in the figure.
Because the cross section is assumed to be rectangular, the area moment of inertia is at point O. The momentum around O exerted by the weight of the beam and an external loading (representative of tides and other random factors) T located at the tip end of the shelf is expressed as
where g is gravitational acceleration (9.8 m s−2) and M is the moment in the positive y direction. In a static state, the resistance moment should in the negative y direction and of the same magnitude. As a result of moment drive, there is potential energy stored around the cross section passing through O,
where r is the curvature of the beam at O. At the yielding condition, is the tensile strength of ice. In Equation (A3), the factor of 2 appears because assumed mass conservation is assumed, so that the cross-sectional area experiencing compression and the area experiencing expansion are the same. The relation between strain and Young’s modulus are applied, producing the factor of 2. The stored potential energy and the moment should have the same value: that is,
Equation (A4) is the master equation for obtaining the limiting length of ice shelf before attrition. For example, assuming the linear ice thickness profile,
where . For the case without external loading,
Although the shelf as a whole is in near hydrostatic balance with the ocean waters, the flow structure inside the ice shelf determines that it is a dynamic scenario of advancing–thinning–breaking, from groundling line toward the calving front. Ice shelves spread under their own weight and the imbalance pushing at the calving front and around the grounding line. According to Reeh (Reeh 1968), as a result of the horizontal compression, the geometry is that of an anvil-shaped outreach at the calving front (Thomas 1973; Hughes 1992; Hughes 2002). The length of this portion (bk) is limited by the tensile strength of ice [fc ~ 2 MPa in Equation (A7)], the ice thickness (H), and the ice creeping speed. In Equation (A7), if we take as the density of seawater, as the density of ice, and CT as a factor taking tidal and sea wind swelling into consideration, assuming the mean ice creeping speed at the calving front is U and its vertical shear is , then
and the average calving period
Suppose the calving front location right after a tabular calving is at the origin point (of a coordinate system) and the ice front advances steadily. After it advances bk length, the next calving event occurs. The term bk is also the dimension of the next tabular iceberg along the (ice) flow direction.
Equation (A8) is a relationship for a steady, systematic background calving rate. It is a physically based relationship between rift opening (and thus tabular iceberg calving) and glaciological stresses (Joughin and MacAyeal 2005). As the ice shelves are dynamic, ice moves constantly toward the calving front. The cracks/rifts upstream (C1 in Figure 2b) caused by stress larger than maximum tensile strength (i.e., at the flanks of an ice stream) or by hydrologic fracturing can eventually be advected toward the calving front (). As from Equation (A8), calving is highly ice thickness dependent. For the hydrologic-caused fracture, the water can deepen cracks and wedge through the ice shelf (Scambos et al. 2004), significantly reducing the local thickness. Downstream advection of the ice cracks makes the calving process less regular and less periodic. For Antarctica, the two largest ice shelves calve more systematically and smaller fringing shelves calve more randomly.
In reality, the ice thickness profile is flow dependent and defies an analytical expression such as Equation (A8). SEGMENT-Ice numerically solves the limiting length and calving frequency. Equation (A8) also explains why some very fast-flowing ice streams do not have extension sections as ice shelves. Although the surface elevation gradients drive the ice flow, fast-flowing ice streams that tend to have banded crevasses often are collocated with sudden bedrock gradient changes. The fast-flowing nature of the ice stream provides limited time for the compress healing of the ice shelves (Hughes 2002) when the ice is on water while still connected to the land ice. These large systematic crevasses greatly reduce the limiting length for systematic calving and discharged ice is removed by ocean water (or colliding with passing icebergs), as icebergs, shortly after leaving the grounding-line transition zone.
According to this hypothesis, accompanying a calving event, both the tabular iceberg and its parent ice shelf will adjust their respectively displaced seawater volumes. This inevitably also displaces the seabed sediments. In this sense, the proposed formula also assists paleoreconstruction of major historical calving events.
The unique frontal profile is, by scale analysis, on the order of 10 km for the major peripheral Antarctic ice shelves. To represent this, kilometer-scale data are needed. Fortunately, Griggs and Bamber (Griggs and Bamber 2011) archived a 1-km-resolution ice thickness and surface elevation dataset for Antarctica ice shelves. The data are retrieved from radar altimeter data from the geodetic phase of the European Remote Sensing Satellite (ERS) during 1994–95 supplemented by Ice, Cloud, and Land Elevation Satellite (ICESat) data (2003–09) for regions south of the ERS-1 polar hole. A new spatial variable firn-density correction scheme supposedly reduced errors in this ice thickness data. This dataset can be looked as the present state of air–ice–water configuration. As remarked above, the frontal section is where the assumption of hydrostatic equilibrium breaks down and the ice can hanging above or submerged under seawater. Figure A2b shows clearly the features anticipated from Reeh’s (Reeh 1968) perspective and our hypothesis. The linear-slanting section caused by vertical shear in ice flow is about 30 km at this section of Ross Ice Shelf. This outspread girder may just had been calved before the data collection period, likely along a crack at the location labeled by the cross sign.