Sensitivity of the ice-shelf/ocean system to the sub-ice-shelf cavity shape measured by NASA IceBridge in Pine Island Glacier, West Antarctica

Abstract

Two high-resolution (1 km grid) numerical model simulations of the Amundsen Sea, West Antarctica, are used to study the role of the ocean in the mass loss and grounding line retreat of Pine Island Glacier. The first simulation uses BEDMAP bathymetry under the Pine Island ice shelf, and the second simulation uses NASA IceBridge-derived bathymetry. The IceBridge data reveal the existence of a trough from the ice-shelf edge to the grounding line, enabling warm Circumpolar Deep Water to penetrate to the grounding line, leading to higher melt rates than previously estimated. The mean melt rate for the simulation with NASA IceBridge data is 28 ma^–1, much higher than previous model estimates but closer to estimates from remote sensing. Although the mean melt rate is 25% higher than in the simulation with BEDMAP bathymetry, the temporal evolution remains unchanged between the two simulations. This indicates that temporal variability of melting is mostly driven by processes outside the cavity. Spatial melt rate patterns of BEDMAP and IceBridge simulations differ significantly, with the latter in closer agreement with satellite-derived melt rate estimates of ~50ma^–1 near the grounding line. Our simulations confirm that knowledge of the cavity shape and its time evolution are essential to accurately capture basal mass loss of Antarctic ice shelves.

Introduction

Large basal melt rates, mass loss and grounding line retreat of ice shelves of the West Antarctic ice sheet (WAIS) have been reported in many studies over the last few years (Jacobs and others, 1996; Hellmer and others, 1998; Rignot, 1998; Thomas and others, 2004; Pritchard and others, 2009). According to Shepherd and others (2004), the WAIS lost ~51±9km³ a^–1 of volume during the 1992–2001 period. Ocean warming is listed as a potential source for increased basal melting and ice-shelf thinning, which in turn leads to grounding line retreat. The influence of the ocean on the melt rate of the Amundsen Sea ice shelves can be attributed to the Circumpolar Deep Water (CDW), which is ~3˚C warmer than the surface freezing-point temperature and which flows onto the shelf through depressions in the depth of the continental shelf break (e.g. Jacobs and others, 1996). To understand the processes that determine the melt rate, Jenkins and others (2010) find that it is important to know the shape of the sub-ice cavity and the distribution of water properties inside and outside the cavity.

Hellmer and others (1998) showed temperature and salinity sections in the eastern part of Pine Island Bay to indicate the onshore flux of CDW and a circulation scheme inside and outside the Pine Island ice shelf cavity. We refer to this easternmost trough as the eastern channel (EC, Fig. 1). Walker and others (2007) reported an onshore flow of CDW through the central channel (CC) and its potential to melt the Pine Island ice shelf. More recently, Wåhlin and others (2010) reported an onshore flow through the western channel (WC). In this study we find that, due to bathymetric constraints, the major path for onshore flow of warm water is the EC.

Fig. 1. (a) Bathymetry (m) of the model domain in the Amundsen Sea. The thick black line depicts the ice-shelf edge. The x- and y-axes show model grid spacing in km. (b) Ice-shelf thicknesses (m) for the Abbot (AB), Cosgrove (CG), Pine Island Glacier (PIG), Thwaites (TH), Crosson (CR), Dotson (DT) and Getz (GZ) ice shelves. SI: Sherman Island; BI: Burke Island; EC: eastern channel; CC: central channel; WC: western channel.

Previous studies of the Amundsen Sea either considered the ocean circulation in front of the cavities (Thoma and others, 2008) or were restricted to the ice-shelf cavity (Payne and others, 2007; Bindschadler and others, 2011). In this study, we combine both aspects, i.e. the onshore flow of CDW and its direct effect on the melting behavior of Amundsen Sea ice shelves, with emphasis on the Pine Island ice shelf. Specifically, we explore the role of the ocean in the recently observed mass loss and grounding line retreat of the Pine Island ice shelf, study the sensitivity of melt rate due to changes in the sub-ice-shelf cavity bathymetry and explore possible links between melting patterns and onshore flow in the Amundsen Sea.

