Figure 1 . Study area. a . Larsen Ice Shelf (A, B, C,
and D) from north to south. The white area is the ice shelf extent in
2015. The color-coded area shows the surface topography of the grounded
region. We divided Larsen C into six flow units (C1, C2, C3, C4, C5, and
C6). b . The ice front locations in 1963, 2000, 2010, and 2017
with the ice shelf surface elevations from REMA DEM (Howat et al.,
2019), twelve flow inlets, and suture zones flowing from peninsulas and
islands.
3 Data and Methods
3.1 Mapping ice front, rifts, and flow velocities
We used multi-source satellite images acquired between 1963 and 2020 to
map the ice front, identify rifts, and derive flow velocities in
different time periods. The images include declassified intelligence
satellite photographs (DISP) acquired by the ARGON mission, optical
images acquired by the Landsat-4/5 TM, Landsat-7 ETM+, and Landsat-8 OLI
sensors, and synthetic aperture radar (SAR) images acquired by the
ERS-1, ERS-2, Radarsat-1, Envisat, and ALOS satellites. We performed
image orthorectification, image co-registration, noise reduction, and
image enhancement on the above set of images (Text S1).
We mapped the ice front locations of Larsen C for the period 1963 to
2020 (Figure S1) and delineated two rifts in the vicinity of Gipps Ice
Rise that were closely related to the 2017 calving event. We derived
flow velocities over the ice shelf for multiple periods using an
image-matching-based feature tracking method (e.g., Heid & Kääb, 2012;
Liu et al., 2012; Scambos et al., 1992). The feature tracking method
tracks the spatial displacements of moving and identical surface
features from two sequential images acquired at different times. The
derived flow velocity represents the average velocity between the
acquisition dates of sequential images. The image-matching algorithm
automatically searches for matching points between two images by
maximizing the cross-correlation coefficient between conjugate areas in
the two images. The cross-correlation-based method can achieve an
accuracy of better than 0.5 pixel (Scambos et al., 1992). For our study,
the image space was examined every 300 m to search for matching points,
with a matching window of 64x64 pixels. Automated directional and local
statistical filters (Liu et al., 2012) plus human inspection were
employed to remove erroneous matches and verify the matching results.
The accuracy of the derived velocity fields is determined by the image
co-registration accuracy, image matching accuracy, and time separation
of sequential images, which was estimated to be 3.7~25
m/year, depending on the image pairs used (Table S1).
3.2 Ice shelf modeling experiments
We used the Ice-sheet and Sea-level System Model (ISSM) (Larour et al.,
2012) to conduct the ice shelf modeling experiments. We used the shallow
shelf approximation (SSA) model (MacAyeal 1989; Weis et al. 1999;
Morland 1987) to examine how ice shelf flow velocities and stress fields
vary with front retreat and rift development (Text S2). For the model
configuration, we used the anisotropic mesh adaptation method to
generate the numerical mesh net nodes. The boundary condition was
specified as kinematic for the interior ice shelf and dynamic for the
ice shelf front. The flow law exponent n was set to 3. The ice
rigidity \(\overset{\overline{}}{B}\) was initialized as
1.5×108 Pa s1/3. For the inversion
of \(\overset{\overline{}}{B}\), we used the M1QN3 optimization
algorithm in ISSM (Morlighem et al., 2015) by minimizing a cost function
defined as the weighted sum of the difference between the modeled and
observed flow velocity fields and the gradient of the optimized\(\overset{\overline{}}{B}\) at each iteration. We further calculated
the deviatoric stresses, including longitudinal deviatoric stresses and
principal stresses (Text S1). We also calculated the backstress across
the ice shelf. Backstress describes the stress transmitted upstream that
opposes spreading of an ice shelf (Cuffey and Paterson, 2010), arising
from side drag by lateral walls and slower flowing ice, and from basal
resistance from grounded ice rises or ice rumples. These stress fields
are tightly related to rift growth and ice shelf instability
development.
4 Observations derived from satellite data
4.1 Retreat history of Larsen C during 1963–2020
During the period from August 1963 to March 2020 (Figure S1), the area
of Larsen C decreased from 55,511 km2 to 43,182
km2, a reduction of ~22%. In the six
decades, Larsen C underwent two large retreat cycles, during 1976–1986,
and during 2004–2017. Each cycle began with northern retreat near the
Bawden Ice Rise (C1–C3 units), followed by southern retreat near the
Gipps Ice Rise (C4–C6 units). For both cycles, the time interval
between the northern and southern retreat was approximately ten years.
In the first cycle, the northern retreat occurred at the C2–C3 units in
1976 due to a collision between an iceberg and the ice shelf (Ferrigno
& Gould, 1987); afterwards, two large icebergs calved off from the
C4–C6 units in 1986 (Skvarca, 1994), causing a retreat of
~100 km in the south. The ice shelf slowly advanced
between 1986 and 2004. In the second cycle, the northern front detached
from the C2–C3 units in 2004/2005, with an area loss of
~1500 km2. The major retreat of the
southern front happened on 12 July 2017, when a
~5800-km2 tabular iceberg detached
from the C3–C6 units.
4.2 Rifting related to the 2017 calving event
The 2017 calving event resulted from the propagation of two rifts
(labeled R1 and R2, Figure 2). Long and transverse fractures (Figure 2a)
tend to form between Kenyon Peninsula and Gipps Ice Rise due to a stress
condition change, where ice shelf flow becomes less constrained by the
lateral shear margins after passing by Kenyon Peninsula. The majority of
the fractures terminate near the Joerg suture zone (Figure 2a). Suture
zones can inhibit fracture growth because of the material properties of
marine ice underneath (Kulessa et al., 2014, 2019; McGrath et al.,
2014). The rift R1, a transverse rift propagating perpendicular to the
flow direction, developed from a preexisting fracture (L1 in Figure 2a)
that is identifiable on the 1963 DISP imagery. The lateral extent of R1
changed little until 2008, when it began to propagate towards R2. By
2010, the southern rift tip of R1 had traveled ~14 km
from its 2006 location (from A’ to A, Figure 2c), to a point less than 1
km from R2. At the northern end, R1 developed a new rift trajectory from
point B’ to B (Figure 2c), and subsequently propagated to the north
along this trajectory. R2 was an oblique rift that formed between Gipps
Ice Rise and a transverse fracture (L2 in Figure 2a). R2 is not present
on the 1963 DISP imagery. Its formation was very likely due to the front
geometry change caused by the 1986 retreat. R2 widened over time, while
its rift tips exhibited little change until 2013. Then, it began
developing towards the south, and the rift tip almost reached Gipps Ice
Rise (Figure 2e). Simultaneously, R1 penetrated through the Joerg suture
zone and started rapid propagation to the north (Figure 2d). R1 cut
through the C5 unit in 2014 and had propagated across the Tonkin suture
zone by January 2016 (Figure 2d). By March 2017, the length of R1
dramatically increased by ~65 km, and its northern rift
tip (‘B’) reached the Cole Peninsula suture zone (connecting the C2 and
C3 units). R1 continued penetrating through the Cole Peninsula suture
zone, but did not extend into the C2 unit. Instead, it was developing
along the C2 unit boundary towards the ice front, resulting in the
complete detachment of the downstream portion in July 2017.
The propagation of R1 appears to
have had a strong association with the development of R2. The extension
of R2 to Gipps Ice Rise was immediately followed by the rapid growth of
R1 towards the north. This suggests that rifting can promote the growth
of another rift by weakening the lateral shear margins near an ice rise.
Ice rises are locally grounded topographic highs in ice shelves. They
often act as pinning points, buttress upstream ice flows and play an
active role in ice shelf stability (Favier & Pattyn, 2015; Matsuoka et
al., 2015). Weakening around ice rises could reduce their buttressing
effect, leading to enhanced longitudinal stretching and creating
favorable conditions for rift propagation. We further explore this
hypothesis through modeling experiments in Section 5.