Hostname: page-component-7ccbd9845f-s2vjv Total loading time: 1.178 Render date: 2023-01-27T09:08:19.483Z Has data issue: true Feature Flags: { "useRatesEcommerce": false } hasContentIssue true
>Journal of Glaciology
>FirstView
>A site for deep ice coring at West Hercules Dome: results...
Journal of Glaciology
Geophysical data and model setup
Geophysical observations
Kinematic model results
A site for deep ice coring at West Hercules Dome: results from ground-based geophysics and modeling
T. J. Fudge [Opens in a new window] ,
Benjamin H. Hills [Opens in a new window] ,
Annika N. Horlings [Opens in a new window] ,
Nicholas Holschuh [Opens in a new window] ,
John Erich Christian ,
Lindsey Davidge ,
Andrew Hoffman ,
Gemma K. O'Connor ,
Knut Christianson and
Eric J. Steig [Opens in a new window]
T. J. Fudge*
Department of Earth and Space Sciences, University of Washington, Seattle, WA, USA
Benjamin H. Hills
Department of Earth and Space Sciences, University of Washington, Seattle, WA, USA Applied Physics Lab, University of Washington, Seattle, WA, USA
Annika N. Horlings
Nicholas Holschuh
Department of Geology, Amherst College, Amherst, MA, USA
John Erich Christian
School of Earth and Atmospheric Sciences, Georgia Tech, Atlanta, GA, USA University of Texas Institute of Geophysics, Austin, TX, USA
Lindsey Davidge
Andrew Hoffman
Gemma K. O'Connor
Knut Christianson
Eric J. Steig
Author for correspondence: T. J. Fudge, E-mail: tjfudge@uw.edu
Save PDF (5 mb) View PDF[Opens in a new window]
Hercules Dome, Antarctica, has long been identified as a prospective deep ice core site due to the undisturbed internal layering, climatic setting and potential to obtain proxy records from the Last Interglacial (LIG) period when the West Antarctic ice sheet may have collapsed. We performed a geophysical survey using multiple ice-penetrating radar systems to identify potential locations for a deep ice core at Hercules Dome. The surface topography, as revealed with recent satellite observations, is more complex than previously recognized. The most prominent dome, which we term 'West Dome', is the most promising region for a deep ice core for the following reasons: (1) bed-conformal radar reflections indicate minimal layer disturbance and extend to within tens of meters of the ice bottom; (2) the bed is likely frozen, as evidenced by both the shape of the measured vertical ice velocity profiles beneath the divide and modeled ice temperature using three remotely sensed estimates of geothermal flux and (3) models of layer thinning have 132 ka old ice at 45–90 m above the bed with an annual layer thickness of ~1 mm, satisfying the resolution and preservation needed for detailed analysis of the LIG period.
Antarctic glaciologyice coreradio-echo sounding
Journal of Glaciology , First View , pp. 1 - 13
DOI: https://doi.org/10.1017/jog.2022.80[Opens in a new window]
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution, and reproduction in any medium, provided the original work is properly cited.
Copyright © The Author(s), 2022. Published by Cambridge University Press
The West Antarctic ice sheet (WAIS) is vulnerable to rapid retreat because it is grounded below sea level (e.g. Mercer, Reference Mercer1978; Bamber and others, Reference Bamber, Riva, Vermeersen and LeBrocq2009), and has collapsed at least once in the Pleistocene (Scherer and others, Reference Scherer1998). During the Last Interglacial (LIG, 116–132 ka), global eustatic sea level was higher than today, though the timing and magnitude of the sea level highstand remain uncertain. Estimates of the upper limit of sea level during the LIG range from 4 m (Dyer and others, Reference Dyer2021) to 9 m higher than today (Dutton and others, Reference Dutton2015). A collapse of the WAIS would help explain the LIG sea level, but there is no unequivocal evidence that this occurred.
Ice core records provide one possible means to obtain information about past changes in the WAIS. Existing ice core records from West Antarctica either do not contain ice dating to the LIG (e.g. WAIS Divide, Buizert and others, Reference Buizert2015) or, if they do, the ice is very near the bed and is difficult to interpret (e.g. Siple Dome, Brook and others, Reference Brook2005; Roosevelt Island, Lee and others, Reference Lee2020; Skytrain Ice Rise, Mulvaney and others, Reference Mulvaney2020). An ice core site from a location that would be sensitive to the climatic changes induced by WAIS collapse, but where the ice flow would not be directly affected, has the potential to elucidate the size and rate of change of the WAIS during the LIG.
Hercules Dome (86° S, 105° W) is located in East Antarctica ~400 km from the South Pole towards the Horlick and Thiel ranges in the Transantarctic mountains. The primary moisture source region for snowfall at Hercules Dome is the Pacific sector, such that the moisture pathway passes over the WAIS; in contrast, the moisture pathways for existing East Antarctic ice cores that preserve LIG records do not pass over the WAIS (Sodeman and Stohl, Reference Sodeman and Stohl2009; Nicolas and Bromwich, Reference Nicolas and Bromwich2011). Hercules Dome is not likely to have changed elevation significantly in elevation from the LIG to present even in scenarios of WAIS collapse (Sutter and others, Reference Sutter, Gierz, Grosfeld, Thoma and Lohmann2016), with changes of <75 m in some models (e.g. Pollard and DeConto, Reference DeConto and Pollard2016; Garbe and others, Reference Garbe, Albrecht, Levermann, Donges and Winkelmann2020). On the other hand, the composition of snowfall at Hercules Dome should be sensitive to changes in the elevation of the WAIS (Steig and others, Reference Steig2015). Simulations with multiple atmospheric models show that a lowered WAIS topography would cause a cyclonic anomaly which changes regional atmospheric circulation, bringing a greater flow of warm, moist air across the WAIS toward the South Pole (Steig and others, Reference Steig2015; Holloway and others, Reference Holloway2016). Based on these simulations, an ice core from Hercules Dome would provide a climate record that is influenced by the size and extent of the WAIS.
Climate records from a deep ice core at Hercules Dome could be used in conjunction with other paleoclimate records to constrain the size and rate of change of the WAIS during the LIG. Indicators of past ice-sheet extent, such as subglacial bedrock exposure (e.g. Hein and others, Reference Hein2016; Spector and others, Reference Spector2018), lack information on the rate of change which is likely preserved in ice core records. A record from Hercules Dome would complement other ice core records both directly affected by WAIS topography changes, such as Skytrain Ice Rise (Mulvaney and others, Reference Mulvaney2020), and indirectly affected, such at Mount Moulton (Korotkikh and others, Reference Korotkikh2011). Indeed, Steig and others (Reference Steig2015) showed if the WAIS collapsed, surface temperatures at Mount Moulton would have warmed less, and Hercules Dome would have warmed more, than sites farther inland on the East Antarctic plateau. Therefore, existing ice core records from the East Antarctic plateau provide independent information against which the Hercules Dome record could be compared.
Hercules Dome was first identified as a promising location for a deep ice core during the 2001/02 US-ITASE traverse from West Antarctica to South Pole (Jacobel and others, Reference Jacobel and Welch2005). Ice-penetrating radar was collected along the main traverse route as well as two perpendicular 50 km long tracks. The radar revealed clear internal layering with ice thicknesses that ranged from ~1500 to 2800 m. A 72 m firn core was also recovered. Analyses of the firn core show that annual layering and volcanic signals are well preserved and the 400 year average accumulation rate was 0.12 m a−1 ice equivalent (Dixon and others, Reference Dixon2013).
We conducted a geophysical survey across the Hercules Dome region during the 2018/19 and 2019/20 austral summers. Recent satellite observations (Cryosat-2, ICESat2 and Maxar Imagery) revealed a complex surface topography (Fig. 1). There are three local surface divides, which we refer to here by their relative positions in polar-stereographic space with respect to grid north (Fig. 1). 'West Dome' (85.8° S, 102.9° W) lies furthest to the true north, and is connected to 'East Dome' by a ~70 km ridge. A secondary ridge extends approximately grid south from East Dome that we term 'South Ridge'. East Dome, where the US-ITASE traverse visited, has large variations in bedrock topography and lacks a clear dome summit; East Dome is the current source region for South Ridge. We therefore focus on West Dome where our geophysical data suggest a simple flow regime and small bedrock topography variations (Fig. 2), and assess the suitability of West Dome as a site for recovering a deep ice core that preserves LIG ice. We focus on a 10 km long radar line (named x135 within our grid) which nearly transects the dome summit, and describe the results of ice-flow modeling constrained by near-surface, deep-profiling and phase-sensitive radar.
Fig. 1. Overview map of the Hercules Dome geophysical survey in polar stereographic (PS) coordinates. The directions in the PS projection are used in the remainder of the paper. Arrow indicates direction of geographic north. The inset shows the two radar lines, x135 and x133, at West Dome. Surface elevation contours from the Reference Elevation Model of Antarctica (REMA, Howat and others, Reference Howat, Porter, Smith, Noh and Morin2019) are shown at 10 m spacing for the overview and 2 m spacing for the inset.
Fig. 2. (a) GNSS derived surface elevation along x135 with ApRES sites shown by colored dots: n0 (yellow), n1000 (orange), n2000 (pink), n5000 (purple). (b) VHF radar of line x135 with arrows pointing to bright layer identifiable throughout survey. Data quality is lower to the north. (c) The depth of the bright layer on x135 agrees well with x133, which does not have the data gap or quality issues. VHF data along x133 is used for the accumulation forcing in the kinematic model. (d) HF radar of x135 showing internal layer and bed reflections. (e) Traced internal layers of x135 radargram shown in (d).
2. Geophysical data and model setup
The three radar systems used at Hercules Dome in the 2018/19 and 2019/20 field seasons as well as the kinematic model setup used to help interpret the observations are described in this section; Section 3 will present the results. The extent of each survey is indicated in Figure 1. We use stereographic grid directions, following the orientation of Figure 1, in the remainder of the paper. This paper focuses on the radar data collected near the West Dome summit. Two 10 km long lines across the divide were collected (Fig. 1, inset); x135 is closest to the dome summit and x133 is 2 km to the east.
2.1 Deep (HF) radio-echo sounding surveys
Radio-echo sounding (RES) has been used to map ice thickness and bed elevation, basal reflectivity and ice-sheet internal stratigraphy along two 10 km lines at West Dome (x135 in Fig. 2). We used a custom high-frequency (HF, 3 MHz) monopulse radar system (Welch and others, Reference Welch and Jacobel2003, Reference Welch, Jacobel and Arcone2009; Christianson and others, Reference Christianson, Jacobel, Horgan, Anandakrishnan and Alley2012, Reference Christianson2016). The RES data were processed using standard techniques, including horizontal stacking, bandpass filtering, correction for antenna separation, geolocation from contemporaneous kinematic Global Navigation Satellite System (GNSS) data, interpolation to standard trace spacing and time-wavenumber migration using the open-source radar-processing and interpretation toolbox ImpDAR (Lilien and others, Reference Lilien, Hills, Driscol, Jacobel and Christianson2020). Englacial reflections were picked through the semi-automated digitization routine in ImpDAR. The two-way travel times (TWTT) were then converted to depth using a spatially uniform firn-density profile (ITASE 02-4, Dixon and others, Reference Dixon2013) and a Looyenga and others (Reference Looyenga1965) mixing model (Fig. 2d). There is a ~100 m data gap at 4.9 km on x135. The noise on the picked layers in Figure 2e is due to sample resolution and the picking algorithm in ImpDAR.
2.2 Shallow (VHF) RES surveys
To image the upper 20–200 m of the ice sheet, we used a Pulse-EKKO Pro radar system with very high-frequency (VHF, 100 MHz) antennas. Our VHF data-processing flow is similar to that described for the HF radar and includes bandpass filtering, correction for antenna separation, geolocation from contemporaneous kinematic GNSS data, interpolation to standard trace spacing and exponential range gain to increase the relative amplitude of deeper reflections. Englacial reflectors were traced using the semi-automated reflection digitization routine in ImpDAR. A distinct high-amplitude reflector at ~75 m depth (Fig. 2b), observed in all the shallow RES profiles, is used to infer the spatial pattern of accumulation at West Dome.
We convert the TWTT to depth with a depth–density profile based on the firn core density measurements; we fit a Herron and Langway (Reference Herron and Langway1980) model to the firn core measurements to produce a smooth depth–density profile. The bright layer is at 73.5 m depth at the closest intersection (within 1 km) of our VHF profiles and the US-ITASE 02-4 core (Dixon and others, Reference Dixon2013; Fig. 1). The core extends to 72 m depth and was collected 17 years prior to our VHF survey. The age at 73.5 m depth is 1600 CE, or 420 years before the VHF radar data were collected, where we have assumed a steady-state firn column and determined the additional age in the bottom 1.5 m from the same model fit used for the TWTT conversion.
For the x133 and x135 radar lines, the conversion of TWTT for the bright layer to depth (Fig. 2c) and the inference of accumulation rate from the layer depth both depend on the density profile. Because there is no firn core at West Dome, we use a Herron and Langway (Reference Herron and Langway1980) firn model to estimate the depth–density profile. We use fixed inputs of 370 kg m−3 for the surface density and −40°C for the surface temperature; we start with an initial value of 0.13 m a−1 for the accumulation rate (all accumulation rates are given in ice equivalent and be converted to a kg m−2 a−2 using a density of 917 kg m−3), which is greater than the value inferred from the ITASE core because the average bright layer depth at West Dome is deeper than at the ITASE core site. We then infer the accumulation rate from the layer depth and the modeled density profile. To find a self-consistent solution, we iterate by calculating a new density profile along the radar line with the inferred accumulation rate at each trace, convert TWTT to depth and infer the updated accumulation rate with the new density profile. The solution becomes self-consistent within three iterations.
2.3 Automated phase-sensitive RES (ApRES)
Acquisitions of autonomous ApRES (Nicholls and others, Reference Nicholls2015) were performed at four locations along the x135 line in the 2018/19 season, and repeated in 2019/20. ApRES measurements are used to determine the relative phase offset of englacial reflectors between the two acquisitions. ApRES processing is implemented in ImpDAR (Lilien and others, Reference Lilien, Hills, Driscol, Jacobel and Christianson2020) and follows the processing flow described by Brennan and others (Reference Brennan, Lok, Nicholls and Corr2014). The wave speed is not corrected for the density variations of the firn and we do not interpret the data in the upper 150 m. The phase offset is converted to a relative range offset which determines a depth–velocity averaged over the period between the two measurements. Since phase is relative, the phase offset profile is unwrapped against a point of known velocity such as the bed or surface which then can shift the velocity profile by a constant.
2.4 Kinematic model setup
We employ a kinematic model to help interpret the internal structure of West Dome. In this section we describe the model setup and the application of the model follows in Section 4. The kinematic model uses inputs of the accumulation rate, bed topography and prescribed shape functions for the horizontal and vertical velocities (Nereson and others, Reference Nereson and Waddington2002; Waddington and others, Reference Waddington, Neumann, Koutnik, Marshall and Morse2007; Koutnik and others, Reference Koutnik and Waddington2012). These model inputs are based on the geophysical observations described in Section 3. The ability to fit the measured layers provides information about the stability of the current flow regime near the divide and about the appropriate forcing for depth–age modeling.
We choose a kinematic model for two primary reasons. First, the existing radar coverage is limited. The primary radar line is 10 km in length, which results in a distance of only ~3 ice thicknesses between the divide and downstream boundary. Using dynamic models requires extending the lower boundary to ~10 ice thicknesses from the divide, so that the area of interest is not affected by the edge boundary (Martín and others, Reference Martín, Hindmarsh and Navarro2006). Second, dynamic models have difficulty reproducing divide flow, particularly in areas with variations in bed topography (Goel and others, Reference Goel, Martin and Matsuoka2018). Therefore, we use a kinematic model where we can specify the shape of velocities and the position of the divide.
The horizontal velocities, u, are calculated as:
(1)$$u( {x, \;\hat{z}, \;t} ) = \bar{u}( {x, \;t} ) \phi ( {x, \;\hat{z}, \;t} ) $$
where $\bar{u}$ is the depth-averaged velocity, x is the distance along the flowline, $\hat{z}$ is the relative height above the bed, t is the time and ϕ is the horizontal velocity shape function. We use the Lliboutry (Reference Lliboutry1979) approximation:
(2)$$\phi ( {\hat{z}} ) = \displaystyle{{\,p + 2} \over {\,p + 1}} \Big ( {1-{( 1-\hat{z}) }^{\,p + 1}} \Big) $$
where p is the shape factor. The corresponding vertical velocity shape function is:
(3)$$\psi ( {\hat{z}} ) = \left({1-\displaystyle{{\,p + 2} \over {\,p + 1}}( {1-\hat{z}} ) + \displaystyle{1 \over {\,p + 1}}{( 1-\hat{z}) }^{\,p + 2}} \right).$$
The depth-averaged horizontal velocity is calculated from the ice flux from upstream, based on the accumulation rate and ice thickness. We assume no divergence in the ice flow (i.e. the flowband width is constant) based on comparisons with the measured surface velocities along x135 (Fig. S1). The vertical velocity is calculated as
(4)$$w = ( {\dot{b}-\dot{S}} ) \psi -\bar{u}\phi \left({\displaystyle{{\partial B} \over {\partial x}} + \hat{z}\displaystyle{{\partial H} \over {\partial x}}} \right)-\bar{u}H\mathop \int \limits_0^{\hat{z}} \displaystyle{{\partial \phi } \over {\partial x}}\partial \hat{\zeta }$$
where w, ψ, ϕ vary in x,$\;\hat{z}$ and t; $\dot{b}, \;\dot{S}, \;\;\;\bar{u}$ vary in x and t; B and H vary in x and $\dot{b}$ is the accumulation rate, $\dot{S}$ is the change in surface elevation, B is the bed elevation and H is the ice thickness.
An example of the modeled velocity field is shown in Figure 3. The details of our model results are presented in Section 4, but we present the individual velocity components, the three terms on the right-hand side of Eqn (4), to aid in visualizing the model construction. By representing the velocity field with these three terms, we can separate the influence of the bed from dynamics within the ice column to better understand the cause of variability in ice-flow along profile.
Fig. 3. (a) Modeled vertical velocity field using Eqn (4) with the divide at 4.7 km (black triangle). (b) First term of Eqn (4), which has the 1-D effects and is set by the accumulation rate and shape function. (c) Second term of Eqn (4), which expresses the vertical velocity for bed-parallel horizontal flow. (d) Third term of Eqn (4), which expresses the impact of the varying horizontal shape function. All panels use the same color scale to show the relative importance of the terms but have contours of different values to emphasize the magnitude of the terms.
We refer to the first term on the right-hand of Eqn (4) side as the 1-D term, which is governed by the surface vertical velocity. The surface vertical velocity (w s) is given by the difference between the accumulation rate ($\dot{b}$ ) and the change in the surface elevation ($\dot{S}$ ): $w_{\rm s} = \dot{b}-\dot{S}$ . We do not consider basal melt (see Sections 3.4 and 4.4).
We refer to the second term on the right-hand side of Eqn (4) as the 2-D term since it ensures bed and surface parallel flow.
The third term on the right-hand side in Eqn (4) accounts for variations in the shape of the horizontal vertical velocity, which in our case is the transition from divide to flank flow.
The bed elevation and ice thickness (Section 3.2) as well as the accumulation rate (Section 3.3) used in the kinematic model are derived from the radar data. The prescribed shape functions are based on vertical velocity profiles fit to the ApRES data (see Section 3.4). The ApRES measurements are at only four locations and we must define the vertical velocity continuously along the flowline. To do so, we define one shape function for the divide and one for the flank, and prescribe a smooth transition between them. We follow Nereson and Waddington (Reference Nereson and Waddington2002) and assume the width of divide flow influence can be approximated as:
(5)$$\beta ( x ) = \exp \left({-\displaystyle{{x^2} \over {2\sigma^2}}} \right)$$
where σ is the half-width of the divide zone for the vertical velocity and the transition for the horizontal velocity is then:
(6)$$\alpha ( x ) = \displaystyle{{\sqrt {2\pi } \sigma } \over {2x}}{\rm erf}\left({\displaystyle{x \over {\sqrt 2 \sigma }}} \right)$$
The transition in the vertical velocity occurs in 1–2 ice thicknesses from the divide, based on modeling studies (Raymond, Reference Raymond1983) and observations of the width of Raymond arches (e.g. Matsuoka and others, Reference Matsuoka2015). The transition in horizontal velocity occurs over a longer distance. The value of p for the shape functions (Eqn (3)) are 0.3 for the divide and 7 for the flank. The divide-zone half-width is 2800 m, ~1.5 ice thicknesses.
3. Geophysical observations
The geophysical observations reveal the internal structure and bed of West Dome and provide inputs and constraints for the kinematic model. We first determine the internal stratigraphy from the HF radar reflections. Second, we use the HF radar to determine the bed elevation and assess ambiguities in the bed picks using the kinematic model and measured internal layers. Third, we find the accumulation rate across the dome using a bright layer in the VHF radar which can be dated to 420 years before 2020, when the data were collected, near the US-ITASE shallow core. Finally, we assess the englacial vertical velocity profiles and identify a change in shape with distance from the present ice divide.
3.1 Internal layer geometry
We traced 22 internal layers from ~400 m below the surface to 200 m above the bed (Figs 2d, e). Layers are consistent and traceable along the full profile to within 400 m of the bed. Below 400 m, there is a difference in the character of the layers on the north and south sides of the divide. To the north, there is a bed-conformal and traceable layer between 100 and 200 m above the bed that extends to the northern end of the profile. Between that layer and the ice bottom, there appears to be coherent layering just north of the divide, but that signal is lost as the ice thickens along the profile. The lack of continuous layering in this interval is not surprising given the challenges of resolving deep layers (e.g. Drews and others, Reference Drews2009): returns from the deep layers are weaker due to greater attenuation; the layers vary more in depth as they drape the irregular bedrock and the relatively long wavelength (~55 m) of the HF system integrates a larger number of small-scale electrical contrasts which can cause more complex reflections. To the south, there is a package of layers that diverge from the traced layer at 400 m, indicating a local thickening of the layer package that sits between 200 and 400 m above the bed at the divide as you move south along the profile. These layers lose coherence as they approach the ice bottom.
The internal layers cannot be dated with existing information from ice cores because there are no continuous radar profiles between the two West Dome radar profiles and the ITASE ground-based or airborne radar datasets. The only layer that can be dated is the bright reflector in the VHF data (Section 2.2, Fig. 2b) because it can be unambiguously matched by correlation to a radar line within 1 km of the dated ITASE core despite the lack of a continuous link. The internal reflectors of the HF radar profiles at West Dome are not distinct enough to be matched to the East Dome and South Ridge HF radar profiles.
3.2 Bed topography
We imaged the bed of West Dome along the x135 and x133 radar lines using the HF radar. The current surface divide is located above a bedrock high (Fig. 2). The ice thickness along x135 ranges from 1510 m below the divide summit to 1750 m at ~8 km along the radar line. The bed topography is smoother than at East Dome (Jacobel and others, Reference Jacobel and Welch2005). The bed reflection is visible along the full profile (Fig. 4), although there are two locations where the bed reflection is ambiguous. The most important of these is at 6 km along the radar line. The brightest reflection is at 1640 m depth; however, there is also a weaker return from 1580 m depth. This smaller amplitude reflection is likely an off-nadir return and may indicate a decameter-scale bump at the bed close to the profile.
Fig. 4. Interpreted bed location is shown (red). There is an ambiguity in the bed reflection at 6 km (arrow).
To evaluate the significance of such a feature on the internal layering, we run the kinematic model with bed topographies based on both the original pick and by picking the decameter-scale bump at ~6 km. We also use two different low-pass filters of 500 and 1000 m, resulting in four different potential boundary conditions. The observed layers (and the associated model-predicted layers) deeper than 1000 m are shown in Figure 5 because they are the most affected by the bed boundary condition. The age of the internal layers is unknown so we cannot compare the modeled and measured layers based on age; instead, we compare the observed layers to a modeled layer with the same average depth. The modeled layers produced using the 500 m filtered beds show greater variability than the observed layers while the 1000 m filtering produces internal layering with smoothness comparable to the measured layers. Using the bed with 1000 m filtering and the decameter-scale bump, the model is able to reproduce the slight dip in the deepest layer at 5.5 km along the radar line. Therefore, we use the 1000 m filtered bed with the bump in all of the model runs used to help interpret the internal layering (Section 4).
Fig. 5. Black lines are observed internal layers. The bottom layers are the four different beds used in the kinematic model. The two brown lines have no bump at 6 km; the two blue lines have the bump added. Low-pass filtering of 500 m is shown as solid lines and 1000 m is shown as dash dot. The modeled internal layers have line color and type that match the bed.
3.3 Accumulation rate
The accumulation rates inferred from the VHF radar profiles along x135 and x133 are shown in Figure 6. There is a ~2 km gap in the x135 line because of poor radar quality. The inferred accumulation rates are similar between the two lines. The highest accumulation rates of 0.14 m a−1 are in the north and are relatively constant from 7 to 10 km along the radar lines. The accumulation rates decrease by ~25% to <0.115 m a−1 at the south end of the radar line. Because the data from x133 is of higher quality and agrees well with the high-quality portions of x135, we use the accumulation rate inferred from x133 as the forcing for the model runs used to help interpret the internal layering (Section 4).
Fig. 6. Accumulation rate inference from the bright layer dated to an age of 420 years (1600 CE). Profile x135 has poor data quality in the center so x133 is used for model forcing given the close agreement between the two. Note that the x133 profile is slightly offset at 5 km (Fig. 1, inset) such that we interpolate linearly from 4.8 to 5.1 km.
3.4 Englacial vertical velocities
Divides with frozen beds have vertical strain rates which are greater near the surface and smaller near the bed compared with locations of larger surface slopes and shear strain rates. We refer to this characteristic shape of the vertical velocity near divides as divide-flow, and the more linear vertical velocity away from divides as flank-flow. The difference in the vertical velocity patterns yields Raymond arches (Raymond, Reference Raymond1983) when the divide position is stable. Englacial vertical velocities were derived using ApRES acquisitions in the 2018/19 and 2019/20 field seasons at four sites shown in Figure 2a with the goal to (1) identify whether the characteristic shape of the vertical velocity is characteristic of those with a frozen bed and (2) provide constraints for kinematic modeling of internal layers.
The vertical velocity profiles are plotted by relative height above the bed in Figure 7a (where 0 represents the bed and 1 the ice surface) and by absolute height above the bed in Figure 7b. A useful way to visualize whether divide flow is occurring is to perform a linear fit of the upper portion of the vertical velocities and to compare the extension of the linear fit to the bed with the measured vertical velocities (Fig. 7a). The fit is calculated by linear regression between a relative height of 0.5–0.9. The near-surface portion is excluded because of the influence of firn compaction on the velocity profile above 0.9. For n5000, ~5 km north of the divide, the linear fit matches well for all velocities except for those below 0.05. This near-linear increase of vertical velocity with height above bed is typical of ice away from a divide (Kingslake and others, Reference Kingslake2014). At sites n0 and n1000, there is a noticeable divergence from the linear fit, as the vertical velocities show little change below 0.3H.
Fig. 7. Modeled and measured vertical velocity profiles for the four ApRES locations labeled by their approximate distance in m from the divide. The horizontal spacing of the vertical lines in both panels is 0.05 m a−1. The positions are shown in Figures 1 and 2, with n0 in yellow, n1000 in orange, n2000 in pink and n5000 in purple. (a) ApRES vertical velocities; black lines are a linear fit to the interval of 0.5H to 0.9H. (b) Lliboutry fits to the vertical velocities where the surface vertical velocity (excluding firn compaction) is required to be within 0.02 m a−1 of the 420 year average accumulation rate, the vertical component of ice flow over uneven bedrock is included and a constant shift for the vertical velocity profile is a free parameter. The transition from divide-flow with low p values to flank-flow with high p values occurs with increasing distance from the divide (see Section 2.4).
The four sites were revisited in the 2021/22 field season. Sites n0, n1000 and n5000 showed nearly identical vertical velocities for the 19/20–21/22 period as the 18/19–19/20 period; however, n2000 did not (Figs S2–S5). Therefore, we focus our discussion on the differences between n0 and n1000 near the present divide and n5000, which is outside the divide region. A discussion of the vertical velocities using the 2021/22 field data are in given in the Supplementary material.
To quantitatively describe the vertical velocity profiles, we fit the vertical velocities using two simplified versions of Eqn (4). First, we follow Kingslake and others (Reference Kingslake2014) and fit with the following 1-D term,
(7)$$w( {\hat{z}} ) = w_{\rm s}\psi{\hat{z}} )$$
The value of p in the calculation of ψ (Eqn (4)) is a useful measure of the shape of the vertical velocity and can be used to characterize the influence of the divide. Small values of p indicate divide-like flow, with increased vertical strain near the surface and reduced vertical strain near the bed. Higher p values indicate more uniform vertical strain with depth.
Since the ApRES measures phase, we can only infer a relative velocity, and we must specify one absolute velocity. We assume no motion at the bed and set the average of the velocities within 100 m of the bed to zero. This allows an estimate of the surface vertical velocity, excluding firn compaction, which should be similar to the ice-equivalent accumulation rate. The best-fit values of p and w s are given in Table 1; the velocities within 100 m of the bed and within 150 m of the surface are not used in the fit. The inferred surface vertical velocity exceeds the inferred 420 year average accumulation rates from the VHF radar data, except at n2000. The inferred surface vertical velocity at n5000 is 0.17 m a−1, nearly 0.04 m a−1 (30%) greater than the 420 year average accumulation rate of 0.13 m a−1. We suspect that this is because n5000 has the greatest horizontal velocity (~0.5 m a−1) and is expected to have the largest component of the vertical velocity due to horizontal flow over the bedrock.
Table 1. Fits to ApRES vertical velocities
Bold values are the parameters inferred.
To account for ice flow over an irregular bed, we also fit the vertical velocities with
(8)$$w( {\hat{z}} ) = w_{\rm s}\psi + \bar{u}\phi \displaystyle{{\partial B} \over {\partial x}}$$
We neglect the component of the vertical velocity due to lateral variation in the shape of the horizontal velocity profile because it is the smallest term (Fig. 3). We estimate the depth-averaged horizontal velocity from the balance velocity. The bed slopes are calculated on a 200 m length scale using the unfiltered bed (Fig. S6).
In the fit using Eqn (8), the surface vertical velocity is constrained to within 0.02 m a−1 of the measured accumulation rate. In addition, the fit finds a uniform shift of the vertical velocity profile since the velocities are relative; this relaxes the constraint we applied when using Eqn (7), where the average velocity in the bottom 100 m was set to 0, and allows us to determine if the vertical velocities are consistent with the glaciological setting. The fit values are also shown in Table 1.
The largest difference in the fit parameters between Eqns (7) and (8) is at n5000, where a reasonable fit with a surface vertical velocity of 0.151 m a−1 is possible because of the negative (upward) vertical velocities near the bed caused by flow up the bedrock slope (e.g. Hills and others, Reference Hills2022). The two measurements closest to the divide, n0 and n1000 are similar to each other with the value of p being indicative of divide-flow; n5000 has a large p value indicative of flank-flow. Because of the uncertainty in the measurements at n2000 and the lack of measurements in between n2000 and n5000, it is difficult to determine the width of the divide zone.
4. Kinematic model results
Our kinematic model and data comparison focuses on five objectives: (1) fit the shallowest layer using the kinematic model to determine how well the modern observations reflect late Holocene dynamics and forcing; (2) assess differences in layer depths in the divide region to evaluate the presence of a Raymond arch; (3) fit all layers across the full radar profile to evaluate the likelihood of divide migration within our domain; (4) determine the likelihood of basal melt and (5) use our best fit model to produce a depth–age scale that captures the annual-layer thickness near the bed.
4.1 Fit to shallowest layer
Before examining the model fit for all of the internal layers, we focus on the shallowest HF radar layer which has a depth of ~430 m and a modeled age of ~4 ka. We assess whether the modern observations of accumulation and vertical velocity are likely to have been constant for the past ~4 ka. The accumulation rate pattern is a 420 year average while the measured vertical velocities are a snapshot over 1 year.
Using the modern patterns of accumulation and vertical velocity for the past 4 ka produces a layer which is too shallow beneath the divide and too deep north of the divide (Fig. 8). This suggests that the modern accumulation and vertical velocity patterns have not been constant for the past 4 ka, which may indicate that the divide has not been stable at its current position. To evaluate this possibility, we test scenarios of divide onset occurring in 1 ka increments. Prior to the divide onset, we assume spatially uniform flank-flow vertical velocities and uniform accumulation. The modeled layers and the misfit are shown in Figure 8. The modeled layers fit the measured layer best with a 2 ka onset of modern conditions although there is still a significant discrepancy on the south edge of the profile.
Fig. 8. Model fit to the shallowest HF (deep) radar layer. The legend gives the spatially averaged misfit to the measured layer (in m) and the age at which the modern vertical velocity and accumulation rate patterns are applied (in ka). The divide zone has a width of 1.5H. The modeled age of the layer is ~4 ka.
For these tests, we have assumed that the divide in the surface topography causes both the vertical velocity pattern and the accumulation pattern. The vertical velocity pattern is primarily caused by the low deviatoric stress beneath the divide although the development of ice fabric may also play a role (Pettit and others, Reference Pettit, Thorsteinsson, Jacobson and Waddington2007; Martín and others, Reference Martín, Gudmundsson, Pritchard and Gagliardini2009); however, it is less clear how closely the surface topography of the divide and the accumulation pattern are coupled. More complicated temporal and spatial variations in the accumulation pattern and divide position can improve the fit; however, without evidence to guide a priori expectations for how the patterns should develop in time, the improvements would not necessarily be meaningful. An improved fit with more complicated variations in the forcings would reinforce the primary conclusion that the modern conditions are unlikely to have been stable for the past ~4 ka. In the following sections, we explore the fit to all of the internal layers to extend our understanding of the divide position and accumulation pattern to older ages.
4.2 Evaluating the presence of a Raymond arch
Although the ApRES vertical velocities indicate that divide-flow is occurring, there is not an obvious Raymond arch beneath the divide. However the measured internal layers are shallower in the ice column at ~6 km than might be expected based on the bed topography (Fig. 2d). In areas of flat bedrock, the amplitude of a Raymond arch can be determined by comparing the depth of an internal layer beneath the divide to the depth on the flank outside of the divide region (e.g. Roosevelt Island, Conway and others, Reference Conway, Hall, Denton, Gades and Waddington1999), but the identification of Raymond arches is challenging in areas with greater ice thickness, lower accumulation rates and without smooth bedrock (e.g. Goel and others, Reference Goel, Martin and Matsuoka2018). Since the irregular bedrock at our site does not allow a simple comparison between the depth of a layer beneath the divide and on the flank, we use the kinematic model.
We compare the height of the measured layer to one from the kinematic model forced with a spatially constant vertical-velocity shape function (p = 4, Fig. 9). The difference between the measured and modeled layers at 5.5 km along the radar line has both similarities and differences to the Raymond arch amplitude at Roosevelt Island (Fig. 9b). The Raymond arch at Roosevelt Island began forming ~3 ka (Conway and others, Reference Conway, Hall, Denton, Gades and Waddington1999; Martín and others, Reference Martín, Hindmarsh and Navarro2006). The onset of a layer mismatch begins deeper in the ice column, at a normalized height of ~0.7 at West Dome compared to 0.9 at Roosevelt Island. There are also large misfits in the layer depths outside of the divide zone at West Dome, particularly to the south (Fig. 9a). Thus, the existence of a Raymond arch is uncertain. In Section 4.3, we model the full set of internal layers to better understand the observed layering.
Fig. 9. (a) Modeled internal layers (blue) assuming uniform flank flow (p = 4) compared to measured layers (black). The bed is shown in brown. (b) The difference between modeled and measured layers (blue circles) at 5.75 km (black vertical line) compared to the Raymond arch amplitude of Roosevelt Island (red squares) (Conway and others, Reference Conway, Hall, Denton, Gades and Waddington1999).
4.3 Model fit to all layers
The kinematic model fit to the shallowest layer (Section 4.1) suggested that the modern divide and accumulation pattern have not been stable for more than a couple of thousands of years while the possibility of a Raymond arch exists (Section 4.2) offset to north from the present divide. Here we evaluate the fit to all the layers for test scenarios of divide migration and variations in the accumulation pattern. The primary goals are to evaluate the likelihood of divide migration within our domain and identify the appropriate forcings for the thermal and depth–age modeling.
The kinematic model fit to the full depth of measured layers is shown in Figure 10 for two scenarios: (A) using the measured accumulation and vertical velocity forcings for all ages and (B) using the measured forcings for only the past 2 ka (i.e. divide onset at 2 ka, see Section 4.1) and spatially uniform patterns of accumulation and vertical velocity for ages older than 2 ka. The RMSE, calculated as fraction of layer depth offset, is improved by a factor of ~2 when the modern patterns are used for only the past the 2 ka. The resulting misfit shows a pattern of the modeled layers being too deep to the north of the divide, and too shallow to south.
Fig. 10. (a, b) Bed (brown), measured internal layers (black) and modeled internal layering (blue) using the modern spatial patterns of vertical velocity and accumulation rate for all ages (a) and only the past 2 ka (b, uniform vertical velocity and accumulation rate before 2 ka). (c, d) Relative misfit (measured depth – modeled depth)/measured depth. Panels (a) and (b) show only 13 of the 22 measured layers for clarity while panels (c) and (d) show all 22 layers.
There are two primary changes in the model forcing that will yield deeper layers to the south, shallower layers to the north and thereby improve the fit: (1) move the region of divide flow to the north and (2) reverse the accumulation gradient such that there is more accumulation to the south and less to the north.
We performed model runs with the divide shifted north to 6 km along the flowline for different ages prior to 2 ka. The divide then migrates from 6 to 4.7 km (the modern divide position) between 2.5 and 2 ka. The best fit occurs with a divide onset at 6 km at 6 ka. In Figure 11, we show the misfit for divide onset at 6 ka as well as for 9 ka; the misfit for the base scenario (Fig. 10b) of divide onset at 2 ka at 4.7 km along the flowline is also shown. In the scenarios with divide onset prior to 2 ka, the modeled layers are shallower to the north because of the increased layer thinning in the upper portion of the ice sheet. The modeled layers are also deeper to the south, better matching the measured layers, because the average depth of the modeled layer is constrained to match the average depth of the measured layer. Divide onset at 6 ka results in better fits from 5 to 7 km and from 0 to 4 km along the flowline throughout most of the depth; however, the misfit increases from 4 to 5 km, the modern divide region. The scenario with divide onset at 9 ka further increases the misfit beneath the modern divide, which drives an overall increase in the misfit for the full flowline. The timing of the divide onset should be interpreted with caution because the layers are not dated and the improvement in fit is modest. These tests with the kinematic model indicate that the divide is unlikely to have been established prior to the mid-Holocene.
Fig. 11. Relative misfit as calculated in Figures 10e, f for kinematic model scenarios with no divide onset prior to 2 ka (a), divide onset at 6 ka at 6 km (b) and divide onset at 9 ka at 6 km (c).
To assess the impact of the spatial pattern of accumulation, we performed two tests to compare with the base scenario in Figure 10b. In both scenarios, the modern pattern of vertical velocities began at 2 ka with uniform flank-flow for older ages. The first scenario kept the modern accumulation gradient fixed for all ages; the second scenario flipped the modern accumulation gradient for ages older than 2 ka (Fig. S7). Using the modern accumulation gradient for all ages increases the relative misfit from 2.1 to 3.2%, while using a flipped accumulation gradient increases the relative misfit from 2.1 to 2.2%. As indicated above, a reversal in the accumulation gradient results in deeper layers to the south and shallower layers to the north, yielding a similar misfit as the base scenario. Different accumulation gradients and/or timing of changes in the past accumulation gradients can yield better fits. This issue of non-uniqueness is inherent in inferring accumulation rates from deep internal layers near ice divides (e.g. Koutnik and others, Reference Koutnik2016), and is more challenging when the layers have not been dated. Scenarios with an onset of the flipped accumulation gradient in the Holocene, and uniform accumulation for older ages, produce a worse fit than when the flipped gradient is applied for all ages. We have no a priori expectation that a reversal in the gradient should have occurred during the early-to-mid-Holocene when West Antarctic climate was relatively stable. While past changes in the accumulation gradient cannot be well constrained, it is unlikely that the modern pattern has persisted for more than a few thousand years.
Overall, the kinematic modeling shows that the divide position was unlikely to be stable before the mid-Holocene and that the pattern of accumulation is difficult to constrain. A key conclusion for the depth–age relationship is that because the divide is unlikely to have been stable, a flank-flow vertical velocity is most appropriate for modeling the depth of the older ice.
4.4 Temperature modeling and basal thermal state
We use a 1-D ice-and-heat flow model to estimate the internal and basal temperature at West Dome. The model is described in Fudge and others (Reference Fudge, Biyani, Clemens-Sewall and Hawley2019) and Buizert and others (Reference Buizert2021). We take two approaches. First, we estimate the internal temperature based on the 65 mW m−2 average of three recent remotely sensed estimates of geothermal flux which are tightly clustered: 62 ± 10 mW m−2 (Shen and others, Reference Shen, Wiens, Lloyd and Nyblade2020), 68 ± 25 mW m−2 (Stal and others, Reference Stal, Reading, Halpin and Whittaker2021) and 66 ± 11 mW m−2 (Martos and others, Reference Martos2017). This provides a prior estimate for the temperature profile for use in planning ice core drilling and analyses. Second, we find the minimum geothermal flux at which melting begins (e.g. Fudge and others, Reference Fudge, Biyani, Clemens-Sewall and Hawley2019).
The ice-and-heat flow model is forced by time-varying surface temperature, accumulation rate and vertical velocity shape factor histories. We use an ice thickness of 1600 m, which is held constant. The histories are based on the two closest deep ice cores, the WAIS Divide ice core (WDC) and the South Pole ice core (SPICEcore); these histories are extended through the LIG using the EPICA Dome C core (EPICA, 2004; Veres and others, Reference Veres2013). The temperature histories are calculated using the below equation:
(9)$$T( t ) _{{\rm Model}} = T( 0 ) _{{\rm Herc}} + ( {T{( t ) }_{{\rm core}}-T{( 0 ) }_{{\rm core}}} ) $$
where T(t)Model is the surface temperature forcing for the thermo-mechanical model, T(0)Herc is the modern temperature at Hercules Dome and is −41°C based on a 10 m firn temperature measurement in 2019, T(t)core is the inferred temperature history from the one of the ice cores, T(0)core is the modern surface temperature of the ice core site and t is the model time. The Last Glacial Maximum (LGM)-modern change is −11°C for WDC (Cuffey and others, Reference Cuffey2016) and −6°C for SPICEcore (Kahle and others, Reference Kahle2021). The accumulation rate histories are calculated as
(10)$$\dot{b}( t ) _{{\rm model}} = \dot{b}( 0 ) _{{\rm Herc}} \times \displaystyle{{\dot{b}{( t ) }_{{\rm core}}} \over {\dot{b}{( 0 ) }_{{\rm core}}}}$$
Using the same notations as for the surface temperature in Eqn (9). The modern accumulation rate at West Hercules Dome, $\dot{b}( 0 ) _{{\rm Herc}}$ , is 0.13 m a−1. Both WDC (Fudge and others, Reference Fudge2016) and SPICEcore (Kahle and others, Reference Kahle2021) have inferred LGM accumulation rates of ~40% of modern.
The englacial temperature profiles are shown in Figure 12. The geothermal flux necessary for melting at present is 85 ± 7 mW m−2 using the SPICEcore scaling and 89 ± 7 mW m−2 using the WDC scaling. The uncertainty was calculated assuming an end-member modern surface temperature of −39°C (2°C warmer) and modern accumulation rate of 0.11 m a−1 (0.02 m a−1 lower). The minimum geothermal flux to produce melting is considerably larger than the remotely sensed estimates inferred by three different studies. Our estimated range of the geothermal flux necessary for basal melting overlaps with the uncertainty of only one of these estimates, Stal and others (Reference Stal, Reading, Halpin and Whittaker2021), and only because of its large (25 mW m−2) uncertainty.
Fig. 12. Temperature profiles for West Dome. Solid lines use the average of the remotely sensed geothermal flux estimates (65 mW m−2); dashed lines assume the bed is at the pressure melting point. The surface temperature and accumulation rate histories are scaled from values at SPICEcore (blue) and WAIS Divide (red) using modern values of −41°C and 0.13 m a−1.
4.5 Depth–age relationship
Because the measured internal layers cannot be dated, preliminary depth–age relationships must be assessed from modeling. We use the onset of the divide at 2 ka (Sections 4.1 and 4.3) and three different choices for the vertical velocity shape factor (4, 7 and 10) for ages prior to 2 ka. We use the two accumulation histories described in Section 4.4 which are based on WDC and SPICEcore. For the kinematic model forcing, the fractional accumulation is multiplied by the modern spatial accumulation-rate pattern.
The resulting depth–age relationships at 5.5 km along the x135 flowline are shown in Figure 13. The choice of vertical-velocity shape factor has a larger influence than the choice of accumulation history. Ice of 132 ka age, the start of the LIG, ranges from 47 to 87 m above the bed. The average LIG annual-layer thickness is ~1 mm.
Fig. 13. Depth–age relationships at 5.5 km along the flowline. Accumulation forcing is scaled to SPC (SPICEcore, solid) or WDC (WAIS Divide, dashed). Three different vertical velocity profiles for flank flow (p = 4, 7 or 10) are shown for each accumulation forcing. The height above the bed for ice of 132 ka age ranges from 47 to 87 m.
5.1 Potential for a deep ice core
One goal of drilling a deep ice core at Hercules Dome is to recover ice that includes the onset of the LIG (132 ka), and preferable Termination II (the transition into the LIG from the previous glacial maximum) as well. To preserve ice of this age, there must be limited or zero basal melt and no more than minor disturbances to the stratigraphy. The current pattern of vertical velocity across the divide at West Dome suggests the ice is frozen to the bed. This is supported by the ice-and-heat flow modeling (Section 4.4). With a geothermal flux of 65 mW m−2, the ice at the bed is ~−10°C, which would allow subglacial bedrock drilling after completion of the ice core.
The presence of a frozen bed today also implies no significant impact from basal melt at older ages, even if the bed is near the pressure melting point. In the model runs where we found the maximum geothermal flux before basal melt occurred, the bed remains frozen at all ages with the WAIS Divide forcing. For the SPICEcore scenario, there is a ~3000 year period of basal melt that occurs at the end of the LGM (i.e. from ~25 to 22 ka). The maximum melt rate is 0.04 mm a−1, and the total predicted melt is 80 mm. While these two scenarios do not encompass the full range of possibilities, they indicate that significant basal melt is unlikely to impact the age of the basal ice.
The impacts of changes in ice thickness are difficult to evaluate because it is not clear if Hercules Dome was thicker or thinner at the LGM. A thicker ice sheet at the LGM would cause a warmer basal temperature and likely increased basal melt, while a thinner ice sheet would have the opposite effect. The modest amount of basal melt would not remove ice from the LIG. In fact, basal melt rates up to ~0.8 mm a−1 would preserve thicker packets of ice from the LIG, although the ice would be closer to the bed. While the evidence for a frozen bed and limited basal melt may not be as definitive as a measured basal temperature, our study suggests that basal melt is unlikely to have removed ice older than 132 ka from the West Dome site.
The preservation of stratigraphic order is more difficult to assess. The layering at West Dome is resolvable until just a few tens of meters above the bed, but the deepest layer that is traceable along the full radar profile is ~400 m above the bed. Ice from the start of the LIG (132 ka) is likely between 47 and 87 m above the bed, and thus continuity of the deep stratigraphy is important for recovering ice of this age. The deep internal reflectors are coherent all the way to the ice bottom just north of the divide, with more bed conformal layering on the north side of the divide generally, suggesting the north is the preferable location for a core.
The internal layers suggest that the current divide position has not been stable beyond the mid-Holocene; however, the internal layers alone are not enough to constrain the extent of divide migration. Modest variations in the divide position of a few kilometers are enough to prevent the development of a clear Raymond arch above an irregular bed. The lack of a clear Raymond arch at West Dome is not surprising because the characteristic time of ~13 ka, ($H/\dot{b}$ ), is substantially larger than for the thinner and higher accumulation coastal domes where Raymond arches are widely observed (Matsuoka and others, Reference Matsuoka2015).
Our best-fit models indicate that average annual layer thicknesses are ~1 mm during the LIG. With continuous-flow analysis methods that are able to make measurements with ~0.5 cm resolution (e.g. Röthlisberger and others, Reference Röthlisberger2000; Osterberg and others, Reference Osterberg, Handley, Sneed, Mayewski and Kreutz2006; Gkinis and others, Reference Gkinis, Popp, Johnsen and Blunier2010; Sigl and others, Reference Sigl2016; Jones and others, Reference Jones2017), measurements of the chemistry and isotopes during the interglacial will have approximately decadal resolution. For gas sampling with a typical sample length of ~10 cm, this translates to a maximum of centennial resolution. The largest challenge from thin annual layering will not be resolution but the limited ice availability for ice of older ages. One approach to overcome the challenge of a thinner packet of LIG ice is to drill one or multiple replicate cores with a deviation from the main borehole.
One limitation of the shallower ice at West Dome is that borehole thermometry (e.g. Cuffey and others, Reference Cuffey2016) to constrain the magnitude of the glacial–interglacial temperature change becomes increasingly difficult. However, new techniques such as water isotope diffusion lengths, and firn thickness and gas-age ice-age difference (Δage) constraints, have been pioneered to provide alternate ways of addressing this question (e.g. Gkinis and others, Reference Gkinis, Simonsen, Buchardt, White and Vinther2014; Holme and others, Reference Holme, Gkinis and Vinther2018; Kahle and others, Reference Kahle2021; Buizert and others, Reference Buizert2021).
The geophysical survey at Hercules Dome indicates that a continuous climate record through the LIG could likely be obtained at West Dome. West Dome is the highest elevation portion of Hercules Dome, and is the most pronounced dome with the simplest ice-flow setting. The ice thickness varies relatively little at West Dome, with a mean depth of ~1600 m, which minimizes the risk of disrupted layering from complicated ice flow. The ice temperature at the bed is likely below freezing, preserving old ice near the bed.
5.2 Future plans
A full-field season focused on the West Dome region is planned for 22/23. Future measurements will include more extensive HF and VHF surveys. Additional ApRES measurements will be made to both more fully characterize the englacial vertical velocities and to constrain ice fabric variations. These future measurements will improve both the spatial resolution and coverage of the data to more precisely define the best ice core site within the West Dome region. Other locations at Hercules Dome are less promising for a deep ice core, and are not prioritized for future research. At East Dome, where the US-ITASE traverse visited (Jacobel and others, Reference Jacobel and Welch2005), the basal topography and ice flow appears prohibitively complex. South Ridge has moderately thicker ice but the ice flow is from East Dome and over complex bedrock topography.
The surface elevation of Hercules Dome is not defined by a single topographic high, but rather three features with local divides. We have focused on West Dome, which is at the highest elevation in the Hercules Dome region and has the simplest ice-flow regime. Coherent internal layers are imaged to within tens of meters of the bed. The ice is of moderate thickness (~1600 m) and the bed is likely frozen; a geothermal flux of ~20 mW m−2 greater than suggested by regional estimates would be required to produce melting. The internal layer pattern cannot be well fit with a kinematic ice-flow model if modern accumulation and vertical velocity patterns are held constant in the past; the lack of fit indicates that the modern divide position has likely been established during the Holocene. As with other inland divides in Greenland and Antarctica, modest variations in divide position would prevent the development of a prominent Raymond arch. Our best-fit modeled depth–age relationship suggests ice from the LIG is between 45 and 90 m above the bed, and that the average annual layer thickness for that time period is ~1 mm. The modeled depth–age relationship is insensitive to the Holocene flow history because the age of the deep ice depends primarily on the average shape of the vertical velocity profile for the Last Glacial–Interglacial cycle. The lack of modern basal melt implies no significant melt in the past and therefore ice from the LIG should not have been removed. In conclusion, West Dome is a promising location for a deep ice core that will preserve a climate record sensitive to the size of WAIS during the LIG.
The supplementary material for this article can be found at https://doi.org/10.1017/jog.2022.80.
Radar data are available at https://www.usap-dc.org/view/project/p0010359.
This study was supported by NSF-OPP awards 1744649 and 1841844. We acknowledge the Antarctic Support Contractor and the Air National Guard for logistical support. We also thank Erika Schreiber for conducting ApRES measurements in the 2021/22 season.
TJF, KC and EJS conceived of the study. AH, ANH, BHH, EJS, GKO, JEC, KC and NH collected the geophysical data. AH, ANH, BHH, KC, LD, NH and TJF contributed to the data processing, analysis and interpretation. BHH and TJF conducted ApRES interpretation in this manuscript. ANH, LD and TJF inferred the snow-accumulation record for this manuscript. NH contributed insight to layer and bed interpretations. TJF performed the kinematic and temperature modeling. All authors contributed to the manuscript.
Bamber, JL, Riva, REM, Vermeersen, BLA and LeBrocq, AM (2009) Reassessment of the potential sea-level rise from a collapse of the West Antarctic ice sheet. Science 324, 900–903.CrossRefGoogle ScholarPubMed
Brennan, PV, Lok, LB, Nicholls, K and Corr, H (2014) Phase-sensitive FMCW radar system for high-precision Antarctic ice shelf profile monitoring. IET Radar, Sonar & Navigation 8, 776–786. doi: 10.1049/iet-rsn.2013.0053CrossRefGoogle Scholar
Brook, EJ and 8 others (2005) Timing of millennial-scale climate change at Siple Dome, West Antarctica, during the Last Glacial period. Quaternary Science Reviews 24(12–13), 1333–1343.CrossRefGoogle Scholar
Buizert, C and 16 others (2015) The WAIS Divide deep ice core WD2014 chronology – part 1: methane synchronization (68–31 ka BP) and the gas age-ice age difference. Climate of the Past 11(2), 153–173.CrossRefGoogle Scholar
Buizert, C and 40 others (2021) Antarctic surface temperature and elevation during the Last Glacial Maximum. Science 372(6546), 1097.CrossRefGoogle ScholarPubMed
Christianson, K and 6 others (2016) Basal conditions at the grounding zone of Whillans Ice Stream, West Antarctica, from ice-penetrating radar. Journal of Geophysical Research-Earth Surface 121(11), 1954–1983.CrossRefGoogle Scholar
Christianson, K, Jacobel, RW, Horgan, HJ, Anandakrishnan, S and Alley, RB (2012) Subglacial Lake Whillans – ice-penetrating radar and GPS observations of a shallow active reservoir beneath a West Antarctic ice stream. Earth and Planetary Science Letters 331–332, 237–245.CrossRefGoogle Scholar
Conway, H, Hall, BL, Denton, GH, Gades, AM and Waddington, ED (1999) Past and future grounding-line retreat of the West Antarctic ice sheet. Science 286(5438), 280–283.CrossRefGoogle ScholarPubMed
Cuffey, KM and 8 others (2016) Deglacial temperature history of West Antarctica. Proceedings of the National Academy of Sciences of the United States of America 113(50), 14249–14254.CrossRefGoogle ScholarPubMed
DeConto, RM and Pollard, D (2016) Contribution of Antarctica to past and future sea-level rise. Nature 531(7596), 591–597. doi: 10.1038/nature17145CrossRefGoogle ScholarPubMed
Dixon, DA and 6 others (2013) Variations in snow and firn chemistry along US ITASE traverses and the effect of surface glazing. The Cryosphere 7(2), 515–535.CrossRefGoogle Scholar
Drews, R and 7 others (2009) Layer disturbances and the radio-echo free zone in ice sheets. The Cryosphere 3, 195–203.CrossRefGoogle Scholar
Dutton, A and 8 others (2015) Sea-level rise due to polar ice-sheet mass loss during past warm periods. Science 349(6244), aaa4019, doi: 10.1126/science.aaa4019CrossRefGoogle ScholarPubMed
Dyer, B and 7 others (2021) Sea-level trends across The Bahamas constrain peak Last Interglacial ice melt. Proceedings of the National Academy of Sciences 118(13), 1–11. doi: 10.1073/pnas.2026839118CrossRefGoogle ScholarPubMed
EPICA Community Members (2004) Eight glacial cycles from an Antarctic ice core. Nature 429, 623–628.CrossRefGoogle Scholar
Fudge, TJ and 8 others (2016) Variable relationship between accumulation and temperature in West Antarctica for the past 31,000 years. Geophysical Research Letters 43(8), 3795–3803.CrossRefGoogle Scholar
Fudge, TJ, Biyani, SC, Clemens-Sewall, D and Hawley, RL (2019) Constraining geothermal flux at coastal domes of the Ross Ice Sheet, Antarctica. Geophysical Research Letters 46(22), 13090–13098.CrossRefGoogle Scholar
Garbe, Julius, Albrecht, Torsten, Levermann, Anders, Donges, Jonathan F and Winkelmann, Ricarda (2020) The hysteresis of the Antarctic Ice Sheet. Nature 585(7826), 538–544. doi: 10.1038/s41586-020-2727-5CrossRefGoogle ScholarPubMed
Gkinis, V, Popp, TJ, Johnsen, SJ and Blunier, T (2010) A continuous stream flash evaporator for the calibration of an IR cavity ring-down spectrometer for the isotopic analysis of water. Isotopes in Environmental and Health Studies 46(4), 463–475.CrossRefGoogle ScholarPubMed
Gkinis, V, Simonsen, SB, Buchardt, SL, White, JWC and Vinther, BM (2014) Water isotope diffusion rates from the NorthGRIP ice core for the last 16,000 years – glaciological paleoclimatic implications. Earth and Planetary Science Letters 405, 132–141.CrossRefGoogle Scholar
Goel, V, Martin, C and Matsuoka, K (2018) Ice-rise stratigraphy reveals changes in surface mass balance over the last millennia in Dronning Maud Land. Journal of Glaciology 64(248), 932–942.CrossRefGoogle Scholar
Hein, A and 9 others (2016) Evidence for the stability of the West Antarctic ice sheet divide for 1.4 million years. Nature Communications 7, 10325. doi: 10.1038/ncomms10325CrossRefGoogle ScholarPubMed
Herron, MM and Langway, CC (1980) Firn densification – an empirical model. Journal of Glaciology 25(93), 373–385.CrossRefGoogle Scholar
Hills, BH and 10 others (2022) Geophysics and thermodynamics at South Pole Lake indicate stability and a regionally thawed bed. Geophysical Research Letters 49, e2021GL096218. doi: 10.1029/2021GL096218CrossRefGoogle Scholar
Holloway, MD and 5 others (2016) Antarctic Last Interglacial isotope peak in response to sea ice retreat not ice-sheet collapse. Nature Communications 7, 1–11. doi: 10.1038/ncomms12293CrossRefGoogle Scholar
Holme, C, Gkinis, V and Vinther, BM (2018) Molecular diffusion of stable water isotopes in polar firn as a proxy for past temperatures. Geochimica et Cosmochimica Acta 225, 128–145.CrossRefGoogle Scholar
Howat, Ian M, Porter, Claire, Smith, Benjamin E., Noh, Myoung-Jong and Morin, Paul (2019) The Reference Elevation Model of Antarctica. The Cryosphere 13(2), 665–674. doi: 10.5194/tc-13-665-2019CrossRefGoogle Scholar
Jacobel, RW and Welch, BC (2005) Glaciological and climatic significance of Hercules Dome, Antarctica: an optimal site for deep ice core drilling. Journal of Geophysical Research 110, F01015.CrossRefGoogle Scholar
Jones, TJ and 7 others (2017) Improved methodologies for continuous-flow analysis of stable water isotopes in ice cores. Atmospheric Measurement Techniques 10(2), 617–632.CrossRefGoogle Scholar
Kahle, EC and 13 others (2021) Reconstruction of temperature, accumulation rate, and layer thinning from an ice core at south pole, using a statistical inverse method. Journal of Geophysical Research-Atmospheres 126(13), 1–20. doi: 10.1029/2020JD033300CrossRefGoogle Scholar
Kingslake, J and 9 others (2014) Full-depth englacial vertical ice sheet velocities measured using phase-sensitive radar. Journal of Geophysical Research-Earth Surface 119(12), 2604–2618.CrossRefGoogle Scholar
Korotkikh, EV and 7 others (2011) The Last Interglacial as represented in the glaciochemical record from Mount Moulton Blue Ice Area, West Antarctica. Quaternary Science Reviews 30(15–16), 1940–1947.CrossRefGoogle Scholar
Koutnik, MR and 7 others (2016) Holocene accumulation and ice flow near the West Antarctic ice sheet divide ice core site. Journal of Geophysical Research: Earth Surface 121, 907–924. doi: 10.1002/2015JF003668CrossRefGoogle Scholar
Koutnik, MR and Waddington, ED (2012) Well-posed boundary conditions for limited-domain models of transient ice flow near an ice divide. Journal of Glaciology 58(211), 1008–1020.CrossRefGoogle Scholar
Lee, JE and 16 others (2020) An 83 000-year-old ice core from Roosevelt Island, Ross Sea, Antarctica. Climate of the Past 16(5), 1691–1713.CrossRefGoogle Scholar
Lilien, DA, Hills, BH, Driscol, J, Jacobel, R and Christianson, K (2020) ImpDAR: an open-source impulse radar processor. Annals of Glaciology 61(81), 114–123.CrossRefGoogle Scholar
Lliboutry, L (1979) A critical review of analytical approximate solutions for steady state velocities and temperatures in cold ice-sheets. Zeitschrift fuer Gletscherkunde und Glazialgeologie 15, 135–148.Google Scholar
Looyenga, H (1965) Dielectric constants of homogeneous mixtures. Molecular Physics 9(6), 501.CrossRefGoogle Scholar
Martín, C, Gudmundsson, GH, Pritchard, HD and Gagliardini, O (2009) On the effects of anisotropic rheology on ice flow, internal structure, and the age-depth relationship at ice divides. Journal of Geophysical Research 114, 1–18. doi: 10.1029/2008JF001204CrossRefGoogle Scholar
Martín, C, Hindmarsh, RCA and Navarro, FJ (2006) Dating ice flow change near the flow divide at Roosevelt Island, Antarctica, by using a thermomechanical model to predict radar stratigraphy. Journal of Geophysical Research 11, 1–15. doi: 10.1029/2005JF000326Google Scholar
Martos, YM and 6 others (2017) Heat flux distribution of Antarctica unveiled. Geophysical Research Letters 44(22), 11417–11426.CrossRefGoogle Scholar
Matsuoka, K and 18 others (2015) Antarctic ice rises and rumples: their properties and significance for ice-sheet dynamics and evolution. Earth-Science Reviews 150, 724–745.CrossRefGoogle Scholar
Mercer, JH (1978) West Antarctic ice sheet and CO2 greenhouse effect: a threat of disaster. Nature 271, 321–325.CrossRefGoogle Scholar
Mulvaney, R and 9 others (2020) WACSWAIN project: does the preliminary timescale suggest we have captured the ice from the Last Interglacial period in the Skytrain ice core? AGU Fall Meeting, Abstracts, C033-01, https://agu.confex.com/agu/fm20/meetingapp.cgi/Paper/744866Google Scholar
Nereson, NA and Waddington, ED (2002) Isochrones and isotherms beneath migrating ice divides. Journal of Glaciology 48(160), 95–108.CrossRefGoogle Scholar
Nicholls, KW and 5 others (2015) A ground-based radar for measuring vertical strain rates and time-varying basal melt rates in ice sheets and shelves. Journal of Glaciology 61(230), 1079–1087.CrossRefGoogle Scholar
Nicolas, JP and Bromwich, DH (2011) Climate of West Antarctica and influence of marine air intrusions. Journal of Climate 24, 49–67. doi: 10.1175/2010JCLI3522.1CrossRefGoogle Scholar
Osterberg, EC, Handley, MJ, Sneed, SB, Mayewski, PA and Kreutz, KJ (2006) Continuous ice core meltersystem with discrete sampling for major ion, trace element, and stable isotope analyses. Environmental Science and Technology 40(1), 3355–3361.CrossRefGoogle ScholarPubMed
Pettit, EC, Thorsteinsson, T, Jacobson, HP and Waddington, ED (2007) The role of crystal fabric in flow near an ice divide. Journal of Glaciology 54(181), 277–288.CrossRefGoogle Scholar
Raymond, CF (1983) Deformation in the vicinity of ice divides. Journal of Glaciology 29(103), 357–373.CrossRefGoogle Scholar
Röthlisberger, R and 6 others (2000) Technique for continuous high-resolution analysis of trace substances in firn and ice cores. Environmental Science and Technology 34(2), 338–342.CrossRefGoogle Scholar
Scherer, RP and 5 others (1998) Pleistocene collapse of the West Antarctic ice sheet. Science 281(5373), 82–85.CrossRefGoogle ScholarPubMed
Shen, WS, Wiens, DA, Lloyd, AJ and Nyblade, AA (2020) A geothermal heat flux map of Antarctica empirically constrained by seismic structure. Geophysical Research Letters 47(14), 1–8. doi: 10.1029/2020GL086955CrossRefGoogle Scholar
Sigl, M and 26 others (2016) The WAIS Divide deep ice core WD2014 chronology – Part 2: annual-layer counting (0–31 ka BP). Climate of the Past 12, 769–786.CrossRefGoogle Scholar
Sodeman, H and Stohl, A (2009) Asymmetries in the moisture origin of Antarctic precipitation. Geophysical Research Letters 36, 1–5. doi: 10.1029/2009GL040242Google Scholar
Spector, P and 5 others (2018) West Antarctic sites for subglacial drilling to test for past ice-sheet collapse. The Cryosphere 12(8), 2741–2757.CrossRefGoogle Scholar
Stal, T, Reading, AM, Halpin, JA and Whittaker, JM (2021) Antarctic geothermal heat flow model: Aq1. Geochemistry, Geophysics, Geosystems 22(2), 1–22. doi: 10.1029/2020GC009428CrossRefGoogle Scholar
Steig, EJ and 7 others (2015) Influence of West Antarctic ice sheet collapse on Antarctic surface climate. Geophysical Research Letters 42(12), 4862–4868.CrossRefGoogle Scholar
Sutter, J, Gierz, P, Grosfeld, K, Thoma, M and Lohmann, G (2016) Ocean temperature thresholds for Last Interglacial West Antarctic Ice Sheet collapse. Geophysical Research Letters 43(6), 2675–2682. doi: 10.1002/grl.v43.6CrossRefGoogle Scholar
Veres, D and 15 others (2013) The Antarctic ice core chronology (AICC2012): an optimized multi-parameter and multi-site dating approach for the last 120 thousand years. Climate of the Past 9(4), 1733–1748.CrossRefGoogle Scholar
Waddington, ED, Neumann, TA, Koutnik, MR, Marshall, HP and Morse, DL (2007) Inference of accumulation-rate patterns from deep layers in glaciers and ice sheets. Journal of Glaciology 53(183), 694–712.CrossRefGoogle Scholar
Welch, BC and Jacobel, RW (2003) Analysis of deep-penetrating radar surveys of West Antarctica, US-ITASE 2001. Geophysical Research Letters 30(8), 1444.CrossRefGoogle Scholar
Welch, BC, Jacobel, RW and Arcone, SA (2009) First results from radar profiles collected along the US-ITASE traverse from Taylor Dome to South Pole (2006–2008). Annals of Glaciology 50(51), 35–41.CrossRefGoogle Scholar
Fig. 1. Overview map of the Hercules Dome geophysical survey in polar stereographic (PS) coordinates. The directions in the PS projection are used in the remainder of the paper. Arrow indicates direction of geographic north. The inset shows the two radar lines, x135 and x133, at West Dome. Surface elevation contours from the Reference Elevation Model of Antarctica (REMA, Howat and others, 2019) are shown at 10 m spacing for the overview and 2 m spacing for the inset.
Fig. 9. (a) Modeled internal layers (blue) assuming uniform flank flow (p = 4) compared to measured layers (black). The bed is shown in brown. (b) The difference between modeled and measured layers (blue circles) at 5.75 km (black vertical line) compared to the Raymond arch amplitude of Roosevelt Island (red squares) (Conway and others, 1999).
Fudge et al. supplementary material
File 1 MB
No CrossRef data available.
T. J. Fudge (a1), Benjamin H. Hills (a1) (a2), Annika N. Horlings (a1), Nicholas Holschuh (a3), John Erich Christian (a4) (a5), Lindsey Davidge (a1), Andrew Hoffman (a1), Gemma K. O'Connor (a1), Knut Christianson (a1) and Eric J. Steig (a1)
DOI: https://doi.org/10.1017/jog.2022.80