Model Description

We use the Massachusetts Institute of Techology general circulation model (MITgcm), which includes a dynamic/ thermodynamic sea-ice model (Losch, 2010). The freezing/ melting process in the sub-ice-shelf cavity is represented by the three-equation thermodynamics of Hellmer and Olbers (1989), with modifications following Jenkins and others (2001) as implemented in MITgcm by Losch (2008). The exchange of heat and fresh water at the base of the ice shelf is parameterized as diffusive flux of temperature and salinity using the constant friction velocity and turbulent exchange coefficients of Holland and Jenkins (1999).

The model domain (Fig. 1a) is located in the Amundsen Sea and is derived from the global cube sphere grid configuration used by the Estimating the Circulation and Climate of the Ocean, Phase II (ECCO2) project (Menemenlis and others, 2008), but with horizontal grid spacing of ~1 km instead of ~18 km in the global configuration. The vertical discretization is the same as that used by the ECCO2 project, i.e. 50 vertical levels of varying thickness (10m at the surface to 450m at the deepest level; 70–90m in the 500–1000m depth range).

Two bathymetry datasets were used in this study. The first is a combination of Nitsche and others (2007) with Smith and Sandwell (1997) topography version 13.1 blended with the BEDMAP bed elevation (Lythe and others, 2001) in the sub-ice-shelf cavities. The blend between these two datasets produces a step at ~100˚W, which is the edge of the Nitsche and others (2007) bathymetry. In the second dataset, the Pine Island sub-ice-shelf cavity bathymetry was extracted from the most recent NASA IceBridge dataset (Fig. 2; Studinger and others, 2010). As shown in the figure, we separate the Pine Island ice shelf into a southern region called PI-South, which is fed by Pine Island Glacier (PIG), and a northern region called PI-North.

Fig. 2. Cavity bathymetry (m) of PIG derived from (a) BEDMAP and (b) IceBridge data. Gray shaded areas represent grounded ice/bedrock. Dashed lines show the separation of the cavity into northern and southern parts as noted in Table 1. The x- and y-axes show model grid spacing in km.

Although all available ice shelves of the Amundsen Sea (Fig. 1b) are included in the model domain, there are large uncertainties in the cavity shape, especially those that do not include IceBridge data. Our choice of bathymetry data and grounding line location in the model domain may explain differences of ice-shelf area sizes (Table 1) compared to the estimates from Griggs and Bamber (2011) (e.g. the Thwaites Glacier ice shelf area is slightly larger due to the inclusion of a grounded iceberg to the model ice shelf).

Initial conditions are derived from the World Ocean Atlas 2005 (Boyer and others, 2006). Inconsistencies in the initial condition temperature fields as a result of the extrapolation are reflected in very high melt rates during the first few months of the integration. Temperature and salinity cavity values adjust quickly, however, as reflected in the stabilization of the melt rates in both model integrations.

Table 1. Ice-shelf area used in the model, derived from Ice, Cloud and land Elevation Satellite (ICESat)/Geoscience Laser Altimeter System (GLAS) data. Time mean (1981–2010), area-averaged freshwater flux (FWF), basal ice mass flux (BMF) and total basal mass loss (BML) using an ice density of 917 kgm^–3. Note the differences compared with Griggs and Bamber (2011): Abbot and Getz Ice Shelf areas are smaller in the model due to the domain size; Pine Island ice shelf area is smaller as Griggs and Bamber (2011) include areas of the ice shelf that appear to be grounded; only Dotson Ice Shelf is referenced in Griggs and Bamber (2011); in the model the ice shelf is separated into Crosson and Dotson Ice Shelves; Thwaites is larger in the model as it includes a grounded iceberg in front of the ice shelf. 1 Sv = 106m³ s^–1

Boundary conditions for hydrography and sea ice were provided by a coarse-resolution (18km grid) circumpolar integration that includes ice-shelf/ocean interaction, which itself obtained boundary conditions from a global ECCO2 ocean state estimate (Menemenlis and others, 2008). Due to the difference in resolution between the 18 km circumpolar integration and the 1 km Amundsen Sea domain, a relaxation (30 gridpoints into the model domain) was applied to temperature and salinity at the boundaries to avoid artifacts such as wave energy radiating into the model interior.

Surface forcing is provided by the Japan Meteorological Agency and Central Research Institute of Electric Power Industry 25 year reanalysis (JRA-25; Onogi and others, 2007), and covers the integration period 1979–2010.

Results

The bathymetry of the Amundsen Sea (Fig. 1) depicts the features controlling the onshore flow of warm CDW onto the shelf. Three channels, the WC (~117˚ W), the CC (~113˚ W) and the EC (~103˚ W), along the shelf break assist CDW penetration onto the shelf, and a deep trough in the east guides it towards PIG. A trough further west guides CDW towards the Getz and Dotson Ice Shelves. The modeled temperature at 552 m depth, the core of the CDW, shows that onshore flow can be observed in the EC and WC (Fig. 3). The warm-water transport in the CC is restricted to the first few kilometers, after which the flow is diverted by bathymetry and by the circulation, which shows a strong northward movement and which effectively blocks the onshore transport of warm CDW in the CC. While Thoma and others (2008) found the main inflow along CC with an eastward deflection before moving south towards the ice shelf, we find that the water entering the CC is mixed with and cooled from waters moving northward. The main heat inflow in our model is through the EC. The bathymetry at the EC favors baroclinic instabilities with warm core eddies and/ or warm filaments transporting heat onshore, feeding the heat towards the ice shelves. The differences of onshore pathways in our model compared to Thoma and others (2008) may be due to differences in bathymetry. Thoma and others (2008) did not use the Nitsche and others (2007) bathymetry, so their model configuration does not contain a depression along 105˚W, which serves as the main pathway for CDW onto the continental shelf in our study.

Fig. 3. Potential temperature (˚C) snapshot at 552m depth (December 1988) depicting the main pathway of warm water along the EC from the continental shelf break towards the PIG cavity. The white line depicts the grounding line, the black line the surface coastline and ice-shelf edge, and the light-gray shaded areas represent the sea-floor at 552 m depth. The x- and y-axes show model grid spacing in km. EC: eastern channel; CC: central channel; WC: western channel.

Simulated temperature higher than 1˚C penetrates into the sub-ice cavity (Fig. 3), as was observed by Jenkins and others (2010). Because the cavity is larger and deeper in the IceBridge configuration, however, the warm waters reach closer to the grounding line. For the IceBridge simulation we obtain a 1981–2010 average melt rate of 28 ma^–1 for the ice shelf in front of PIG, that is, a mass loss of 118 Gt a^–1 (Table 1). The southern part of the Pine Island cavity has an area-averaged melt rate of 33 ma^–1 for 1994, 50% higher than the 1994 estimate of Jacobs and others (2011) and 37% higher than the 1992–96 melt rate of 24 ±4ma^–1 estimated by Rignot (1998). In contrast to this, the simulation with the BEDMAP cavity shape yields a 1981–2010 area-averaged melt rate of 21ma^–1 and a 1994 melt rate of 24 ma^–1 for the southern cavity. The 1994, BEDMAP-based estimate compares well with the estimates of Rignot (1998) and Jacobs and others (2011) for the same period. While for the early part of the simulation period the BEDMAP results are closer to observations, the reverse is true during the later part of the simulation period. In 2009 the BEDMAP simulation yields a melt rate of 20 ma^–1, which is low compared to the 33 ma^–1 melt rate estimate of Jacobs and others (2011) and the 37±4ma^–1 melt rate estimate from interferometric synthetic aperture radar (InSAR)/Regional Atmospheric Climate Model (RACMO) data (Rignot and others, 2008). The 2009 melt rate based on the IceBridge simulation is 31 ma^–1 and much closer to observations.

The differences in melt rate between the BEDMAP and IceBridge simulations are primarily caused by differences in the southern part of the cavity, where major changes to the bathymetry are found (Fig. 4c). The northern part of the cavity has a higher variability and is the principal contributor to the multiple peaks in the melt rate curve. These peaks result from surface forcing and, in particular, from stronger surface winds (e.g. during the 1993 event), which lead to cooler water masses entering the northern cavity and thus reducing the melt rate. This is similar to what Holland and others (2010) found for the melting behavior of the George VI Ice Shelf in the Bellingshausen Sea. Both the BEDMAP and IceBridge simulations show the same interannual variability, with increasing melt rates until 1995, followed by a decrease until 1999, after which they increase again (Fig. 4).

Fig. 4. (a) Area-averaged melt rate (m a^–1) of Pine Island ice shelf from 1979 to 2010 with BEDMAP (gray) and IceBridge (black) bathymetry in the sub-ice-shelf cavity. Estimates from observations are included: triangle and square after Jacobs and others (2011); diamond with error bars after Rignot (1998). (b, c) Area-averaged melt rate for BEDMAP (b) and IceBridge (c) simulation divided into northern (gray) and southern (black) parts of the Pine Island ice shelf.

To complete the picture of the Amundsen Sea, the melt rates of the remaining ice shelves are listed in Table 1. These melt rates are highly uncertain because the cavity bathymetry and other empirical model parameters are unknown. For example, the simulated basal mass loss of the Cosgrove Ice Shelf is 32 Gt a^–1, which is high compared to the 15 Gt a^–1 combined influx of Smith and Pole Glaciers, which feed the Cosgrove Ice Shelf (Rignot and others, 2008).

The circulation on the Amundsen shelf is dominated by a cyclonic gyre with a mean transport of 3±0.5 Sv (Fig. 5a). The ridge that extends from Bear Peninsula towards the shelf break separates the flow of Pine Island Bay from a second cyclonic gyre that bathes the Dotson and Getz Ice Shelves. The strength of the cyclonic gyre that flows into Pine Island Bay increases until 1995, then decreases until 2001 and increases again thereafter (Fig. 5c). The horizontal stream-function and melting show a correlation of r = 0.42, which is significant at the 95% confidence level.

Fig. 5. (a) Mean horizontal streamfunction (1996 as typical year). Bathymetry contours are depicted in 500m intervals. In addition, thick solid lines are included at the 500 and 600m contours, and dashed lines at the 550 and 650m contours. The red line shows the ice-shelf edge. (b) Difference of horizontal streamfunction 2001– 1992, weaker vs stronger gyre circulation. The thin black line indicates the bathymetry, and the solid black line the coastline and ice-shelf edge. (c) Streamfunction time series (black) on the continental shelf of Pine Island Bay just outside the cavity, and 30 day mean, volume-averaged temperature of the PIG ice shelf cavity (gray) for the IceBridge simulation. See Figure 1a for annotations. The x- and y-axes show model grid spacing in km.

Discussion

Mean melt rate

The mean melt rate for the IceBridge simulation (southern Pine Island ice shelf) of 30.5 ma^–1 is higher than estimates from hydrography and numerical models. Jacobs and others (1996) estimated a melt rate of 10–12ma^–1 and a mass loss of 28 Gt a^–1, but they used a different position for the grounding line and therefore a biased estimate of the ice-shelf area. Model estimates by Hellmer and others (1998) of 6–12.5ma^–1 depend on cavity configuration and on source water characteristics and represent a conservative estimate of the basal melt rate.

Melt rate estimates have uncertainties due to water column thickness, grounding line location and ice-shelf area. For example, part of the discrepancy between our results and those of Jacobs and others (1996) and Hellmer and others (1998) may be due to differences in the ice-shelf areas. There are also uncertainties in water column thickness due to the limitations of deriving bathymetry from gravity observations. For example, comparing the IceBridge bathymetry with the lines from Autosub (Jenkins and others, 2010) reveals differences of ±96m (~10% of water column depth) in seabed topography on length scales that are shorter than the spatial resolution of the airborne gravity data. The IceBridge data, however, provide much better spatial coverage than the Autosub data. Ideally, Autosub data should be used to adjust IceBridge inversion parameters (e.g. bedrock density and sediment thickness).

Temporal variability

The total melt rates of the BEDMAP and IceBridge simulations differ significantly (by ~7ma^–1), but the temporal variability remains almost unchanged (Fig. 4a). The main difference is the deeper inner part of the cavity that allows warmer water to reach the now extended grounding line and thus increases the mean melt rate (Fig. 6). These model results indicate that most of the temporal variability is driven by processes outside the cavity. Our results show that the melt rate is correlated with the heat transport, which is determined by the horizontal streamfunction on the continental shelf (Fig. 5). The BEDMAP and IceBridge simulations both show an increase in melting until 1995. This coincides with Rignot’s (2002) finding of an acceleration of PIG in the 1990s, which he attributed to a thinning of PIG due to a possible increase of intrusions of warm CDW on the continental shelf. Our simulations show that these intrusions occur regularly via eddies and filaments and that their properties are set at the boundaries of the model domain.

Fig. 6. Horizontal distribution of mean melt rate (m a^–1) for BEDMAP (a) and IceBridge (b) bathymetries for the period 1981–2010. The dashed line depicts the partition into northern and southern cavities. The x- and y-axes show model grid spacing in km.

The significant correlation (r = 0.42) between the melt rate and the strength of the Pine Island Bay streamfunction indicates that the increased basal melt rates of the PIG ice shelf may be due, in part, to local circulation changes. For example, Figure 5b shows the streamfunction difference between periods of low and high melt rate conditions. On the one hand, there is little change in the main part of the gyre, which is adjacent to (but does not cross) the shelf break, indicating that the flux across the shelf break is not coupled to the continental shelf circulation. On the other hand, the strength of the smaller gyre, just outside the PIG cavity, is amplified during the period of high melt rate, indicating that circulation on the shelf is more important than the initial provision of CDW across the shelf break. This model result is consistent with the argument of Jacobs (2006) and Holland and others (2008) that atmospheric variability, rather than warming of the CDW, can have the larger impact on ice-shelf melting.

Our study finds maximum heat transport to the sub-ice-shelf cavity at the end of summer and beginning of autumn, which is associated with a maximum in the streamfunction of the small gyre, just outside the PIG cavity. This result differs from the Thoma and others (2008) study, which reports that the inflow of warm CDW on the shelf is strongest during the winter and spring. Thoma and others (2008) did not have coupled ice-shelf/ocean thermodynamics in their model, so they could not quantify the influence of CDW on glaciological changes. In addition, the Thoma and others (2008) study was carried out on a coarser horizontal grid (13.5–28 km), which did not resolve eddies, and it used an older bathymetry, which did not contain the EC. Thoma and others (2008) also report that there is some correlation between onshore transports and Southern Annular Mode (SAM) and El Niño Southern Oscillation (ENSO) in their modeling efforts, which we have not investigated in this study. The mechanism that drives CDW waters on the continental shelf, and the relationship between observed melting rate variability and SAM and ENSO atmospheric conditions require further investigation.

Melting pattern

Although the melt rate temporal variability of the BEDMAP and IceBridge simulations is the same, the spatial melt rate pattern is significantly different (Fig. 4). In their numerical model, Payne and others (2007) found that the maximum melt occurs in a narrow zone within 20km of the grounding line, where the melt rates exceed 100ma^–1. In both of our simulations this extreme melting value only occurs at the beginning of the integration, during the initialization phase, before the cavity temperature and salinity properties have adjusted to basal melting. In the simulation with IceBridge bathymetry, where warm waters can penetrate deeper into the cavity, the mean melt rates at the grounding line are higher by ~20ma^–1 than in the BEDMAP simulation (a mean of 47 ma^–1 vs 24 ma^–1). This value is in agreement with the grounding line melt rates of ~50ma^–1 observed by Rignot (1998), and with the most recent estimates by Bindschadler and others (2011). Contrary to Payne and others (2007), the maximum melt rate at all model gridcells near the grounding line, for both of our simulations, never exceeds 75 ma^–1 after the model spin-up period.

Dynamic ice shelf

With remote sensing, area-averaged melt rates of 24±4ma^–1 for 1994 (Rignot, 1998) and 37±4ma^–1 for 2009 are estimated. Using hydrographic data, Jacobs and others (2011) estimated similar melt rates for these two years: 22 ma^–1 for 1994 and 33 ma^–1 for 2009. The melt rate for 1994 is closer to the estimates of the BEDMAP simulation, whereas the melt rate for 2009 is closer to that of the IceBridge simulation. While there are many potential reasons for this result, we suggest that transient changes in the shape of the cavity and grounding line retreat, which are due to the dynamics of the ice shelf as influenced by calving and basal melting, are not captured with our model. A dynamic ice-shelf model (e.g. Joughin and others, 2010) would allow changes in cavity configuration, which would alter the melting rate pattern, possibly leading to thinning of the ice shelf, and hence to acceleration of the ice stream as more water could penetrate the cavity and increase the melt rate. Joughin and others (2003) found that the largest acceleration of PIG occurred after 1996, indicating that large changes in cavity shape must have occurred around that time. Starting the numerical integration from a given grounding line, in this case the 1996 grounding line, means that previous changes to the shape of the cavity were neglected. For example, the IceBridge simulation overestimates melt rates relative to the 1994 observations (Rignot, 1998; Jacobs and others, 2011) because the water column thickness is overestimated in the years prior to 1996. To address this model deficiency, we are in the process of coupling the MITgcm to the Ice Sheet System Model (ISSM; Morlighem, 2010).

Summary and Concluding Remarks

Numerical simulations carried out using the BEDMAP and IceBridge bathymetries showed that a trough detected by the IceBridge campaign enables warm CDW to penetrate to the grounding line of the Pine Island ice shelf. The area-averaged melt rate for the IceBridge simulation (28ma^–1) for the period 1981–2010 is 25% higher than for the BEDMAP simulation (21ma^–1). The BEDMAP simulation may be closer to the true cavity shape in 1979, but does not reflect changes in water column thickness due to accelerated melting and associated dynamic thinning of the ice shelf in the 1980s and 1990s and therefore underestimates the melting in the later period of the simulation. The IceBridge simulation, on the other hand, starts from a cavity shape that exaggerates the water column thickness and overestimates the melt rates in the earlier period of the simulation, when compared to estimates from hydrography and remote sensing. This interpretation relies on the accuracy of our parameterization of basal melt rate and other aspects of the model set-up, including atmospheric forcing, open boundary conditions and sea ice. However, the close agreement between modeled basal melt rates and those estimated from satellite measurements is encouraging.

Spatial patterns of time-averaged (1981–2010) melt rates for the BEDMAP and IceBridge simulations differ significantly, with the latter in closer agreement with satellite-derived melt rate estimates (~50ma^–1) near the grounding line (Rignot, 1998). Despite large differences between the BEDMAP and IceBridge bathymetries, mean melt rates and melt rate patterns, the temporal evolution of the melt rate remains unchanged (Fig. 4), indicating that temporal variability is mostly driven by processes outside the cavity. The two simulations show that knowledge of the cavity shape and its temporal evolution is essential to accurately capture the basal mass loss of Antarctic ice shelves. In turn, this mass loss can lead to glacier acceleration and possibly contribute to sea-level change. Therefore, measuring the shape of the sub-ice-shelf cavity from remote-sensing instruments, as is done by the NASA IceBridge gravity-derived bathymetry under the Pine Island ice shelf, is critical to enabling accurate simulations of sea-ice, ice-shelf and ocean interactions.